# Problem with solving differential equation using NDSolve

I am trying to solve differential equations subject to condition z[0] == 0 and z[tmax] == 10 Pi.

a = 1; k = 0.5; b = 1/(30*Pi); c = 3;
Λ = (a^2 (Cosh[3] - 1))/(Cosh[3] -  a^2 + (a^2 - 1) Cosh[z 3]);

L = NDSolve[{D[z[t],t] ==
(Sqrt[c + 1])/(Λ Sqrt[c Exp[-k Exp[b t] Sin[0.1 z[t]]] + Exp[-k Sin[0.1 z[t]]]]),
(z[0] == 0)}, z, {t, 0, 1000}];
SOLPLOT1[t_] := z[t] /. L;
Tend = Solve[SOLPLOT1[t] == 10 Pi , t][[1, 1, 2]];
P1 = Plot[SOLPLOT1[t], {t, 0, Tend}, PlotStyle -> {Black}]


However if I add accuracy goal to improve accuracy like below I get completely different picture

a = 1; k = 0.5; b = 1/(30*Pi); c = 3;
Λ = (a^2 (Cosh[3] - 1))/(Cosh[3] - a^2 + (a^2 - 1) Cosh[z 3]);

L = NDSolve[{D[z[t], t] ==
(Sqrt[c + 1])/(Λ Sqrt[c Exp[-k Exp[b t] Sin[0.1 z[t]]] +
Exp[-k Sin[0.1 z[t]]]]), (z[0] == 0)}, z, {t, 0, 1000},
AccuracyGoal -> ∞];
SOLPLOT1[t_] := z[t] /. L;
Tend = Solve[SOLPLOT1[t] == 10 Pi , t][[1, 1, 2]];
P1 = Plot[SOLPLOT1[t], {t, 0, Tend}, PlotStyle -> {Black}]


I tried different accuracy goal because if I choose different a for example a=1.3 will give following graph.

However this graph should not be the correct one. In case I increase accuracy goal it will give me a straight line. Is there any other way to improve solving such differential equation?

• Here should be z[t]? Λ = (a^2 (Cosh[3] - 1))/(Cosh[3] - a^2 + (a^2 - 1) Cosh[z 3]) – Alex Trounev Nov 1 '18 at 21:19
• a=1; k=1/2; b=1/(30*Pi); c=3; Λ=(a^2 (Cosh[3]-1))/(Cosh[3]-a^2+(a^2-1) Cosh[z[t] 3]); q=z[t]/.NDSolve[{D[z[t],t] == Sqrt[c+1]/(Λ Sqrt[c Exp[-k Exp[b t] Sin[1/10 z[t]]] + Exp[-k Sin[1/10 z[t]]]]), z[0]==0}, z[t], {t,0,1000}, AccuracyGoal->32][[1]]; Plot[{q,10 Pi}, {t,0,1000}, PlotRange->All] Solve[q==10 Pi, t] works just fine. Changing to a=13/11 blows up for small values of t. – Bill Nov 1 '18 at 22:15
• Thank you for your answer, the idea is to change a for different values > 1. – Alexander Nov 2 '18 at 6:30

Using WhenEvent, Rationalize, ParametricNDSolve there is no need to set AccuracyGoal:

k = 0.5; b = 1/(30*Pi); c = 3;\[CapitalLambda] = (a^2 (Cosh[3] -1))/(Cosh[3] -a^2 + (a^2 - 1) Cosh[z[t] 3]);
dgl = D[z[t],t] == (Sqrt[c + 1])/(\[CapitalLambda] Sqrt[c Exp[-k Exp[b t]Sin[0.1 z[t]]] + Exp[-k Sin[0.1 z[t]]]]) // Rationalize[#, 0] &
sol = ParametricNDSolveValue[{dgl, (z[0] == 0),WhenEvent[z [t] == 10    Pi , "StopIntegration"]},z, {t, 0, 1000}, {a} , PrecisionGoal -> 10 ,WorkingPrecision -> 100];


With Option WorkingPrecision now the Eventvalue 10Pi is reached.

For example a=11/3

te = sol[11/3]["Domain"][[1]][[2]];
sol[11/3][te] // Round[#, .0001] &
(* 31.4159 *)

Plot[sol[#][t], {t, 0, sol[#]["Domain"][[1]][[2]]}] &[11/3]


• Thank you for the answer, however z[tmax] must be equal 10Pi which is not possible to see from the above graph. – Alexander Nov 2 '18 at 15:49
• @Alexander Interesting, NDSolve stops without warning, so I conclude, the WhenEvent happened. But the solution doesn't include 10Pi ??? Perhaps you can ìmprove the results with option Method->"StiffnessSwitching" – Ulrich Neumann Nov 2 '18 at 16:18
• @Alexander I edited my answer, now the EventValue 10Pi is reached. – Ulrich Neumann Nov 2 '18 at 16:39
• Thank you for your answer – Alexander Nov 3 '18 at 6:43
• @Alexander You're welcome. Without success I tried to update the plot, but I don't know how to plot the WorkingPrecision->50 results ... – Ulrich Neumann Nov 3 '18 at 8:38