# Series and Recurrences

This tutorial contains many *Mathematica* scripts. You, as the user, are free to use all codes for your needs, and have the
right to distribute this tutorial and refer to this tutorial as long as
this tutorial is accredited appropriately. Any comments and/or contributions for this tutorial are welcome; you can send your remarks to <Vladimir_Dobrushkin@brown.edu>

Return to computing page for the second course APMA0330

Return to Mathematica tutorial page for the first course

Return to Mathematica tutorial page for the second course

Return to the main page for the course APMA0330

Return to the main page for the course APMA0340

Part V. Series and Recurrences

Recurrences

Sturm-Liouville problem

Heat transfer equations

Wave equations

Laplace equation

Applications

# 5.1. Recurrences

## I. Case Sensitivity

To solve the linear constant coefficient recurrence (also called the difference equation) of the second order

y_{n+2} =6y_{n+1} -y_n subject to y_0 = 0, y_1 = 1

you may try some options:

Here the first pair represents the coefficients of the difference equation of order 2; next pair tells *Mathematica* the first two initial values of the recurrence. The final number (10) is the number of terms the user wants to see as output. Another option:

�� NestList[Flatten[{Rest[#], ListCorrelate[{-1, 6}, #]}] &, {0, 1},

���� 10]][[1]]

If you know that its solution is represented through the Chebyshev polynomial of the second kind at the point x=3, then type:

If the generating function Sum[y[n]*z^n , {n,0,infinity}] is known, then (to see the first 10 coefficients)

Consider another recurrence that generates the Pell numbers:

y_{n+2} = 2 y_{n+1} + y_n , y_0 =0, y_1 =1, n=0,1,2,.... .

Then these numbers can be generated as follows:

�� a = b; b = c, {n, 20}];

lst

Recurrences, although a very tedious computation method by hand, is very simple to do in *Mathematica*. The best way to learn how to do recurrences in *Mathematica* are by examples, and a perfect example for this topic is the Fibonacci sequence.

# 5.2. Generating Functions

Recurrences, although a very tedious computation method by hand, is very simple to do in *Mathematica*. The best way to learn how to do recurrences in *Mathematica* are by examples, and a perfect example for this topic is the Fibonacci sequence.

The process of shifting (or slipping) indices can be obtained by the following *Mathematica* code:

**Example**. Consider the function

# 5.3. Series Solutions for the first Order Equations

## I. How to

Consider the initial value problem for the Abel equation

y' = y^3 +x, y(0) =0

We find its power series expansion using the following *Mathematica* commands:

c[5] x^27 + c[6] x^32 + c[7] x^37 + c[8] x^42 + c[9] x^47 +O[x]^48}

coefeq = LogicalExpand[de];

coeffs = Solve[coefeq, Table[c[i], {0, 1, 17}]];

y = y /. coeffs

x^22 / 262721536 + 879 x^27 / 2101772288 + }}

This problem can be solved in slightly different way. First, we define the differential operator (non-linear because the given equation is non-linear)

Then we set n to be 13, the number of terms in the power series:

The next cell says to create a sum of n terms and effectively turn it into a series saying that the terms beyond n are indefinite. O[x]^(n+1) indicates that we know nothing about terms of order n+1 and beyond. So we define the (finite) series

a[6] x^6 + a[7] x^7 + a[8] x^8 + a[9] x^9 + a[10] x^10 + a[11] x^11 +

a[12] x^12 + a[13] x^13 + O[x]^14

Substitute this expression into the differential operator

Equate the coefficients to 0. Recall that the operator && means

**And**.

Find the coefficients a[i] in the poer series of y[x] using
**Reverse[Table[a[i],{i,0,n}]]**. The latter makes the form of the

coefficients agree with what you would find by hand: the next term is
expressed through the previous terms.

To get the solution series, set **a[0]** equals to zero:

# 3. Series Solutions

Example: solve the IVP:

a[n_]:=a[n]=((n-2) a[n-2] -a[n-1])/4/n

a[0] = a0

a[1] = -a0/4

TableForm[Table[{n,a[n]},{n,0,5}]]

0 a0

1 -a0/4

2 a0/32

3 -3 a0/128

4 11 a0/2048

5 -31 a0/8192

When we apply the initial condition, we set a0=1 in the general solution.

Now we solve the initial value problem exactly:

exactsol=DSolve[{(4-x^2)y'[x]+y[x]==0,y[0]==1},y[x],x]

formula=Simplify[(2-x)^(1/4)/(2+x)^(1/4)]

Out[64]= (2 - x)^(1/4)/(2 + x)^(1/4)

var=Table[D[y[x],i]/.x->0,{i,1,11}]

????sols=Solve[sysofeqs,vars]

sersol=Series[y[x],{x,0,11}/.sols[[1]]]

Legendre's equation

genform =

Solve[(n + 2)*(n + 1)*a[n + 2] + (-n*(n - 1) - 2 n + k*(k + 1))*

a[n] == 0, a[n + 2]]

Out[8]= {{a[2 + n] -> -(((k + k^2 - n - n^2) a[n])/((1 + n) (2 + n)))}}

Factor[genform[[1, 1, 2]]]

Out[9]= ((-k + n) (1 + k + n) a[n])/((1 + n) (2 + n))

We obtain a formula for a[n] by replacing each occurrence of n in a[n+2] by n-2.

genform[[1, 1, 2]] /. n -> n - 2

Out[10]= -(((2 + k + k^2 - (-2 + n)^2 - n) a[-2 + n])/((-1 + n) n))

a[n_] := -(((2 + k + k^2 - (-2 + n)^2 - n) a[-2 + n])/((-1 + n) n))

DSolve[(1 - x^2)*y''[x] - 2 x y'[x] + k*(k + 1)*y[x] == 0, y[x], x]

Out[12]= {{y[x] -> C[1] LegendreP[k, x] + C[2] LegendreQ[k, x]}}