# Topic 5 Series Solutions

## 5.1 Power series solutions

Recall that a point \(x_0\) is a regular point of the equation \(y''+p(x)y'+q(x)y=0\), where \(p\) and \(q\) are rational functions, if \(x_0\) is not a discontinuity of both \(p(x)\) and \(q(x)\). At a regular point \(x_0\), we can expect a power series solution \(y(x)=\sum\limits_{n=0}^\infty a_nx^n\). The coefficients are be determined recurrence relations that can be obtained by the method of undetermined coefficients. In maple, we can find a power series solution by the method described above, or using the command `powsolve`

supported by the package `powerseries`

, or using `dsolve`

with the option `series`

. The `powsolve`

returns a procedure which can be expressed explicitly but truncated using the command `tpsform`

.

**Example 5.1 **Find the general solution in power series of the equation
\[y''-xy'-y=0.\]

*Solution*. We first defined the differential equation.

`ode41:=diff(y(x), x, x)-x*diff(y(x), x)-y(x)=0;`

**Method 1:** Follow the manual procedure.

Suppose \(y(x)=\sum\limits_{n=0}^\infty a_nx^n\). In Maple, we define \(y\) using the command `sum`

.

`y:=sum(a[n]*x^n, n=0..infinity);`

Now plug \(y\) in to the equation and simplify.

`pode41 := combine(subs[eval](y(x) = y, ode41));`

The output is \[\color{blue}{\sum\limits_{n=0}^\infty\left(-n a_{n} x^{n -2}+a_{n} x^{n -2} n^{2}-n a_{n} x^{n}-a_{n} x^{n}\right) = 0}.\]

For polynomials, Maple can combine like terms. Unfortunately, for series, Maple couldn’t do the shifting or combine like terms for us. Since we are interested in the coefficient of the general term in \(x^n\), we can extract the terms that contains \(x^n\). The trick it to substitute `n=0..infinity`

by `n=k..k+2`

. Because, only when \(n=k\), \(k+2\), there are terms in \(x^k\).

`kterms:=simplify(subs((n=0..infinity)=(n=k..k+2), pode41));`

Now we can extract the coefficient of \(x^k\) using the command `coeff`

. To keep the equation form, we use `map`

to apply the command `coeff`

to \(x^k\) in both sides of the equation.

`coeffk:=map(coeff, kterms, x^k);`

The output is \[\color{blue}{\left(1+k \right) \left(k a_{k +2}-a_{k}+2 a_{k +2}\right) = 0}.\]

We can then get the recurrence relations by isolate \(a_{k+2}\) and replace \(k\) by \(n-2\).

`recurrence_relations:=subs[eval](k=n-2, isolate(coeffk, a[k+2]));`

The output is \[\color{blue}{a_{n} = \frac{a_{n -2}}{n}}.\]

To get an explicit formula for \(a_n\), we use the Maple command `rsolve`

`coeffn:=rsolve(recurrence_relation, a[n]);`

We can get a polynomial approximation up to \(x^{10}\) for the solution \(y\).

`y:=a[0] + a[1]*x + sum(eval(coeffn, n = k)*x^k, k = 2 .. 10);`

**Method 2:** Using the command `powsolve`

.

Load the package `powseries`

`with(powseries):`

Solve the equation using `powsolve(equation)`

. Note that the output is a procedure that that produces in each step the coefficient of the degree `_k`

term.

`ypowsol:=powsolve(ode41);`

The coefficient of the \(n\)-th term is given by

`an:= eval(ypowsol(_k), _k=n);`

To get a truncated power series, we use the Maple command `tpsform`

. For example, the truncated power series for `ypowsol`

up to 10th order can be obtained by the following command.

`tpsform(ypowsol, x, 10);`

**Method 3:** Using `dsolve(equation, y(x), type=‘series’).

Using this method, we get a power series solution directly.

`dsolve({ode41, y(0) = a[0], D(y)(0) = a[1]}, y(x), 'series');`

Running this command gives the solution \[y\left(x \right) = a_{0}+a_{1} x +\frac{1}{2} a_{0} x^{2}+\frac{1}{3} a_{1} x^{3}+\frac{1}{8} a_{0} x^{4}+\frac{1}{15} a_{1} x^{5}+\mathrm{O}\left(x^{6}\right)\]

Each method has advantages and disadvantages.

- The first method works in general, but not efficient.
- The
`powsolve`

command can only work with polynomial coefficient equations and the power series solution is always at 0. - The
`dsolve`

command cannot produce recurrence relations.

In all three methods, initial conditions can be added.

- In the first method, initial conditions can be imposed with
`rsolve({recurrence_relation, a[0]=a0, a[1]=a1})`

. - In the second and the third method, initial conditions can be imposed by
`{equation, y(0)=a0, y'(0)=a1}`

.

**Exercise 5.1 **Find the general solution in power series form of the equation
\[y''+ xy'+y=0.\]

## 5.2 Euler Equations

For differential equations with regular singular point, the first method and the the third method still work. If the equation is an Euler equation, that is, \(ax^2y''+bxy'+c=0\). The solution is easier to find even by hands. Since once can always use `dsolve`

to solve a differential equation. We will focus on how to solve Euler equations manually. Recall the following theorem.

**Theorem 5.1 **Suppose the solutions of the indicial equation
\[ar(r-1)+br+c=0\]
are \(r_1\) and \(r_2\). Then the general solution of the Euler equation
\[ax^2y''+bxy'+cy=0\]
on \((0,\infty)\) is

\[y= c_1x^{r_1}+c_2x^{r_2}\] if \(r_1\) and \(r_2\) are distinct real numbers;

\[y= x^{r}(c_1+c_2\ln x)\] if \(r_1=r_2=r\);

\[y=x^{\lambda}\left[c_1\cos\left(\omega\ln x\right)+ c_2\sin\left(\omega\ln x \right)\right]\] if \(r_1,r_2=\lambda\pm \mathrm{i}\omega\) with \(\omega>0\).

From the theorem, we see that the key to find the general solution is to solve the indicial equation \(ar(r-1)+br+c=0\). This equation can be obtained substituting \(y''\) by \(r(r-1)\), \(y'\) by \(r\), and \(x\) and \(y\) by \(1\). It can also be obtained by the Maple command `indicialeq(ode, indepvar, regsp, depvar)`

supported by `DETools`

.

**Example 5.2 **Find the general solution of the Euler equations
\[x^2y''-xy'+5y=0\]

*Solution*. Define the equation

`ode42:=x^2*diff(y(x), x$2)-x*diff(y(x), x)+5*y(x)=0;`

The indicial equation can be obtained using one of the following command

```
indeq1:=subs[eval]({diff(y(x), x, x)=r*(r-1), diff(y(x), x)=r, y(x)=1, x=1}, ode42);
indeq2:=DETools[indicialeq](ode42, x, 0, y(x));
```

Solve the indicial equation, we will use the equation `indeq1`

,

`rcomplex:=solve(indeq1, r);`

The solutions are \[\color{blue}{1+2I, 1-2I}\]

Construct the solution using the theorem. First find the real and imaginary parts.

```
alpha:=Re(rcomplex[1]);
beta:=Im(rcomplex[2]);
```

So the general solution is \[y=e^{alpha x}(c_1\cos(\beta\ln x)+c[2]\sin(\beta\ln x)).\]

You can check that using the Maple command `odetest`

.

```
yg := x -> x^alpha*(c[1]*cos(beta*ln(x)) + c[2]*sin(beta*ln(x))):
odetest(y(x)=yg(x), ode42);
```

The output \(\color{blue}{0}\) means that `yg`

is a solution of the equation `ode42`

.

**Exercise 5.2 **Find the general solution of the Euler equation
\[x^2y''-2y'-4y=0.\]

**Exercise 5.3 **Find the general solution of the Euler equation
\[x^2y''-3y'+5y=0.\]

**Exercise 5.4 **Find the general solution of the Euler equation
\[x^2y''+3y'+y=0.\]