Sorry but cannot fsolve this and don't know why
4 次查看(过去 30 天)
显示 更早的评论
Good evening MATLAB community,
I am conceding defeat after working at this for the last 3 days and I hope somebody can tell me what I have done incorrectly.
I am trying to extract two values from a set of data, the real and imaginary part of the index of refraction: N=n+ik=n1(1)+i*n1(2). The two equations describing the reflection and transmission coefficients are nonlinear in N so I thought using fsolve would be a good approach. The paper I am working from states that using the Newton method should result in solutions quickly; since fsolve is a variation on the Newton method I thought it was perfectly suited.
There are actually 3 materials in this system where the light goes from air into one medium and then into a second medium. These have indices n0, n1, and n2.
clear, clc
n0=[1 0]; %index of refraction for air/first layer
n10=[1.8 0];%initial guess for index of refraction for second layer
n2=[3.4 0]; %index of refraction for the third layer
lambda=10 %wavelength
d=2.3; %thickness of layer 2
r=0.297520661157025; %reflection coefficient for Si
t=0.702479338842975; %transmission coefficient for Si
next(w,:)=extract(r,t,n0,n2,d, lambda(w),n10); %call for fsolve function
My function extract is
function y = extract(r, t, n0, n2, d, lambda, n10)
y = fsolve(@call, n10);
function y = call(n1) ;
y = [
-(1+r)...
+t/(4*n0(1)*n2(1)*(n1(1)^2+n1(2)))*...
((n0(1)^2+n1(1)^2+n1(2)^2)*...
((n1(1)^2+n2(1)^2+n1(2)^2+n2(2)^2)*cosh(2*2*pi*n1(2)*d/lambda)+...
2*(n1(1)*n2(1)+n1(2)*n2(2))*sinh(2*2*pi*n1(2)*d/lambda))+...
(n0(1)^2-n1(1)^2-n1(2)^2)*...
((n1(1)^2-n2(1)^2+n1(2)^2-n2(2)^2)*cos(2*2*pi*n1(1)*d/lambda)+...
2*(n1(1)*n2(2)-n2(1)*n1(2))*sin(2*2*pi*n1(1)*d/lambda))),...
-(1-r)...
+t/(2*n2(1)*(n1(1)^2+n1(2)))*...
(n1(1)*...
((n1(1)^2+n2(1)^2+n1(2)^2+n2(2)^2)*sinh(2*2*pi*n1(2)*d/lambda)+...
2*(n1(1)*n2(1)+n1(2)*n2(2))*cosh(2*2*pi*n1(2)*d/lambda))+...
n1(2)*...
((n1(1)^2-n2(1)^2+n1(2)^2-n2(2)^2)*sin(2*2*pi*n1(1)*d/lambda)-...
2*(n1(1)*n2(2)-n2(1)*n1(2))*cos(2*2*pi*n1(1)*d/lambda)))
];
end
end
Now I have verified my equations are correct by solving for r and t for a dimple junction between air and silicon n0=[1 0] and n1=[3.4 0] and n2=[3.4 0]. That is how I obtained r and t in my code above and it follows the simple normal incidence Fresnel coefficients. So my simple test is to insert r and t obtained from the equations already back into the fsolve program to obtain n1=[3.4 0] but it does not work. The warning I get is:
Warning: Trust-region-dogleg algorithm of FSOLVE cannot handle non-square systems; using
Levenberg-Marquardt algorithm instead.
> In fsolve at 324
In extract at 4
In myfun_3 at 14
Optimizer appears to be converging to a minimum that is not a root:
Sum of squares of the function values exceeds the square root of
options.TolFun. Try again with a new starting point.
I have tried to use other peoples forum queries to solve my problem but I have not been able to fix it.
If anybody can take a look at it and give me some insight I would greatly appreciate it and I am sorry for asking something other people have asked about in the past.
Thank you all for your consideration,
Justin
0 个评论
回答(2 个)
Sargondjani
2012-5-18
fsolve can only find local solutions. like the exit message says: it finds a local minimum but it is not 0. So you should try a different starting value... Did you try this?
also, it seems that you can calculate the gradient analytically. providing the analytical gradient is always a good idea (it should be the second output argument of the objective function and then you set 'Jacobian' to 'on' in options)
0 个评论
Teja Muppirala
2012-5-18
Ok, you have a very subtle, but pathological problem with your code.
Compare this:
y = [ (1)...
+2]
To this:
y = [ (1)+...
2]
The first gives you y = [1 2], while the second gives you y = 3. This is exactly the problem in your code. As it is written, y will be a 4-element vector, when you really mean for it to be a 2-element vector. You are essentially trying to solve 4 equations with 2 unknowns, which is why you are getting that message.
This can be fixed by changing this:
-(1+r)...
to this:
-(1+r)+...
in the two places that it occurs. When I do that, FSOLVE converges instantly to
y = 2.173587481230821 -0.000000000000000
0 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Systems of Nonlinear Equations 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!