# How to implement the second loop

Here I have a problem that probably needs two loops, but I am not sure how to implement them together. The code calculates M for various values of parameter G. Now I want to find M which crosses the real axis (Imaginary part is zero) for each α1 = (1.1, 2.2). So the code is working fine and for example for α1 = 1.1, I got the curve which crosses the real axis, so I have one real solution (for chosen segment G), now I want to find these real points M for various α1 and to plot M[α1].

 \[Alpha]1 = 1.1;
\[Nu] = 0.001;
Tt = 0;
poly = -g^2 + 0.25 x^4 - (Tt + \[Alpha]1^2) x^2 + 2 \[Alpha]1 g x +
I \[Nu] g - I \[Nu] \[Alpha]1 x + 1;

LM = {};
Do[Lroots = x /. NSolve[(poly /. {g -> G}) == 0, x];
res = 0;
Do[root = Lroots[[i]];
If[Im[root] >= 0, L1 = Drop[Lroots, {i}];
res = res + 1/(Times @@ (root - L1));];, {i, 1, Length[Lroots]}];
M = 1/((4*I*G^2)*res);
AppendTo[LM, M];, {G, 0.01, 1.5, 0.001}];

a1 = Min@Select[Table[Im[LM[[i]]], {i, 1, Length[LM]}], Positive];
a2 = Max@Select[Table[Im[LM[[i]]], {i, 1, Length[LM]}], Negative];

b1 = Table[LM[[i]], {i, 1, Length[LM]}] - a1*I;
b2 = Table[LM[[i]], {i, 1, Length[LM]}] - a2*I;

c1 = Select[b1, Im[#] == 0 &][];
c2 = Select[b2, Im[#] == 0 &][];

d1=Re[(c1 + c2)/2]


The curve crosses the real axis, but just for α1 = 1.1; now I need to extract for various α1 = 1.1, ..., 2.2, all points which cross the real axis.

 p1 = ListPlot[Table[{Re[LM[[i]]], Im[LM[[i]]]}, {i, 1, Length[LM]}],
PlotRange -> {{-1, 18}, {-10, 10}}]

• Yes Öskå, from the code I need just one point when the Imaginary part is zero and then to calculate it for various alpha – Pipe May 6 '14 at 11:58
• not this point, just one which cross the real axis, if doesn't exist it is infinity – Pipe May 6 '14 at 12:02
• Yes and maybe I should change the G and change the step to see the crossing real axis – Pipe May 6 '14 at 12:04
• Öskå please take a look now. I made a mistake because of segment. Now for each alpha the plotted curve cross the real axis, I need to extract these points for each alpha and to plot the M(alpha) – Pipe May 6 '14 at 17:26
• Öskå it is not important to be exact, this is because of chosen step, it is enough just to choose the point with smallest imaginary part and for this point to extract the real value for each alpha on the segment. For all values alpha which I changed the curve cross the real axis for sure, but I just need to do this automatically, not to calculate for each alpha manually – Pipe May 6 '14 at 22:36

Tell me if the following fits you. Note that I haven't tried to optimize your code. The idea is simply to put everything under the same function getNearZero which depends on α1. Then you are free to do every operations you want and do your second loop over α1.

getNearZero[α1_] :=
Module[{ν, Tt, Poly, g, G, x, M, LM, poly, res, root, Lroots, L1},
ν = 0.001; Tt = 0;
poly = -g^2 + 0.25 x^4 - (Tt + α1^2) x^2 + 2 α1 g x + I ν g - I ν α1 x + 1;
LM = {};
Do[Lroots = x /. NSolve[(poly /. {g -> G}) == 0, x];
res = 0;
Do[root = Lroots[[i]];
If[Im[root] >= 0, L1 = Drop[Lroots, {i}];
res = res + 1/(Times @@ (root - L1));];, {i, 1, Length[Lroots]}];
M = 1/((4*I*G^2)*res);
AppendTo[LM, M];, {G, 0.01, 1.5, 0.001}];
(* return the Nearest to zero and LM *)
{Extract[LM, Position[Im /@ LM, First@Nearest[Im /@ LM, 0]]], LM}]


Checking if that's the closest point to zero:

With[{LM = Last@getNearZero@1.3, nearZero = First@getNearZero@1.3},
ListPlot[Table[{Re[LM[[i]]], Im[LM[[i]]]}, {i, 1, Length[LM]}],
PlotRange -> {{-1, 18}, {-10, 10}},
Epilog -> {PointSize@0.02, Red, Point[{First@Re@nearZero, First@Im@nearZero}]}]] Then you can run the getNearZero from 1.1 to 2.2 every di = 0.1 and ListLinePlot the result (note the Tooltip showing the corresponding α1).

data = {#, First@getNearZero@#} & /@ Range[1.1, 2.2, .1];
ListLinePlot[
Tooltip[#, "α1=" <> ToString@#2] & @@@
Frame -> True, PlotMarkers -> Automatic] If you want to plot the real part of the points nearest to zero in terms of α1:

ListLinePlot[Thread[{First /@ data, Flatten[Re /@ Last /@ data]}],
PlotMarkers -> Automatic, AxesLabel -> {"α1"}] And finally, if you want to work with your newly edited ai, bi, ci and d1 the following is shorter:

d1 = First @ Mean[Re /@ Extract[LM, Position[Im /@ LM, #]] & /@ Nearest[Im /@ LM, 0, 2]]


where Nearest[Im /@ LM, 0, 2] replaces your ai, the bi are not needed anymore, and Re /@ Extract[LM, Position[Im /@ LM, #]] & replaces your ci.

• Thank you very much Öskå. This is it. I should learn through your code – Pipe May 7 '14 at 10:20