Solving Highly Irregular Non-Linear Equations

Hey all, I'm trying to solve a system of 3 highly irregular nonlinear equations using Matlab, but the solutions I get are not the expected values. Here are the equations in the most basic form, and the variables I am currently using:

n0=1;
n21=3.861;
n22=3.8593;
n23=3.8576;
k21=.015;
k22=.015;
k23=.015;
y1=640;
y2=641;
y3=642;
R1=.11175;
R2=.112907;
R3=.114064;

fcns(1)=((((n0-(n-1i.*k))./(n0+(n-1i.*k)))+(((n-1i.*k)-(n21-1i.*k21))./((n-1i.*k)+(n21-1i.*k21))).*(exp(-2.*1i.*(2.*pi./y1).*(n-1i.*k).*d)))./(1+((n0-(n-1i.*k))./(n0+(n-1i.*k))).*(((n-1i.*k)-(n21-1i.*k21))./((n-1i.*k)+(n21-1i.*k21))).*(exp(-2.*1i.*(2.*pi./y1).*(n-1i.*k).*d))))-R1;

fcns(2)=((((n0-(n-1i.*k))./(n0+(n-1i.*k)))+(((n-1i.*k)-(n22-1i.*k22))./((n-1i.*k)+(n22-1i.*k22))).*(exp(-2.*1i.*(2.*pi./y2).*(n-1i.*k).*d)))./(1+((n0-(n-1i.*k))./(n0+(n-1i.*k))).*(((n-1i.*k)-(n22-1i.*k22))./((n-1i.*k)+(n22-1i.*k22))).*(exp(-2.*1i.*(2.*pi./y2).*(n-1i.*k).*d))))-R2;

fcns(3)=((((n0-(n-1i.*k))./(n0+(n-1i.*k)))+(((n-1i.*k)-(n23-1i.*k23))./((n-1i.*k)+(n23-1i.*k23))).*(exp(-2.*1i.*(2.*pi./y3).*(n-1i.*k).*d)))./(1+((n0-(n-1i.*k))./(n0+(n-1i.*k))).*(((n-1i.*k)-(n23-1i.*k23))./((n-1i.*k)+(n23-1i.*k23))).*(exp(-2.*1i.*(2.*pi./y3).*(n-1i.*k).*d))))-R3;

The answers are never close to the expected values, which are n=1.5, k=0, and d=300. I'm using the fsolve command to perform the calculations. Any tips or suggestions?

Answer 1

[Walter Roberson]

I used the equivalent expressions in Maple,

fcns[1] := ((n0-n+I*k)/(n0+n-I*k)+(n-I*k-n21+I*k21)*exp(-2*(2*I)*Pi*(n-I*k)*d/y1)/(n-I*k+n21-I*k21))/(1+(n0-n+I*k)*(n-I*k-n21+I*k21)*exp(-2*(2*I)*Pi*(n-I*k)*d/y1)/((n0+n-I*k)*(n-I*k+n21-I*k21)))-R1;
fcns[2] := ((n0-n+I*k)/(n0+n-I*k)+(n-I*k-n22+I*k22)*exp(-2*(2*I)*Pi*(n-I*k)*d/y2)/(n-I*k+n22-I*k22))/(1+(n0-n+I*k)*(n-I*k-n22+I*k22)*exp(-2*(2*I)*Pi*(n-I*k)*d/y2)/((n0+n-I*k)*(n-I*k+n22-I*k22)))-R2;
fcns[3] := ((n0-n+I*k)/(n0+n-I*k)+(n-I*k-n23+I*k23)*exp(-2*(2*I)*Pi*(n-I*k)*d/y3)/(n-I*k+n23-I*k23))/(1+(n0-n+I*k)*(n-I*k-n23+I*k23)*exp(-2*(2*I)*Pi*(n-I*k)*d/y3)/((n0+n-I*k)*(n-I*k+n23-I*k23)))-R3;

When I substitute your proposed n, k, d values in to that, the expressions come out non zero and not equal:

[0.5428451377e-1+.2713792967*I, 0.4862923664e-1+.2764299574*I, 0.4291778332e-1+.2813827716*I]

If I did not make an editing error in doing the conversion (and that is a possibility) then your proposed values are not solutions ?

[Alex]

That is a possibility. Since this is a test of the program that I am creating for such an analysis of the data, all the values being used are measured and known values representing the physical properties of the materials being tested. The equations being used happen to be a Fresnel representation of the relationships between three materials using reflectance. We generate 3 equations by analyzing the equation at various wavelengths of light, and we want to approximate the valued of n, k, and d with this method, which are known for the 3 wavelengths. The problem as I see it could arise from 3 possible issues (there may be more, of course):

1) The data for n2 and k2 may be wrong. This could lead to a solution that doesn't match the rest of the data, causing the discrepancy.

2) The estimation of a single value for the three actual values using this method might not be a valid result for the data.

3) Matlab simply cannot compute this.

If it helps in trouble shooting, when I change any of the values in the guess for the fsolve, I get really different answers, and I don't think that should happen. If I alter the guess, shouldn't that just affect the number of iterations before a zero is reached?

Here is a fresh data set if you want to try again. I will include the actual values of n, k, and d for this set for each wavelength:

y,n2,k2,n0,k0,R,n,k,d

500,4.293,.045,1,0,.348033,1.548756,0,300

501,4.2876,.0444,1,0,.344843,1.548689,0,300

502,4.2822,.0438,1,0,.341654,1.548622,0,300

[Walter Roberson]

With those values, your equations cannot produce anywhere near the desired result.

You have given 3 lines of values, with a header line that does not directly correspond to many of the variables in the equations. For example there is no "y" in the equation, so what are we to do with the y value from the table? The seemingly obvious thing to do would be to guess that the first line corresponds to y1, the second line to y2, the third to y3, and so on for the other values. But then we get to the variables "n", "k" and "d", with there being no k1, k2, k3, d1, d2, d3. But the three lines have identical k and d values so we can speculate that we should use that constant value for k, and the constant value for d. But now what do we do with "n" ? The values of "n" are different for the three lines and there is no series n1/n2/n3, and your search is looking for "n" and "k" and "d". So we hold all three n, k, and d as being variables to be solved for and hope we get the values from the table out... a plausibility for k and d, but only one single "n" value can be calculated at best, and there are three different "n" values in the table.

If you do the variable splitting as I describe and use n=1.548756 from the first line of the table, and use k = 0, and then individually solve each of the equations for d to determine whether you can get something "close to" the expected 300, the answers for k that come out are values that are not so far apart -- but they are all complex values, and all have magnitude on the order of 80.5 which is a long way from 300.

If you substitute in k=0 and d=300 and try to solve the individual equations for n, then you get complex values. In the _approximate_ area of the solution you are looking for, at roughly n = 1.25, there is a point at which the imaginary part is 0... but then the real part of the expression is non-zero which should not be the case if n was a solution to the equation. At least it is positive there, though.

Roughly speaking, the imaginary portion for n follows the pattern of a sine wave of increasing frequency, but soon after n=1.3, the real portion goes negative and becomes a decaying sine-like wave approximately centered around -0.7, probably never crossing 0 again.

Your equations appear to be, to use the technical term, "borked".

Answer 2

[Siamak Lotfi]

hello

if its does work , you may consider give a try with blocks in Simulink

[Walter Roberson]

Simulink cannot calculate things that have no solution. Simulink could at most be used to numerically find a point of local minima of error relative to the expected values, which is something that could be done directly in MATLAB if that were the task.