1
$\begingroup$

I need to get a phase portrait for a non-linear oscillator, for this I wrote down the corresponding equations. 3 equations for different "a" which can take values from minus to plus infinity. In each equation, the last term is an expression for the nonlinearity of our oscillator. In general, the problem (error message NDSolve::ndsz) arises when I try to solve equations for "a" which value is more than zero. For example, if i change last term for another like that: (a*x)/(a+x), it works. Can't anderstand where is my mistake in the first case? My code is:

δ = 0.1;
ω0 = 1;
α1 = 0.1;
α2 = 0;
α3 = -0.5;
Subscript[t, 0] = 0;
Subscript[t, 1] = 100;

x1 = DSolve[{x''[t] + 
  2δ*x'[t] + (ω0^2)*x[t] + α1/x[t]^2 == 0, 
x[0] == 1, x'[0] == 0}, x, {t, Subscript[t, 0], Subscript[t, 1]}];

x2 = DSolve[{x''[t] + 
  2δ*x'[t] + (ω0^2)*x[t] + α2/x[t]^2 == 0, 
x[0] == 1, x'[0] == 0}, x, {t, Subscript[t, 0], Subscript[t, 1]}];

x3 = DSolve[{x''[t] + 
  2δ*x'[t] + (ω0^2)*x[t] + α3/x[t]^2 == 0, 
x[0] == 1, x'[0] == 0}, x, {t, Subscript[t, 0], Subscript[t, 1]}];

Show[ParametricPlot[{x[t], x'[t]} /. x1, {t, 0, 60}, 
PlotRange -> {{-1, 1}, {-2, 2}}, PlotStyle -> Blue, 
PlotLegends -> {"α = 0.1"}], 
ParametricPlot[{x[t], x'[t]} /. x2, {t, 0, 60}, 
PlotRange -> {{-1, 1}, {-2, 2}}, PlotStyle -> Black, 
PlotLegends -> {"α = 0"}], 
ParametricPlot[{x[t], x'[t]} /. x3, {t, 0, 60}, 
PlotRange -> {{-1, 1}, {-2, 2}}, PlotStyle -> Red, 
PlotLegends -> {"α = -0.5"}]]
$\endgroup$
5
  • 1
    $\begingroup$ Only NDSolve[] and ilk can produce the NDSolve::ndsz error, yet you are using DSolve[] in your code. What am I missing? $\endgroup$ Commented Mar 21, 2018 at 23:52
  • 1
    $\begingroup$ You might try tracking down which command give the error and restrict your code to what is necessary. If DSolve doesn't work, try NDSolve, if a numeric solution would be acceptable. (It seems x1 and x3 fail.) $\endgroup$
    – Michael E2
    Commented Mar 21, 2018 at 23:53
  • 1
    $\begingroup$ @J.M. I think inside ParametricPlot, N[DSolve[..]] is tried, which now calls NDSolve, when DSolve fails. $\endgroup$
    – Michael E2
    Commented Mar 21, 2018 at 23:54
  • 1
    $\begingroup$ Yes, DSolve is unable to solve the x1 and x3 ODEs. Incidentally, because the ODEs are autonomous, they can be reduced to first order, although it is unclear whether that would help. Also, I would guess that, if NDSolve eventually is called, it would fail at the turning points of the nonlinear oscillators or when x[t] == 0. $\endgroup$
    – bbgodfrey
    Commented Mar 22, 2018 at 2:15
  • $\begingroup$ I'm sorry, I used NDsolve and it doesn't work too, than I tried to use Dsolve to see what changes and forgot to change it back $\endgroup$
    – John
    Commented Mar 22, 2018 at 8:31

1 Answer 1

2
$\begingroup$

Just use ParametricNDSolveValue to solve your problem:

sol = ParametricNDSolveValue[{x''[t] +2 \[Delta]*x'[t] + \[Omega]0^2*x[t] +\[Alpha] x[t]^2 == 0,x[0] == 1, x'[0] == 0}, x, {t, t0, t1}, {\[Alpha]}]
(*ParametricFunction[ <> ] *)

ParametricPlot[Evaluate[Table[{sol[\[Alpha]][t], 
sol[\[Alpha]]'[t]}, {\[Alpha], {\[Alpha]1, \[Alpha]2, \[Alpha]3}}]], {t, t0,t1}, PlotRange -> All]

enter image description here

$\endgroup$
4
  • $\begingroup$ Everything works, but you wrote last term as α*x[t]^2 but i need α/x[t]^2 and if change last it doesn't work. $\endgroup$
    – John
    Commented Mar 22, 2018 at 22:13
  • 1
    $\begingroup$ The problem seems to be the singularity x==0. If you set t1=1 everthing works! $\endgroup$ Commented Mar 23, 2018 at 7:06
  • 1
    $\begingroup$ Check the equilibrium of your ode 0 == \[Omega]0^2*x0 + \[Alpha] /x0^2: Only for \[Alpha] <0 a real solution x0>0exists . $\endgroup$ Commented Mar 23, 2018 at 7:22
  • $\begingroup$ Thank you for your answer! $\endgroup$
    – John
    Commented Mar 23, 2018 at 20:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.