0
$\begingroup$

I have to solve a set of coupled differential equations. First I have created a table containing all the equations and then solved using NDsolve.

 equations = Table[
  I D[M1n[[i]][[j]][[k]], t] == 
   Part[(m2/(b[[i]]*(*10**)10^6*1.98*10^-10) + 
        10^5/t^4*Sum[If[l == i, 0, M1n[[l]]], {l, 1, 3, 1}]).M1n[[
       i]] - M1n[[
       i]].(m2/(b[[i]]*(*10**)10^6*1.98*10^-10) + 
        10^5/t^4*Sum[If[l == i, 0, M1n[[l]]], {l, 1, 3, 1}]), j, 
    k], {i, 1, 3}, {j, 1, 2}, {k, 1, 2}]
equations1 = Flatten[equations];
s1 = NDSolve[{equations1, ci, di, 
   ei}, {Flatten[{Variables[M1n[[1]]], Variables[M1n[[2]]], 
     Variables[M1n[[3]]]}]}, {t, 2.2, 10^4},
  PrecisionGoal -> 5, WorkingPrecision -> 7]

where,

  b = Array[En, 3];
En[1] = 10;
En[2] = 11;
En[3] = 15;
M1n = Table[
  Array[Subscript[Subscript[\[Rho], i], #1, #2][t] &, {2, 2}], {i, 1, 
   3}]
\[Theta] =11.53*\[Pi]/180;
m2 = 1/(2*2)*2*10^-3*{{Cos[
     2 \[Theta]], -Sin[2 \[Theta]]}, {-Sin[2 \[Theta]], -Cos[
      2 \[Theta]]}}
M0 = {{10, 10^-5}, {10^-5, 5}};
ci = Thread[Flatten[M1n[[1]]] == Flatten[M0]] /. {t -> 2};
M20 = {{10, 10^-5}, {10^-5, 4}};
di = Thread[Flatten[M1n[[2]]] == Flatten[M20]] /. {t -> 2};
M30 = {{10, 10^-5}, {10^-5, 3}};
ei = Thread[Flatten[M1n[[3]]] == Flatten[M30]] /. {t -> 2};

But this is showing an error

NDSolve::precw: The precision of the differential equation ({{I (TemporaryVariable$13055^\[Prime])[t]==TemporaryVariable$13057[t] (-0.0989128+100000 Power[<<2>>] Plus[<<2>>])-TemporaryVariable$13056[t] (-0.0989128+100000 Power[<<2>>] Plus[<<2>>]),I (TemporaryVariable$13056^\[Prime])[t]==TemporaryVariable$13056[t] (0.232347 +100000 Power[<<2>>] Plus[<<2>>])-TemporaryVariable$13055[t] (-0.0989128+100000 Power[<<2>>] Plus[<<2>>])+TemporaryVariable$13058[t] (-0.0989128+100000 Power[<<2>>] Plus[<<2>>])-TemporaryVariable$13056[t] (-0.232347+100000 Power[<<2>>] Plus[<<2>>]),<<20>>,TemporaryVariable$13065[2]==1/100000,TemporaryVariable$13066[2]==3},{},{},{},{}}) is less than WorkingPrecision (7.`). >>
NDSolve::ndsz: At t == 2.`7., step size is effectively zero; singularity or stiff system suspected. >>
NDSolve::ndsz: At t == 2.`7., step size is effectively zero; singularity or stiff system suspected. >>

How to resolve this?

$\endgroup$
3
  • 2
    $\begingroup$ Delete: " WorkingPrecision -> 7" and add: "MaxSteps -> Infinity" then your code will execute. $\endgroup$ Sep 29 '20 at 8:58
  • $\begingroup$ Yes I tried that but now it is showing $\endgroup$ Sep 29 '20 at 11:19
  • $\begingroup$ NDSolve::dsfun: {Subscript[Subscript[[Rho], 1], 1,1][t],Subscript[Subscript[[Rho], 1], 1,2][t],Subscript[Subscript[[Rho], 1], 2,1][t],Subscript[Subscript[[Rho], 1], 2,2][t],Subscript[Subscript[[Rho], 2], 1,1][t],Subscript[Subscript[[Rho], 2], 1,2][t],Subscript[Subscript[[Rho], 2], 2,1][t],Subscript[Subscript[[Rho], 2], 2,2][t],Subscript[Subscript[[Rho], 3], 1,1][t],Subscript[Subscript[[Rho], 3], 1,2][t],Subscript[Subscript[[Rho], 3], 2,1][t],Subscript[Subscript[[Rho], 3], 2,2][t]} cannot be used as a function. >> $\endgroup$ Sep 29 '20 at 11:21
1
$\begingroup$

Maybe we have a version problem. On Mathematica 12.1 I get 12 numerical functions, some of them with complex output. I post a plot of the first one. The rest with real output look similar. enter image description here

$\endgroup$
10
  • $\begingroup$ Ok. I am using Mathematica 10. $\endgroup$ Sep 30 '20 at 4:03
  • $\begingroup$ I am using Mathematica 10.4 $\endgroup$ Sep 30 '20 at 5:01
  • $\begingroup$ As you see, the result is highly oscillatory and it is well possible that version 10 is not able to do it. But I can not claim this with certainty. $\endgroup$ Sep 30 '20 at 10:46
  • $\begingroup$ I have solved the same set of equations by writing each equation separately and it is giving some solution. $\endgroup$ Oct 1 '20 at 17:30
  • $\begingroup$ s = NDSolve[{ I D[M1n[[1]], t] == (m2/(10*10^6*1.98*10^-10) + 10^5/t^4*2*(M1n[[2]] + M1n[[3]])).M1n[[1]] - M1n[[1]].(m2/(10*10^6*1.98*10^-10) + 10^5/t^4*2*(M1n[[2]] + M1n[[3]])), I D[M1n[[2]], t] == (m2/(11*10^6*1.98*10^-10) + 10^5/t^4*2*(M1n[[1]] + M1n[[3]])).M1n[[2]] - M1n[[2]].(m2/(11*10^6*1.98*10^-10) + 10^5/t^4*2*(M1n[[1]] + M1n[[3]])), I D[M1n[[3]], t] == (m2/(15*10^6*1.98*10^-10) + $\endgroup$ Oct 1 '20 at 17:38

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.