Outline
1. Source
2. Linear multistep methods
2.1 Absolute stability
2.2 Stability regions
References
RJL Chapter: 6 & 7
Notes
The central concepts in the analysis of linear multistep
methods, and indeed any numerical method for differential
equations, are
convergence, order /consistency , and stability .
Consider the initial value problem
u ′ = f (u, t), t ∈ [a, b], u(a) = α
and the corresponding implicit 2-step method
1
[3vn+1 − 4vn + vn−1 ] = f (vn+1 , tn+1 ), n = 0, 1, ...
2k
Exercise 1
Show that the above implicit 2-step scheme is a consistent
second order approximation for the ordinary differential
equation u ′ = f (u, t).
Note!
Consistency does not automatically guarantee convergence,
however. . .
Exercise 2
Consider the differential equation: u ′ = f (u, t) where
f = 2t(1 + u 2 ) and t ∈ (0, 1], with the following scheme.
3 1
vn+1 = − vn + 3vn−1 − vn−2 + 3kf (vn , tn ).
2 2
The exact solution of the differential equation is
u(t) = tan(t 2 ).
▶ Show that the scheme consistent of O(k 3 ).
▶ Write a Python code to solve the equation using the above
scheme. Comment on the behavior of the error.
Notes
▶ All the methods discussed above fall under
linear-multistep-methods. Their general structure is
r
X r
X
αj vn+j = k βj f (vn+j , tn+j )
j=0 j=0
▶ One can distinguish between explicit and implicit methods.
If βr = 0, then the method is explicit.
▶ On the other hand, if βr ̸= 0, the method is implicit.
▶ Iterative methods such as Newton’s method or fixed-point
iterations are often used to solve the implicit formula.
▶ In this case, we need to solve
r −1
X
vn+r = kβr f (vn+r , tn+r ) − [αj vn+1 − kβj f (vn+j , tn+j )]
j=0
for vn+r .
Notes
▶ If k is sufficiently small, implicit linear multistep
methods also have unique solutions given v0 , v1 , . . . , vr −1 .
▶ To see this, let L be the Lipschitz constant for f . Given
vn , . . . , vn+r −1 , the value of vn+r is obtained by solving the
equation
vn+r = kβr f (vn+r , tn+r ) + hn ,
where
r −1
X
hn = [kβj f (vn+j , tn+j ) − αj vn+1 ] , which is a constant...
j=0
▶ i.e., we are looking for a fixed point of
g (x) = kβr f (x, tn+r ) + hn
Notes
▶ Notice that if k |βr | L < 1, then g is a contraction:
|g (x) − g (y )| ≤ k |βr | |f (x, tn+r ) − f (y , tn+r )| ≤ k |βr | L |x − y |
▶ So by the Contraction Mapping Fixed-Point-Theorem, g has a
(0)
unique fixed point. Any initial guess vn+r will converge
to the fixed point of
(ℓ+1) (ℓ)
vn+r = kβr f (vn+r , tn+r ) + hn ,
▶ In practice, one chooses either
▶ iterate to convergence, or
▶ a fixed number of iterations.
Stability of linear multistep methods
We introduce the so-called characteristic polynomials to check
the stability of linear multistep methods.
r
X r
X
ρ(ξ) = αj ξ j , and σ(ξ) = βj ξ j
j=0 j=0
i.e., the 1st and 2nd characteristic polynomials resp.
Definition 1 (The Dahlquist root condition)
An r −step linear multistep method is said to be zero-stable
if the roots of the characteristic polynomial ρ(ξ) satisfy the
following conditions
▶ |ξj | ≤ 1 for all j = 1, 2, · · · , r ,
▶ If ξj is a repeated root, then |ξj | < 1.
Remark 1
Methods that do not satisfy the root condition are called
unstable.
Notes
▶ The numerical solution of a one-step method depends on the
initial condition v0 , but the numerical solution of an
r −step method depend on all the r starting values,
v0 , v1 , . . . , vr −1 .
▶ It is thus of interest whether the numerical solution is
stable with respect to perturbations in the starting
values.
▶ A linear multistep method is zero-stable for a certain
differential equation on a given time interval, if a
perturbation in the starting values of size ε causes the
numerical solution over that time interval to change by no
more than K ε for some value of K which does not depend on
the step size k.
▶ This is called zero-stability because it is enough to
check the condition for the differential equation u ′ = 0.
Notes
▶ Consider the general homogeneous linear difference
equation
Xr
αj vn+j = 0
j=0
▶ We assume this equation has a solution of the form vn = ξ n ,
which suggests
Xr
αj ξ j = 0
j=0
▶ That is, ξ n is solution of the homogeneous equation
provided ξ is a solution of the characteristic equation
r
X
ρ(ξ) = αj ξ j
j=0
Example 1
Consider the 2-step Adams-Moulton scheme
k
vn+2 = vn+1 + {−f (vn ) + 8f (vn+1 ) + 5f (vn+1 )}
12
The characteristic polynomials are
1
ρ(ξ) = ξ 2 − ξ, σ(ξ) = (−1 + 8ξ + 5ξ 2 ).
12
Exercise 3
Consider the scheme considered earlier
3 1
vn+1 = − vn + 3vn−1 − vn−2 + 3kf (vn ).
2 2
Use the characteristic polynomial of the scheme to discuss the
stability of the scheme.
Theorem 1 (Dahlquist equivalence theorem)
For linear multistep methods applied to the initial value
problem for u ′ = f (u, t),
consistency + zero-stability ⇐⇒ convergence
Definition 2
A linear-multistep-method is called zero-stable if it
satisfies the Dahlquist root condition.
Remark 2
A consistent one-step linear-multistep-method is always
zero-stable.
Exercise 4
Consider the following 4th order Miline’s method
4
vn+1 = vn−3 + k {2f (vn , tn ) − f (vn−1 , tn−1 ) + 2f (vn−2 , tn−2 )}
3
Use the Dahlquist Equivalence Theorem to investigate the
convergence of the scheme.
Notes
▶ To determine whether a numerical method will produce
reasonable results with a given value of k > 0, we need a
notion of stability that is different from zero-stability.
▶ There are a wide variety of other forms of stability that
have been studied in various contexts.
▶ The one that is most basic and suggests itself from the
above examples is absolute stability.
▶ We can look at the simplest case of the test problem in
which g (t) = 0 and we simply have
u ′ = λu(t)
▶ Euler’s method applied to this problem gives
vn+1 = (1 + kλ)vn
Notes
▶ We say that this method is absolutely stable when
|1 + kλ| ≤ 1
otherwise it is unstable.
▶ Note that there are two parameters k and λ, but only their
product z = kλ matters.
▶ The method is stable whenever z ∈ [−2, 0], and we say that
the interval of absolute stability for Euler’s method is
[−2, 0].
▶ It is more common to speak of the region of absolute
stability as a region in the complex z plane, allowing the
possibility that λ is complex (of course the time step k
should be real and positive).
Notes
▶ For a general linear multistep method the region of
absolute stability is found by applying the method to
u ′ = λu, obtaining
r
X r
X
αj vn+j = k βj λvn+j
j=0 j=0
▶ Which can be rewritten as
r
X
(αj − zβj )vn+j = 0.
j=0
▶ The stability polynomial can be expressed in terms of the
characteristic polynomials for the linear multistep method
as
π(ξ; z) = ρ(ξ) − zρ(ξ).
Definition 3
The region of absolute stability of a linear-multistep method
is the set of points z in the complex plane for which the
polynomial π(ξ; z) satisfies the root condition.
Remark 3
Zero-stability tells us that a method will converge if we take
k small enough, but tells us nothing about what the solution
will look like for a fixed k.
Exercise 5
Consider the trapezoidal method
un+1 − un 1
= [f (un ) + f (un+1 )].
k 2
Show that the region of absolute stability of the scheme is
the left half plane.
Notes
▶ For the forward Euler’s method, π(ξ, z) = ξ − (1 + z),
▶ With a single root ξ = 1 + z,
▶ Now |1 + z| ≤ 1 if −2 ≤ z ≤ 0
Plotting the stability region
import matplotlib . pyplot as plt
import numpy as np
N = 100
a = -np.pi
b = np.pi
x = np. linspace (a, b, N)
y = np. linspace (a, b, N)
X,Y = np. meshgrid (x, y)
Z = X + 1j*Y # defin z
xi = 1 + Z # the stability polynomial
zlevel = np.abs(xi)
h = plt. contourf (X,Y,zlevel ,[0 ,1])
plt.show ()
Exercise 6
To approximate the initial value problem: u ′ = λu, for t > 0,
consider a multistep method
5 1
vn+1 = 2vn−1 − vn + k f (vn , tn ) + f (vn−1 , tn−1 ) .
2 2
What is the order of consistence of the scheme? Find and
sketch the region of absolute stability of the scheme?
Exercise 7
Derive the Adams-Bashforth three-step method by the following
method. Set
u(tn+1 ) = u(tn ) + akf (tn , u(tn )) + bkf (tn−1 , u(tn−1 )) + ckf (tn−2 , u(tn−2 )).
Expand u(tn+1 ), f (tn−1 , u(tn−1 )) and f (tn−2 , u(tn−2 )) in Taylor
series abou (tn , u(tn )) and equate coefficients of k, k 2 and k 3
to obtain a, b and c.