# Ordinary differential equation with variable (polynomial) coefficients -form of solution

If you have a differential equation with polynomial coefficients then quite remarkably Mathematica can find a solution. For example for the second order differential equation

$$y''(t)+b y'(t)+\left(\text{c0}+\text{c1} t+\text{c2} t^2+\text{c3} t^3\right) y(t)$$

where the coefficient of y(t) is a polynomial in t we can get a solution from DSolve

eqn = y''[t] + b y'[t] + (c0 + c1 t + c2 t^2 + c3 t^3) y[t] == 0;

sol = y[t] /. First@DSolve[eqn, y[t], t]


If we put in values we can plot

vals = {c0 -> 10, c1 -> 0.5, c2 -> 0.2, c3 -> 0.1, b -> 0.1,
C[1] -> 0, C[2] -> 1};
Plot[Evaluate[sol /. vals], {t, 0, 10}]


The solution may be used to do algebraic calculation. Here I take the derivative and plot on the phase plane.

dsol = D[sol, t];
ParametricPlot[Evaluate[{sol, dsol} /. vals], {t, 0, 10},
AspectRatio -> 1]


However, I would really like to know what the function looks like. FunctionExpand and DifferentialRootReduce do nothing.

sol // FunctionExpand
DifferentialRootReduce[sol, t]


I know that the function behind DifferentialRoot is probably so complicated and long that it is difficult to understand and do anything useful with. However, the general form may be useful. I am curious and would like to see the function behind the icon. Can this be done?

• You are right, the analytic solution of such an ODE [homogeneous second order differential equation with non-constant coefficients] is too complicated. However, it exists and the algorithm is known due to Kovacic. For higher orders it is known that linear ODEs do not have a general solution Picard–Vessiot theory. Oct 9, 2021 at 21:10
• If in your equation c3=0, the solution can be expressed in terms of linear combination of BesselJ and BesselY of order 1/3. Oct 9, 2021 at 21:16
• Thanks for all the comments. We can get a glimpse of the solution by working out the series around some point, This shows some structure but is not quite what I was hoping for.
– Hugh
Oct 11, 2021 at 8:39
• Your may also find relevant Nikiforov-Uvarov formalism arxiv.org/pdf/quant-ph/0604021.pdf .
– Acus
Nov 9, 2021 at 6:25
– Hugh
Nov 9, 2021 at 16:46

The manual/help for DifferentialRoot says " FunctionExpand will attempt to expand DifferentialRoot objects in terms of ordinary special and elementary functions."

Test it for various combinations with Manipulate. E.g. for all c4 == 0 or for (c1 and c2 and c3 and b) == 0 you get explicit functions.

eqn = y''[t] + b y'[t] + (c0 + c1 t + c2 t^2 + c3 t^3) y[t] == 0;

ysol[d1_, d2_, b_, c0_, c1_, c2_, c3_] =
y /. First@DSolve[eqn && y[0] == d1 && y'[0] == d2, y, t]

Manipulate[{Plot[
Evaluate[ysol[d1, d2, b, c1, c2, c3, c4][t]], {t, 0, 7},
PlotRange -> All],
ysol[d1, d2, b, c1, c2, c3, c4][t] // FunctionExpand[#, Assumptions -> t > 0] &}, {{d1, 1}, 0,
4, 1/100, Appearance -> "Labeled"}, {{d2, 2}, 0, 3, 1/100,
Appearance -> "Labeled"}, {{b, 1}, 0, 3, 1/100,
Appearance -> "Labeled"}, {{c1, 1}, 0, 4, 1/100,
Appearance -> "Labeled"}, {{c2, 1}, 0, 4, 1/100,
Appearance -> "Labeled"}, {{c3, 1}, 0, 4, 1/100,
Appearance -> "Labeled"}, {{c4, 1}, 0, 4, 1/100,
Appearance -> "Labeled"}, ContinuousAction -> False]

• Thanks. I had looked at this. It is not clear how the simpler forms progress to the one I want to see. Any ideas?
– Hugh
Oct 11, 2021 at 8:31