differential equation of order2

I want to solve the following differential equation

ode = x*y''[x] + y'[x] + 4*x^2*(x^2 - 1) y[x] == 0


But do not know how to actually solve it. Any suggestion?

• Change the last y to y[x] and try feeding it again into DSolve[]. May 2, 2015 at 14:01
• I obtain answer to form "DifferentialRoot", and I can not obtain y[x]. May 2, 2015 at 14:03
• and I add conditions: DSolve[ode && y[1] == 0 && y'[1] == 1, y[x], x] May 2, 2015 at 14:11
• DifferentialRoot actually gives you the symbolic solution of the equation. Inside Mma you may operatre with it as with any regular function, i.e., plot, integrate etc. Have a look in Help. May 2, 2015 at 16:46

At first, I was going to answer as bbgodfrey did, except I did one little test that stumped me. It also occurs to me that this ODE might be solvable in terms of hypergeometric functions, but that Mathematica either does not find the solution or does not know how to find it. (This is an area of active research, but it is not my field. I don't know if anyone on this site will know the answer, but the OP might ask on math.stackexchange.com or mathoverflow.com)

Normal DifferentialRoot usage

Normally DifferentialRoot works just like a function, including evaluating it and its derivatives. For instance:

f = DifferentialRoot[
Function[{y, x}, {y'''[x] + 8 y[x] == 0, y[0] == 0,
y'[0] == Sin[2 Pi/3], y''[0] == Sin[4 Pi/3]}]]


Plot[{f[x], f'[x], f''[x]}, {x, 0, 5}, PlotRange -> 20]


The OP's case: A bug?

But with the OP's differential equation we seem to run into a bug:

ode = x*y''[x] + y'[x] + 4*x^2*(x^2 - 1) y[x] == 0;
{ysol} = DSolve[{ode, y[1] == 0, y'[1] == 1}, y, x]


The function y[x] evaluates fine:

Table[y[x], {x, 1., 4.}] /. ysol
(*  {0., 0.0405724, -0.122343, -0.0857817}  *)


However, the derivative at the initial point results in nonsense and remains unevaluated at other values of x:

y'[1.] /. ysol
(*  ReplaceAll[1.]  *)

y'[2.] /. ysol


It seems likely that the problem is that the differential equation in the DifferentialRoot for y'[x] is singular at x == 1:

 y'[x] /. ysol


Nonetheless, ReplaceAll[1.] hardly seems an acceptable return value.

Other approaches

One might:

• Use a differential initial point, such as x == 2, avoiding x == 0 and x == 1. Things seem to work better (and faster).

• One can develop (a finite portion of) a Frobenius series at x == 0, if local analysis is desired.

• It's possible to use the solution returned by DSolve to find another independent solution and to construct the general solution (albeit in terms of the troublesome DifferentialRoot).

Example: Initial conditions at x == 2:

ode = x*y''[x] + y'[x] + 4*x^2*(x^2 - 1) y[x] == 0;
{ysol} = DSolve[{ode, y[2] == C[1], y'[2] == C[2]}, y, x]


Note than the constant C[1] was replace by an expression in terms of C[2] and a new C[3] for some reason. Now we can find initial conditions at x == 2 that give y[1] == 0 and y'[x] == 1. The DifferentialRoot in this case is well behaved, even at x == 1. Note also that DifferentialRoot requires exact coefficients; hence the SetPrecision[_, Infinity].

obj[a_?NumericQ, b_?NumericQ] := {y[1.], y'[1.]} - {0, 1} /. ysol /.
{C[3] -> SetPrecision[a, Infinity], C[2] -> SetPrecision[b, Infinity]};
ic0 = FindRoot[obj[a, b], {{a, 1.}, {b, 1.}}] /. {a -> C[3], b -> C[2]}
(*
{C[3] -> -0.286103, C[2] -> -1.37527}
*)

Plot[{y[x], y'[x]} /. ysol /. SetPrecision[ic0, Infinity] // Evaluate,
{x, 0.2, 5}]


• Nice analysis. Solving the ODE with conditions y[1] == 1, y'[1] == 0 gives the same strange y[x] results you obtained for y[1] == 0, y'[1] == 1. However, the case I plotted in my answer does not give these strange results. May 4, 2015 at 2:29
• @bbgodfrey Thanks. Yeah, I tried it both ways (ics at x == 1), but it's a long answer already. Your method avoids trouble by starting close to the singular point, instead on directly on it (which I guess you realize). May 4, 2015 at 3:11
• You correctly point out that y'[x] /. ysol yields a DifferentialRoot expression that is singular at x == 1. Yet, x == 1 clearly is an ordinary point of ode, and it is straightforward to construct a power series for y[x] there (e.g., using Series[ode[[1]], {x, 1, 3}]). It would appear, therefore, that this instance of DifferentialRoot is incorrect, reinforcing your observation that a bug exists. May 4, 2015 at 13:20
• Dear Michael, can you help me for this question: mathematica.stackexchange.com/questions/86302/… Jun 19, 2015 at 7:47

The point made by Alexei Boulbitch can be illustrated by computing and plotting a solution well behaved near x = 0.

sol = y /. First@DSolve[ode && y[1/1000] == 1 && y'[1/1000] == 0, y, x];
Plot[sol[x], {x, 0, 6}]


The solution decreases both in amplitude and wavelength roughly as 1/x. Like any such second order linear ODE, it admits a second solution, singular at x = 0.