FindRoot::lstol error

Question

I'm getting the message:

FindRoot::lstol: The line search decreased the step size to within tolerance specified by AccuracyGoal and PrecisionGoal but was unable to find a sufficient decrease in the merit function. You may need more than MachinePrecision digits of working precision to meet these tolerances.

I am also getting results and the plot that I get is acceptable, but I want to know: does the problem described in the message influence the results? Or can I continue?

 e1=1; e2=-1; e3=3.0625; e4=-1; e5=1; lam=632.8*10^-9; c=3*10^8; w=2*3.14*c/lam; m1=1; m2=-1; m3=1; m4=-1; m5=1; d2=60*10^-9; d4=60*10^-9; f=0.5; j=2;
        k0=w/c;
m = 0;
Wi = 150 ;                 Wf = 500;      Wstep = 10;                      

nf = (Wf - Wi)/Wstep + 1

        q3 = Sqrt[(k0^2*e3*m3) - B^2];
        q1 = Sqrt[B^2 - (k0^2*e1*m1)];
    q2 = Sqrt[B^2 - (k0^2*e2*m2)];
    q4 = Sqrt[B^2 - (k0^2*e4*m4)];
    q5 = Sqrt[B^2 - (k0^2*e5*m5)];
    wValues = Table[i, {i, Wi, Wf, Wstep}];
    R1 = (1 - (((q1*m2*k0)/(q2*m1))*f))*Exp[-q2*d2];
    R2 = (1 + (((q1*m2*k0)/(q2*m1))*f))*Exp[q2*d2];
    R3 = (1 - (((q4*m5)/(q5*m4*k0))*j))*Exp[-q4*d4];
    R4 = (1 + (((q4*m5)/(q5*m4*k0))*j))*Exp[q4*d4];
    R = q3*wValues*10^(-9);
    X = ((q2*m3)/(q3*m2))*((R2 - R1)/(R1 + R2));
    Y = ((q4*m3)/(q3*m4))*((R3 + R4)/(R4 - R3));
    disp = R - ArcTan[X] - ArcTan[Y] - m*Pi;

    n = Table[j, {j, 1, nf, 1}];
    BValues = 
     Table[FindRoot[disp[[n]], {B, 1.5578*10^7}, 
       MaxIterations -> 10^5], {n, 1, nf, 1}]
    BValues = Re[B /. BValues]
    Nvalues = BValues/k0

I tried changing the starting value, but the error still appears. The values that I get from programme are complex. Is it necessary for the starting value that I choose be complex?

Answer 1

This computation fails, because there are no solutions for the parameters chosen. An additional complication is the branch cuts that occupy much of the real axis. Both can be demonstrated by plotting a typical disp expression (with B renormalized by k0 for convenience).

ComplexPlot3D[disp[[nf/2]], {B, -4 - I, 4 + I}]

Incidentally, the renormalization is accomplished by redefining q1 as

q1 = k0 Sqrt[B^2 - (e1*m1)];

and similarly for {q2, q3, q4, q5}. For future reference, applying Simplify to disp simplifies these expressions considerably, especially if 3.14 is replaced by the exact expression Pi in the definition of w.

Addendum: Solutions for modified parameters

An alternative way to display the disp functions is to plot them in the region between the two branch points that occur at 1 and 7/4.

Plot[disp, {B, 1, 7/4}]

For none of the nf expressions does a zero occur in the region. However, changing some of the parameters may yield solutions. For instance, with m = 2, the curves decrease in value by 2 Pi, and those for n >= 18 cross the axis.

Table[{n, FindRoot[disp[[n]], {B, 1.5578*10^7/k0}] // Values}, {n, 18, nf}]

(* {{18, {1.00403}}, {19, {1.08856}}, {20, {1.15744}}, {21, {1.21482}}, 
    {22, {1.26341}}, {23, {1.30511}}, {24, {1.34127}}, {25, {1.3729}}, 
    {26, {1.40079}}, {27, {1.42552}}, {28, {1.4476}}, {29, {1.46739}}, 
    {30, {1.48523}}, {31, {1.50137}}, {32, {1.51603}}, {33, {1.52939}}, 
    {34, {1.54161}}, {35, {1.55281}}, {36, {1.5631}}} *)

Varying other parameters also may yield solutions (and also may move the branch points).

At the request of the OP, here is the entire code.

e1 = 1; e2 = -1; e3 = 3.0625; e4 = -1; e5 = 1; lam = 632.8*10^-9; 
c = 3*10^8; w = 2*Pi*c/lam; m1 = 1; m2 = -1; m3 = 1; m4 = -1; m5 = 1; d2 = 60*10^-9; 
d4 = 60*10^-9; f = 1/2; j = 2; k0 = w/c; m = 2; Wi = 150; Wf = 500; Wstep = 10;

nf = (Wf - Wi)/Wstep + 1

q3 = k0 Sqrt[(e3*m3) - B^2];
q1 = k0 Sqrt[B^2 - (e1*m1)];
q2 = k0 Sqrt[B^2 - (e2*m2)];
q4 = k0 Sqrt[B^2 - (e4*m4)];
q5 = k0 Sqrt[B^2 - (e5*m5)];
wValues = Table[i, {i, Wi, Wf, Wstep}];
R1 = (1 - (((q1*m2*k0)/(q2*m1))*f))*Exp[-q2*d2];
R2 = (1 + (((q1*m2*k0)/(q2*m1))*f))*Exp[q2*d2];
R3 = (1 - (((q4*m5)/(q5*m4*k0))*j))*Exp[-q4*d4];
R4 = (1 + (((q4*m5)/(q5*m4*k0))*j))*Exp[q4*d4];
R = q3*wValues*10^(-9);
X = ((q2*m3)/(q3*m2))*((R2 - R1)/(R1 + R2));
Y = ((q4*m3)/(q3*m4))*((R3 + R4)/(R4 - R3));
disp = Simplify[R - ArcTan[X] - ArcTan[Y] - m*Pi];

Table[{n, FindRoot[disp[[n]], {B, 1.5578*10^7/k0}] // Values}, {n, 18, nf}]