Frobenius and power series solutions to ODEs in Wolfram

Not all ordinary differential equations, or ODEs, can be solved analytically, and power series or Frobenius series solutions can then be extremely helpful. The current essay provides Wolfram Language function to help automate power and Frobenius series solutions, a manual work that is done much better by computer than by human. The essay looks at linear homogenous ODEs with polynomial coefficients. A major disadvantage of series solutions is that generally, they only converge within a limited region, which severely restricts their usage; in this essay, a method for piecewise Frobenius/power series is given, that allows to overcome the challenge of limited convergence radius by sewing Frobenius series solution with several power series solutions.
Here is a newer post on this topic: https://community.wolfram.com/groups/-/m/t/2182346. It presents a method for finding piecewise solutions that is more efficient and is easier to use.
In[]:=
SaveDefinitions  True;

Quick demonstration of the results

Solve the equation
2
2
x
y''+7x(x+1)y'-3y=0
near the singular point
x=0
step-by-step using Frobenius’s method:
In[]:=
frobeniusDSolve[2
2
x
y''[x]+7x(x+1)y'[x]-3y[x],y,x,0,"StepByStep"True]
ODE to solve:
-3y[x]+7x(1+x)
′
y
[x]+2
2
x
′′
y
[x] == 0
Try solution of the form y[x]
r
t
∞
∑
n=0
n
t
a
n
, where t≡x:
-3
∞
∑
n=0
n+r
t
a
n
+7
∞
∑
n=0
(n+r)
n+r
t
a
n
+2
∞
∑
n=0
(-1+n+r)(n+r)
n+r
t
a
n
+7
∞
∑
n=0
(n+r)
1+n+r
t
a
n
 0
Reindex:
∞
∑
n=1
7(-1+n+r)
n+r
t
a
-1+n
+
∞
∑
n=0
-3
n+r
t
a
n
+
∞
∑
n=0
7(n+r)
n+r
t
a
n
+
∞
∑
n=0
2(-1+n+r)(n+r)
n+r
t
a
n
 0
Equate the coefficients before all powers of t to zero:
-3
a
0
+5r
a
0
+2
2
r
a
0
0
a
n
-
7(-1+n+r)
a
-1+n
-3+7(n+r)+2(-1+n+r)(n+r)
Possible values of r such that
a
0
≠ 0:
r ∈ -3,
1
2

Out[]=
-3,
1
2
,{},a[n]-
7(-1+n+r)a[-1+n]
-3+7(n+r)+2(-1+n+r)(n+r)

The output tells us that
r
is either
-3
or
1/2
and that
a
n
=-
7(-1+n+r)
a
n-1
-3+7(n+r)+2(-1+n+r)(n+r)
. The function can also solve equations that don’t have a singular point; in that case, the solution reduces to a power series solution:
In[]:=
frobeniusDSolve[f''[x]+f[x],f,x,0]
Out[]=
{0,1},{r(1+r)a[1]0},a[n]-
a[-2+n]
(-1+n+r)(n+r)

r
being either 0 or 1 corresponds to the first term being proportional to either
0
x
or
1
x
, which is the good old power series. Numerically compute the two linearly independent solutions of
2
2
x
y''+7x(x+1)y'-3y=0
up to the first
10
terms:
In[]:=
frobeniusNDSolve[2
2
x
y''[x]+7x(x+1)y'[x]-3y[x],y,x,0,10]
Out[]=

1-
21x
5
+
49
2
x
5
-
343
3
x
15
3
x
,
x
1-
7x
18
+
49
2
x
264
-
1715
3
x
20592
+
16807
4
x
494208
-
117649
5
x
9335040
+
823543
6
x
193489920
-
117649
7
x
89303040
+
823543
8
x
2190901248
-
5764801
9
x
57994444800
+
40353607
10
x
1648263168000

In[]:=
{Plot[%[[1]],{x,0,2},ImageSizeMedium],Plot[%[[2]],{x,0,2},ImageSizeMedium]}
Out[]=

0.5
1.0
1.5
2.0
-20
-15
-10
-5
5
10
,
0.5
1.0
1.5
2.0
0.2
0.4
0.6
0.8

The more terms we compute, the more accurate the solution becomes. Here is what the second solution (
x
1-
7x
18
+
49
2
x
264
+...
) looks like for different ncoefs:
Out[]=
ncoefs = 4
ncoefs = 8
ncoefs = 100
Error, |2
2
x
y'' + 7x(x + 1) y' - 3y|
With ncoefs = 100, it is possible to reduce the error to be of order
-15
10
. Obtain a piecewise Frobenius solution to an equation whose Frobenius series solution has a limited convergence radius:
Out[]=
Piecewise Frobenius
'Regular' Frobenius
Solutions to 2
2
x
​(x + 3
2
)
y''[x] + 7x(x + 3) y'[x] - 3 y[x] = 0
While the ‘regular’ Frobenius solution seems to start diverging at
x=3
(since
x=-3
is a singular point of the equation, from Fuchs’ theorem we would expect the radius of convergence to be around 3), the piecewise solution continues to be accurate. We can assess the accuracy of the piecewise solution by plotting the absolute value of the right-hand-side of the equation,
2
2
x
2
(x+3)
y''+7x(x+3)y'-3y
, which is zero for an exact solution.
Out[]=
|2
2
x
​(x + 3
2
)
y'' + 7x(x + 3) y' - 3 y|
Solutions to 2
2
x
​(x + 3
2
)
y''[x] + 7x(x + 3) y'[x] - 3 y[x] = 0

Step-by-step solutions

Numerically computing the coefficients

Piecewise solutions

The overall principle

Placing the nodes

Gluing two solutions together

Assembling the piecewise solution

Results

Comparison to frobeniusNDSolve

Options

Open questions

Bibliography

http://dslavsk.sites.luc.edu/courses/other/classnotes/frobenius.pdf
​
https://www.math.mcgill.ca/gantumur/math315w14/downloads/frobenius.pdf
​
https://ir.library.oregonstate.edu/downloads/6395w964x