7
$\begingroup$

I am trying to solve the following ODE, $$ \frac{1}{2}\sigma^2 x^2 g'' + \frac{\nu x^2 g'^2}{\eta^2} + \frac{(\alpha\eta+2\nu) x g'}{\eta} - \beta g + \nu = 0\ , $$ where the Greek letters are constants. I have so far only managed to locate the steady-state solution $g(x) = \nu /\beta$. Due to a broken laptop, I'm currently forced to use Wolfram Cloud, where the computation times out. So my question: is Mathematica in fact able to solve this ODE? If not, I would be grateful for any other solution ideas/hints, e.g. in terms of implicit functions.

The code is the following:

pdeg = \[Beta] g[x] == \[Nu] + (x (\[Alpha] \[Eta]+2 \[Nu]) g'[x])/\[Eta] + (x^2 \[Nu] g'[x]^2)/\[Eta]^2 + 1/2 x^2 \[Sigma]^2 g''[x]

DSolve[pdeg, g, x]
$\endgroup$
10
  • $\begingroup$ It does not seem to be solvable., Even with numerical values for the parameters. Mathematica can't solve and neither can Maple. So it looks like lost case. Try numerical solution or series solution. It is because of the $(g')^2$ term in there. $\endgroup$
    – Nasser
    Commented Jan 12 at 14:11
  • $\begingroup$ @Nasser Thank you for running the code! Would you be able to explain how to construct a series solution in Mathematica, or provide a reference? Thanks! $\endgroup$
    – Anthony
    Commented Jan 12 at 14:14
  • 1
    $\begingroup$ @Anthony What is known about the x-range? 0<x<Infinity? $\endgroup$ Commented Jan 12 at 15:16
  • $\begingroup$ @UlrichNeumann Yes, that is correct, $0 < x < \infty$. For my application, I am interested in $x$ rather large, say $\mathcal{O}(10^4)$. $\endgroup$
    – Anthony
    Commented Jan 12 at 15:24
  • 1
    $\begingroup$ It might be of interest to note that Mathematica can solve the equation analytically when $\beta = 0$. In this case, the solution generally has a logarithmic singularity as $x \to 0$, and diverges proportionally to $\ln x$ as $x \to \infty$. $\endgroup$ Commented Jan 12 at 20:34

2 Answers 2

5
$\begingroup$

To illustrate my comment about $x=e^t$ turning the ODE into an autonomous one. Autonomous just means the coefficients are constants or equivalently the vector field $dY/dt = F(Y)$ is time-independent. See https://mathworld.wolfram.com/Autonomous.html, for instance. It reduces the dimension of the system by 1. In a second-order ODE, the reduction permits the visualization of phase portrait (2D) or the contact manifold (surface in 3D).

Phase portrait example (note g[x] = h[Log[x]]):

hode = DSolveChangeVariables[Inactive[DSolve][pdeg, g, x], h, t, 
     x == Exp[t]] // First // Simplify;
With[{hA = -3, hB = 5, hpA = -4, hpB = 2,
  hpp = First@SolveValues[hode, h''[t]]},
 Manipulate[
  StreamPlot[{h'[t], hpp}, {h[t], hA, hB}, {h'[t], hpA, hpB},
   Epilog -> {RGBColor[0, 0.8, 0.6], PointSize@Large, 
     Point[{\[Nu]/\[Beta], 0}]}],
  {{\[Alpha], 1}, 0, 2}, {{\[Beta], 1}, 0, 2}, {{\[Eta], 1}, 0, 2},
  {{\[Nu], 1}, 0 + $MachineEpsilon, 2}, {{\[Sigma], 1}, 0, 2},
  TrackedSymbols :> {\[Beta], \[Nu], \[Alpha], \[Sigma], \[Eta]}
  ]
 ]

Note on symmetry: In an autonomous ODE, if $h_1$ is a solution, then so is $h_1(t+C)$ for any constant $C$. Correspondingly, under the substitution $x=e^t$, if $g_1$ is a solution to the OP's ODE, then so is $g_1(C\,x)$ (for $C \ne 0$).

$\endgroup$
2
  • $\begingroup$ Clever idea x->Exp[t]. But unfortunately the stationary point doesn't look stable. Perhaps that might be the reason why series solution can't be derived? $\endgroup$ Commented Jan 13 at 21:39
  • $\begingroup$ @UlrichNeumann Thanks, but why "unfortunately"? Shouldn't we say fortunately we are able to determine the stability of the stationary point? In fact, the determinant of the Jacobian is a constant $-{2 \beta }/{\sigma^2}$. Of course, with different parameter settings, we can get stationary points of different types — but not for $0 < \beta < 1$, which is the range that the OP is interested in. It's a saddle if $\beta$ is positive. $\endgroup$
    – Goofy
    Commented Jan 13 at 22:22
9
$\begingroup$

how to construct a series solution,

It can solve it by series, but expansion point can not be zero. (the reason is, at $x=0$ the whole ode goes away, and you end up with the equation

 β  g[0] == ν

Here is your ode

pdeg=β g[x]==ν+(x (α η+2 ν) g'[x])/η+(x^2 ν g'[x]^2)/η^2+1/2 x^2 σ^2 g''[x]

Mathematica graphics

AsymptoticDSolveValue[pdeg, g[x], {x, 0, 4}]
(*no solution*)

but using 4 terms, the solution around $x=1$ is

AsymptoticDSolveValue[pdeg, g[x], {x, 1, 4}]

Mathematica graphics

$\endgroup$

Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.