# Too slow using Dsolve

I am referring to the accepted answer to the question: Solving n simultaneous differential equation The issue is that the code is running too slow for even N=5! Is there any way to make it faster? I tried NDSolve as well. It’s also too slow! I am using the following code:

ClearAll["Global*"]
ClearAll[t, h, b, n, M];
NN = 5;
h = 21; b = 400;(*some made up values*)odes =
Table[I ToExpression["M" <> ToString[n]]'[t] ==
b Sqrt[NN + 3 + n]*ToExpression["M" <> ToString[n + 1]][t] +
h*Sqrt[n*(2*NN + 5)]*ToExpression["M" <> ToString[n - 1]][t], {n,
0, NN}]
deps = Table[ToExpression["M" <> ToString[n]][t], {n, 0, NN}]
M6[t_] := 2;(*some function for the last one,which has no ODE*)ic = \
{M0[0] == 1, M1[1] == 2, M2[0] == 3, M3[0] == 2, M4[0] == 1,
M5[0] == 4};(*some IC*)
NDSolve[{odes, ic}, deps, {t, 0, 3000}]
lot[{Evaluate[(M0[t]*Conjugate[M0[t]])] /. MH,
Evaluate[(M1[t]*Conjugate[M1[t]]) /. MH],
Evaluate[(M2[t]*Conjugate[M2[t]]) /. MH],
Evaluate[(M3[t]*Conjugate[M3[t]]) /. MH],
Evaluate[(M4[t]*Conjugate[M4[t]]) /. MH],
Evaluate[(M5[t]*Conjugate[M5[t]]) /. MH]}, {t, 0, 3000},
PlotRange -> All]

• I tried NDSolve as well. It’s also too slow! then you should post the code you tried. For me NDSolve finishes instantly for N=5 and for N=10 – Nasser May 21 at 4:02
• @Nasser Please see the edited question with code. – Jasmine May 21 at 4:48
• It only took 136 seconds for NDSolve to complete on my computer using V 12.2 and 12.3. This is not too slow considering your t goes up to 3000 now. Try Timing@NDSolve[{odes, ic}, deps, {t, 0, 3000}] and see what you get for the time used. screen shot !Mathematica graphics – Nasser May 21 at 5:34
• Jasmine, what do you mean by not getting any result? Can you, please, be more specific about what you expected to occur, and what it was that actually did occur? – CA Trevillian May 21 at 6:22
• In ic is it supposed to be M1[1] == 2 or M1[0] == 2 -- that is, is the problem supposed to be a BVP or an IVP? (In the current set up, it uses the shooting method to solve the BVP, which means it is solving many IVPs searching for the solution. That's why it's so slow, I think.) – Michael E2 May 21 at 15:16

Hint.

Try

sols = Solve[odes, D[deps, t]][[1]]
A = Grad[(D[deps, t] /. sols), deps]
MatrixExp[A t]


or also use

{Lambda,T} = Eigensystem[A]


then calling

udeps = Table[u[k][t],{k,1,6}];


we solve

deps = Inverse[T].udeps/.DSolve[Thread[D[udeps, t] == Lambda.udeps],udeps, t]
`

following with the initial conditions.