HIGHER ORDER METHODS There are two principal means to derive higher order methods = aj + ...

28 Pages·2001·0.08 MB·English

HIGHER ORDER METHODS
There are two principal means to derive higher order
methods
(cid:1)p (cid:1)p
y = a y + h b f(x , y )
− − −
n+1 j n j j n j n j
−
j=0 j= 1
(a) Method of Undetermined Coefficients
(b) Numerical Integration
UNDETERMINED COEFFICIENTS
Pick p > 0 and m > 1. Solve the linear system
(cid:1)p
a = 1
j
j=0
(cid:1)p (cid:1)p
−
− i − i 1
( 1) a + i ( 1) b = 1
j j
−
j=0 j= 1
for i = 1, ..., m.
EXAMPLE
For p = 1, the general two-step formula is
y = a y + a y
n+1 0 n 1 n−1
(1)
+h [b y + b y + b y ]
−1 n+1 0 n 1 n−1
For order m = 2, we must have the coefficients satisfy
the system
a + a = 1
0 1
−
a + b + b + b = 1
−
1 1 0 1
−
a + 2 [b b ] = 1
−
1 1 1
These three equations have two degrees of freedom;
and their general solution is
−
a = 1 a
1 0
− 1 − 1
b− = 1 a b (2)
1 0 0
4 2
− 3 − 1
b = 1 a b
1 0 0
4 2
SOLVING THE LINEAR SYSTEM
The linear system has the matrix
 
1 1 0 0 0 1
 
−
 0 1 1 1 1 1 
−
0 1 2 0 2 1
By elementary row operations, this can be reduced to
 
1 1 0 0 0 1
 
 1 1 
0 1 0 1
 
4 2
3 1
0 0 1 1
4 2
From this, we can solve for a , b , b .
−
1 1 1
For the truncation error in the resulting method, recall
the formula
(cid:8) (cid:9)
c
m+1 m+1 (m+1) m+2
T (Y ) = h Y (x ) + O h
n n
(m + 1)!
 
(cid:1)p (cid:1)p
−  − m+1 − m 
c = 1 ( j) a + (m + 1) ( j) b
m+1 j j
−
j=0 j= 1
For our case with m = 2, we have
c
3 3 (cid:2)(cid:2)(cid:2) 4
T (Y ) = h y (x ) + O(h ) (3)
n n
3!
−
c = 4 + 2a + 3b (4)
3 0 0
How might we choose a and b ?
0 0
§
To apply the earlier convergence theory of 6.3, we
would want to require
≤ ≤
0 a , a 1
0 1
Subject to this, we could choose a and b to
0 0
(cid:10) (cid:10)
(cid:10)− (cid:10)
4 + 2a + 3b
(cid:10) 0 0(cid:10)
Minimize
(cid:10) (cid:10)
6
In fact, for any a in [0, 1], we can choose b so as to
0 0
make this quantity zero, namely let
−
4 2a
0
b =
0
3
Then
4
T (Y ) = O(h )
n
and we could actually calculate a formula
c
4 4 (4) 5
T (Y ) = h y (x ) + O(h )
n n
4!
In fact, we can then go ahead and choose a = 0,
0
4
b = , leading to the method
0
3
(cid:11) (cid:12)
h
(cid:2) (cid:2) (cid:2)
y = y + y + 4y + y (5)
n+1 n−1 n+1 n n−1
3
This method has a truncation error of
5
h
− (5)
T (Y ) = y (ξ )
n n
90
It arises from applying Simpson’s numerical integra-
tion rule to
(cid:13)
x
n+1
Y (x ) = Y (x ) + f(x, Y (x)) dx
−
n+1 n 1
xn−1
See pages 256-258 for Simpson’s quadrature rule.
The above formula (5) is called Milne’s method. It
appears to be quite a good formula. It is convergent
4
with an error of O(h ); and it is stable in the original
§ §
sense of 6.2 and 6.3. Unfortunately, it is only weakly
stable, and this appears when x is sufficiently large.
n
EXAMPLE - CONTINUED
Recall (1)
y = a y + a y
n+1 0 n 1 n−1
+h [b y + b y + b y ]
−1 n+1 0 n 1 n−1
with the choice of coefficients from (2),
−
a = 1 a
1 0
− 1 − 1
b = 1 a b
−
1 0 0
4 2
− 3 − 1
b = 1 a b
1 0 0
4 2
If we want an explicit method, then we must force
b = 0. Doing so yields the choice of coefficients
−
1
1 1
− − −
a = 1 a , b = 2 a , b = a
1 0 0 0 1 0
2 2
Again, our earlier convergence and stability theory
≤ ≤
works if we choose 0 a 1. The truncation error
0
from (3)-(4) becomes
(cid:14) (cid:15)
a 1
0 3 (cid:2)(cid:2)(cid:2) 4
T (Y ) = + h Y (x ) + O(h )
n n
12 3
≤ ≤
For 0 a 1, the leading coefficient is minimized
0
by choosing a = 0. This yields the numerical method
0
y = y + 2hf(x , y )
n+1 n−1 n n
which is simply the midpoint method. Thus choosing
a so as to just minimize the truncation error is not
0
an adequate basis for choosing a method.
In contrast, choosing a = 1 leads to the method
0
h
−
y = y + [3f(x , y ) f(x , y )]
n+1 n n n n−1 n−1
2
with a truncation error
5
(cid:2)(cid:2)(cid:2)
3 4
T (Y ) = h Y (x ) + O(h )
n n
12
This is the Adams-Bashforth method of order 2, and
it does not have the weak stability property.
NUMERICAL INTEGRATION
The most popular multistep formulas in use today are
based on numerical integration, by means of suitably
chosen polynomial interpolation formulas. Integrate
the equation
(cid:2)
Y (x) = f (x, Y (x))
over some interval [x − , x ]. This yields
n r n+1
(cid:13)
x
n+1
Y (x ) = Y (x − ) + f(x, Y (x)) dx
n+1 n r
xn−r
(cid:2)
Approximate Y (x) = f (x, Y (x)) with a polynomial
interpolant, and then integrate this polynomial to ap-
proximate the integral. The most satisfactory formu-
las have been found by doing the case
(cid:13)
x
n+1
Y (x ) = Y (x ) + f(x, Y (x)) dx
n+1 n
x
n
and these are called the Adams family of methods.
POLYNOMIAL INTERPOLATION
Let t , ..., t be p + 1 distinct node points, and let
0 p
z , ..., z be corresponding function values. The unique
0 p
≤
polynomial of degree p which interpolates these val-
ues at the given node points is given by
(cid:1)p
P
(t) = z (cid:24) (t) (6)
p j j
j=0
(cid:17) (cid:18)
(cid:16)p −
t t
i
(cid:24) (t) =
j
−
t t
j i
i=0
(cid:5)
i=j
With this definition, deg((cid:24) ) = p and
j
(cid:19)
1, i = j
≡
(cid:24) (t ) δ =
j i i,j (cid:5)
0, i = j
Formula (6) is called the Lagrange form of the inter-
polation polynomial.