1
$\begingroup$

I have got a differntial equation which I have tried to solve with DSolve

sol = DSolve[p'[t] - (2 p[t]^2) - 3 p[t] + 4 p[t]^-9 == 0, p, t]

The answer implies RootSum. The question is how to use and plot the solution. I have tried

Plot[Evaluate[p[t] /. sol /. {C[1] -> 1}], {t, -7, 7}]

which doesn't work. Moreover, this case is not documented in the documentation "Plotting the Solution".

Need help; Thanks

$\endgroup$
7
  • $\begingroup$ Do you necessarily need an analytical solution ? It may be best to go for NDSolve. $\endgroup$
    – Lotus
    Nov 3 '16 at 11:24
  • $\begingroup$ When using NDSolve[] it appears that for any initial condition, the system quickly reaches a singularity. $\endgroup$
    – Feyre
    Nov 3 '16 at 11:26
  • $\begingroup$ @Feyre: Not when I tried to solve it in the range {0,7} $\endgroup$
    – Lotus
    Nov 3 '16 at 11:35
  • $\begingroup$ I am getting some kind of solution and also a plot of it but I am unable to make anything out if it specifically the cone like shape in the plot, NDSolve[{p'[t] - (2 p[t]^2) - 3 p[t] + 4 p[t]^-9 == 0, p'[0] == 2}, p, {t, 0, 10}] $\endgroup$ Nov 3 '16 at 11:40
  • $\begingroup$ @Feyre and @Lotus, The plot code I have used is Plot[Evaluate[p[t] /. %], {t, 0, 10}] but I am unable to make anything out if it specifically the cone like shape in the plot. $\endgroup$ Nov 3 '16 at 11:51
2
$\begingroup$

The problem is not RootSum but Solve:

sol = DSolve[p'[t] - (2 p[t]^2) - 3 p[t] + 4 p[t]^-9 == 0, p, t]
(*
  Solve[1/2 RootSum[-4 + 3 #1^10 + 2 #1^11 &, 
      Log[p[t] - #1]/(15 + 11 #1) &] == t + C[1], p[t]]
*)

Workarounds for plotting include replacing Solve with FindRoot or # &, or replacing the variable to be solved for by t:

Block[{Solve = fr, fr, t},
 fr[eq_, x_] /; NumericQ[t] :=(*x/.*)FindRoot[eq, {x, 1}];
 Plot[Evaluate@Quiet[p[t] /. (sol /. {C[1] -> 1})], {t, -1.45, -1}, 
  PlotRange -> {0, 50}]
 ]

Mathematica graphics

Block[{Solve = # &, p},
 ContourPlot[sol /. {p[t] -> p, C[1] -> 1} // Evaluate,
  {t, -1.4, -1}, {p, 0, 50}, Exclusions -> None]
 ]

(* plot like that above *)

Solving for t gives us the inverse function of the solution, which we can plot with ParametricPlot:

Block[{eq, t},
 eq = Block[{Solve = #1 &}, sol];
 With[{tt = t /. First@Solve[eq /. {C[1] -> 1, p[t] -> p}, t]},
  ParametricPlot[{tt, p}, {p, 0, 50}, PlotRange -> All, 
   AspectRatio -> 0.6, Exclusions -> None]
  ]]

(* plot like that above *)

Update: Alternatives for visualizing the solutions of a first-order ODE

One can plot the direction field with StreamPlot. Here is a utility to construct a direction field.

Clear[dirfield];
dirfield[de_, v_, t_] := Module[{y, df},
   y = First@Flatten[{v /. (y1_)[t] :> y1}];
   df = {1, y'[t]} /. Solve[de, y'[t]];
   First[df /. y[t] -> y] /; Length[df] == 1
   ];

Expand@dirfield[p'[t] - (2 p[t]^2) - 3 p[t] + 4 p[t]^-9 == 0, p, t]
(*  {1, -(4/p^9) + 3 p + 2 p^2}  *)

Since the OP's ODE is autonomous, it is natural to consider the equilibria.

equilibria = NSolve[-(4/p^9) + 3 p + 2 p^2 == 0, Reals]
(*  {{p -> -1.45198}, {p -> -1.21521}, {p -> 0.978767}}  *)

There is also a singularity at p == 0. We can incorporate these features into StreamPlot[].

StreamPlot[
 dirfield[p'[t] - (2 p[t]^2) - 3 p[t] + 4 p[t]^-9 == 0, p, t],
 {t, -1, 1}, {p, -2, 1.5},
 StreamPoints -> {Join[{{0, p}, Red} /. equilibria, {Automatic}]} // Evaluate,
 GridLines -> {None, {0}}, GridLinesStyle -> Green,
 PlotRangePadding -> 0
 ]

Mathematica graphics

The behavior near p == 0.978767 can be inferred, but a better picture can be obtained by zooming in. With StreamPlot[], I would advise keeping the intervals of t and p approximately equal in length. If you want a much larger range of t than p or vice versa, I would advise rescaling so that the ranges are about the same.

StreamPlot[
 dirfield[p'[t] - (2 p[t]^2) - 3 p[t] + 4 p[t]^-9 == 0, p, t],
 {t, -0.1, 0.1}, {p, 0.9, 1.1},
 StreamPoints -> {{Last[{{0, p}, Red} /. equilibria], Automatic}},
 PlotRangePadding -> 0
 ]

Mathematica graphics

$\endgroup$
3
  • $\begingroup$ I am unable to understand the negative values of t, is there a discontinuity issue. I also tried to use the manipulate the plot. I used the following code,P[t_, a_, b_, c_, d_] = p[t] /.First@NDSolve[{p'[t] - (a*p[t]^2) - b*p[t] + c*(p[t]^(1 - d)) == 0, p[0] == 1}, p, {t, 0, 0.335}]; Manipulate[Plot[P[t, a, b, c, d], {t, 0, 0.335}, PlotLabel -> P[x, a, b]], {{a, 1, "a"}, -10, 10, Appearance -> "Labeled"}, {{b, 1, "b"}, -5, 5, Appearance -> "Labeled"}, {{c, 1, "c"}, -5, 5, Appearance -> "Labeled"}, {{d, 1, "d"}, -5, 5, Appearance -> "Labeled"}] but I am not getting any graph $\endgroup$ Nov 4 '16 at 6:51
  • $\begingroup$ @SudiptaSen There does seem to be a discontinuity. The ODE is autonomous, so I don't see what's not to understand about pos./neg. t, except in light of an initial condition. What condition corresponds to C[1] -> 1? -- As for your code, doesn't it give errors? $\endgroup$
    – Michael E2
    Nov 4 '16 at 10:37
  • $\begingroup$ @ Michael E2 There maybe a problem with the code, as NDSolve is encountering problem with calculating the derivative at t=0 and I chose the value of C1 arbitarily, anyways I got the solution for this problem. By the way, thanks for showing how to bypass the rootsum issue while plotting. $\endgroup$ Nov 5 '16 at 6:04

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.