0
$\begingroup$

Some sample data

data0 = {{0.00005, -0.99994}, {0.0001, -0.999882}, {0.0004, -0.99953},
         {0.0104, -0.988117}, {0.0204, -1.00472}, {0.0304, -1.01358},
         {0.0404, -1.01729}, {0.0504, -1.01776}};
NICs = Length[data0];

and some definitions of derivatives

x1 = m3*Sqrt[3];
y1 = 0;
x2 = Sqrt[3]/2*(2*m3 - 1);
y2 = 1/2;
x3 = x2;
y3 = -y2;

r1 = Sqrt[(x - x1)^2 + (y - y1)^2];
r2 = Sqrt[(x - x2)^2 + (y - y2)^2];
r3 = Sqrt[(x - x3)^2 + (y - y3)^2];

Ω = m1/r1 + m2/r2 + m3/r3 + 1/2*(x^2 + y^2);

m2 = m3;
m1 = 1 - 2*m3;

Ωx = D[Ω, x];
Ωy = D[Ω, y];

Ωxx = D[Ωx, x]; 
Ωxy = D[Ωx, y];  
Ωyx = D[Ωy, x];
Ωyy = D[Ωy, y];

Now let's setup a FOR loop

data = {};
For[j = 1, j <= NICs, j++, 
    m30 = data0[[j, 1]];
    x0 = data0[[j, 2]]; 
    eq = λ^4 + (4 - Ωxx - Ωyy)*λ^2 + Ωxx*Ωyy - Ωxy^2 /. {m3 -> m30, 
         x -> x0, y -> 0};
    sol = Solve[eq == 0, λ];
    AppendTo[data, {m30, x0, k}];
   ]

As we can see, inside the new list data there is a undefined parameter $k$. The value of this parameter is related with the solutions of the equation eq==0. This equation has always four roots. Now, if all four of them are purely imaginary (their real part is equal to zero) then $k = 1$, otherwise $k = -1$.

So my question is: how can I add a condition inside the FOR loop so as to have either $k = 1$ or $k = -1$ in data, according to the particular initial conditions?

The new list should be as follows

data = {{0.00005, -0.99994, 1}, {0.0001, -0.999882, 1}, 
        {0.0004, -0.9995295, 1}, {0.0104, -0.9881174, -1},  
        {0.0204, -1.0047199, -1}, {0.0304, -1.0135774, -1}, 
        {0.0404, -1.0172852, -1}, {0.0504, -1.0177593, -1}}

EDIT

Following @george2079 suggestion we can define $k$ as

k = If[Re[λ /. sol[[1]]] == 0 && Re[λ /. sol[[2]]] == 0 && 
       Re[λ /. sol[[3]]] == 0 && Re[λ /. sol[[4]]] == 0, 1, -1];

Is there any other, more elegant solution? Many thanks in advance!

$\endgroup$
  • 1
    $\begingroup$ k=If[condition,1,-1] $\endgroup$ – george2079 Sep 23 '16 at 11:20
  • $\begingroup$ @george2079 Is there a more elegant way of writing k = If[Re[λ /. sol[[1]]] == 0 && Re[λ /. sol[[2]]] == 0 && Re[λ /. sol[[3]]] == 0 && Re[λ /. sol[[4]]] == 0, 1, -1];? $\endgroup$ – Vaggelis_Z Sep 23 '16 at 11:32
  • $\begingroup$ you can use AllTrue or VectorQ $\endgroup$ – george2079 Sep 23 '16 at 11:48
  • 1
    $\begingroup$ a bit aside but i think you can get that result more elegantly if you write the expression in terms of I lambda and check the discriminant. $\endgroup$ – george2079 Sep 23 '16 at 12:03
  • $\begingroup$ @george2079 Could you post a quick answer so as to accept it? Perhaps showing the use of VectorQ... $\endgroup$ – Vaggelis_Z Sep 23 '16 at 12:06
1
$\begingroup$

using VectorQ :

sol = \[Lambda] /. Solve[eq == 0, \[Lambda]];
k = If[VectorQ[sol, Chop[Re[#] == 0] &], 1, -1];
AppendTo[data, {m30, x0, k}];

I added Chop in case some values might be machine precision small but not quite exactly zero.

note you really should use Map here instead of For, it will be much faster and cleaner.

data = (
     m30 = #[[1]];
     x0 = #[[2]];
     eq = \[Lambda]^4 + (4 - \[CapitalOmega]xx - \[CapitalOmega]yy)*\
     \[Lambda]^2 + \[CapitalOmega]xx*\[CapitalOmega]yy - \
     \[CapitalOmega]xy^2 /. {m3 -> m30, x -> x0, y -> 0};
     sol = \[Lambda] /. Solve[eq == 0, \[Lambda]];
     k = If[VectorQ[sol, Chop[Re[#] == 0] &], 1, -1];
     {m30, x0, k}) & /@ data0;

Note this can be written even more compactly but I left it looking like your For construct to illustrate Map

eq = \[Lambda]^4 + (4 - \[CapitalOmega]xx - \[CapitalOmega]yy)*\
      \[Lambda]^2 + \[CapitalOmega]xx*\[CapitalOmega]yy - \
       \[CapitalOmega]xy^2;
data =
  {#[[1]], #[[2]],
     If[VectorQ[\[Lambda] /.
        Solve[(eq /.
            {m3 -> #[[1]], x -> #[[2]], y -> 0}) == 0, \[Lambda]],
       Chop[Re[#] == 0] &], 1, -1]} & /@ data0;
$\endgroup$
  • $\begingroup$ I'm not familiar with Map. How can Map replace the existing For loop? $\endgroup$ – Vaggelis_Z Sep 23 '16 at 15:50
  • $\begingroup$ You are right! When the original data0 list contains about 10000 initial conditions the compact forms are much more faster then the traditional For loop! Many thanks again! $\endgroup$ – Vaggelis_Z Sep 23 '16 at 16:18

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.