Linear multistep method
2090865
222149219
2008-06-27T20:10:59Z
Giftlite
37986
/* Analysis */ -sp
'''Linear multistep methods''' are methods used in mathematics for the [[numerical ordinary differential equations|numerical solution of ordinary differential equations]]. One-step methods (such as [[Euler's method]] and [[Runge–Kutta methods]]) refer only to one previous value to determine the current value. Multistep methods refer to several previous function values in an effort to achieve greater accuracy. In the case of ''linear'' multistep methods, a [[linear combination]] of the previous function values is used.
==Introductory example==
Linear multistep methods solve initial value problems of the form
: <math> y' = f(t,y), \quad y(t_0) = y_0. </math>
Consider for example the problem
: <math> y' = y, \quad y(0) = 1. </math>
The exact solution is <math> y(t) = \mathrm{e}^t </math>.
A simple numerical method is Euler's method:
: <math> y_{n+1} = y_n + hf(t_n, y_n). \, </math>
This method, applied with step size <math> h = \tfrac12 </math> on the problem <math> y' = y </math>, gives the following results:
: <math> \begin{align}
y_1 &= y_0 + hf(y_0) = 1 + \tfrac12\cdot1 = 1.5, \\
y_2 &= y_1 + hf(y_1) = 1.5 + \tfrac12\cdot1.5 = 2.25, \\
y_3 &= y_2 + hf(y_2) = 2.25 + \tfrac12\cdot2.25 = 3.375, \\
y_4 &= y_3 + hf(y_3) = 3.375 + \tfrac12\cdot3.375 = 5.0625.
\end{align} </math>
Euler's method is a one-step method. A simple multistep method is the two-step Adams–Bashforth method
: <math> y_{n+2} = y_{n+1} + \tfrac32 hf(t_{n+1},y_{n+1}) - \tfrac12 hf(t_n,y_n). </math>
This method needs two values, <math> y_{n+1} </math> and <math> y_n </math>, to compute the next value, <math> y_{n+2} </math>. However, the initial value problem provides only one value, <math> y_0 = 1 </math>. One possibility to resolve this issue is to use the <math> y_1 </math> computed by Euler's method as the second value. With this choice, the Adams–Bashforth method yields (rounded to four digits):
: <math> \begin{align}
y_2 &= y_1 + \tfrac32 hf(y_1) - \tfrac12 hf(y_0) = 1.5 + \tfrac32\cdot\tfrac12\cdot1.5 - \tfrac12\cdot\tfrac12\cdot1 = 2.375, \\
y_3 &= y_2 + \tfrac32 hf(y_2) - \tfrac12 hf(y_1) = 2.375 + \tfrac32\cdot\tfrac12\cdot2.375 - \tfrac12\cdot\tfrac12\cdot1.5 = 3.7812, \\
y_4 &= y_3 + \tfrac32 hf(y_3) - \tfrac12 hf(y_2) = 3.7812 + \tfrac32\cdot\tfrac12\cdot3.7812 - \tfrac12\cdot\tfrac12\cdot2.375 = 6.0234.
\end{align} </math>
The exact solution at <math> t = t_4 = 2 </math> is <math> \mathrm{e}^2 = 7.3891\ldots </math>, so the two-step Adams–Bashforth method is more accurate than Euler's method. This is always the case if the step size is small enough.
==Definition==
A linear multistep method is a method of the form
: <math> \begin{align}
& y_{n+s} + a_{s-1} y_{n+s-1} + a_{s-2} y_{n+s-2} + \cdots + a_0 y_n \\
& \qquad {} = h \bigl( b_s f(t_{n+s},y_{n+s}) + b_{s-1} f(t_{n+s-1},y_{n+s-1}) + \cdots + b_0 f(t_n,y_n) \bigr),
\end{align} </math>
where ''h'' denotes the step size and ''f'' the right-hand side of the differential equation. The coefficents <math> a_0, \ldots, a_{s-1} </math> and <math> b_0, \ldots, b_s </math> determine the method.
The method is explicit if <math> b_s = 0 </math>. In that case, the recurrence relation can be used directly to compute <math> y_{n+s} </math>. If <math> b_s \ne 0 </math> then the method is implicit and the recurrence relation is an equation for <math> y_{n+s} </math> which has to be solved.
==Examples==
Three families of linear multistep methods are commonly used: Adams–Bashforth methods, Adams–Moulton methods, and the [[backward differentiation formula]]s (BDFs).
=== Adams–Bashforth methods ===
The Adams–Bashforth methods are explicit methods. The coefficients are <math>a_{s-1} = -1</math> and <math>a_{s-2} = \cdots = a_0 = 0</math>, while the <math>b_j</math> are chosen such that the methods has order ''s'' (this determines the methods uniquely).
The Adams–Bashforth methods with ''s'' = 1, 2, 3, 4 are {{harv|Hairer|Nørsett|Wanner|1993|loc=§III.1}}:
* <math> y_{n+1} = y_n + hf(t_n, y_n) \, </math> — this is simply the Euler method;
* <math> y_{n+2} = y_{n+1} + h\left( \tfrac32 f(t_{n+1}, y_{n+1}) - \tfrac12 f(t_n, y_n)\right); </math>
* <math> y_{n+3} = y_{n+2} + h\left( \tfrac{23}{12} f(t_{n+2}, y_{n+2}) - \tfrac43 f(t_{n+1}, y_{n+1}) + \tfrac{5}{12}f(t_n, y_n)\right); </math>
* <math> y_{n+4} = y_{n+3} + h\left( \tfrac{55}{24} f(t_{n+3}, y_{n+3}) - \tfrac{59}{24} f(t_{n+2}, y_{n+2}) + \tfrac{37}{24} f(t_{n+1}, y_{n+1}) - \tfrac38 f(t_n, y_n) \right). </math>
The coefficients <math>b_j</math> can be determined as follows. Use [[polynomial interpolation]] to find the polynomial ''p'' of degree <math>s-1</math> such that
:<math> p(t_{n+i}) = f(t_{n+i}, y_{n+i}), \qquad \text{for } i=0,\ldots,s-1. </math>
The [[Lagrange polynomial|Lagrange formula]] for polynomial interpolation yields
:<math> p(t) = \sum_{j=0}^{s-1} \frac{(-1)^{s-j-1}f(t_{n+j}, y_{n+j})}{j!(s-j-1)!h^{s-1}} \prod_{i=0 \atop i\ne j}^{s-1} (t-t_{n+i}). </math>
The polynomial ''p'' is locally a good approximation of the right-hand side of the differential equation <math>y' = f(t,y)</math> that is to be solved, so consider the equation <math>y' = p(t)</math> instead. This equation can be solved exactly; the solution is simply the integral of ''p''. This suggests taking
:<math> y_{n+s} = y_{n+s-1} + \int_{t_{n+s}}^{t_{n+s-1}} p(t)\,dt. </math>
The Adams–Bashforth method arises when the formula for ''p'' is substituted. The coefficients <math>b_j</math> turn out to be given by
:<math> b_{s-j-1} = \frac{(-1)^j}{j!(s-j-1)!} \int_0^1 \prod_{i=0 \atop i\ne j}^{s-1} (u+i) \,du, \qquad \text{for } j=0,\ldots,s-1. </math>
Replacing ''f''(''t'', ''y'') by its interpolant ''p'' incurs an error of order ''h''<sup>''s''</sup>, and it follows that the ''s''-step Adams–Bashforth method has indeed order ''s'' {{harv|Iserles|1996|loc=§2.1}}
The Adams–Bashforth methods were designed by [[John Couch Adams]] to solve a differential equation modelling [[capillary action]] due to [[Francis Bashforth]]. {{harvtxt|Bashforth|1883}} published his theory and Adams' numerical method {{harv|Goldstine|1977}}.
=== Adams–Moulton methods ===
The Adams–Moulton methods are similar to the Adams–Bashforth methods in that they also have <math>a_{s-1} = -1</math> and <math>a_{s-2} = \cdots = a_0 = 0</math>. Again the ''b'' coefficients are chosen to obtain the highest order possible. However, the Adams–Moulton methods are implicit methods. By removing the restriction that <math>b_s = 0</math>, an ''s''-step Adams–Moulton method can reach order <math>s+1</math>, while an ''s''-step Adams–Bashforth methods has only order ''s''.
The Adams–Moulton methods with ''s'' = 0, 1, 2, 3 are {{harv|Hairer|Nørsett|Wanner|1993|loc=§III.1}}:
* <math> y_n = y_{n-1} + h f(t_n,y_n) \, </math> — this is the [[backward Euler method]];
* <math> y_{n+1} = y_n + \tfrac12 h \big( f(t_{n+1},y_{n+1}) + f(t_n,y_n) \big) </math> — this is the [[trapezoidal rule (numerical differential equations)|trapezoidal rule]];
* <math> y_{n+2} = y_{n+1} + h \big( \tfrac5{12} f(t_{n+2},y_{n+2}) + \tfrac23 f(t_{n+1},y_{n+1}) - \tfrac1{12} f(t_n,y_n) \big); </math>
* <math> y_{n+3} = y_{n+2} + h \big( \tfrac38 f(t_{n+3},y_{n+3}) + \tfrac{19}{24} f(t_{n+2},y_{n+2}) - \tfrac5{24} f(t_{n+1},y_{n+1}) + \tfrac1{24} f(t_n,y_n) \big). </math>
The derivation of the Adams–Moulton methods is similar to that of the Adams–Bashforth method; however, the interpolating polynomial uses not only the points ''t''<sub>''n''−1</sub>, … ''t''<sub>''n''−''s''</sub>, as above, but also <math>t_n</math>. The coefficients are given by
:<math> b_{s-j} = \frac{(-1)^j}{j!(s-j)!} \int_0^1 \prod_{i=0 \atop i\ne j}^{s} (u+i-1) \,du, \qquad \text{for } j=0,\ldots,s. </math>
The Adams–Moulton methods are solely due to [[John Couch Adams]], like the Adams–Bashforth methods. The name of [[Forest Ray Moulton]] became associated with these methods because he realized that they could be used in tandem with the Adams–Bashforth methods as a [[predictor–corrector]] pair {{harv|Moulton|1926}}; {{harvtxt|Milne|1926}} had the same idea. Adams used [[Newton's method]] to solve the implicit equation {{harv|Hairer|Nørsett|Wanner|1993|loc=§III.1}}.
== Analysis ==
The central concepts in the analysis of linear multistep methods, and indeed any numerical method for differential equations, are [[Numerical ordinary differential equations#Analysis|convergence, order, and stability]].
The first question is whether the method is consistent: is the difference equation
:<math> \begin{align}
& y_{n+s} + a_{s-1} y_{n+s-1} + a_{s-2} y_{n+s-2} + \cdots + a_0 y_n \\
& \qquad {} = h \bigl( b_s f(t_{n+s},y_{n+s}) + b_{s-1} f(t_{n+s-1},y_{n+s-1}) + \cdots + b_0 f(t_n,y_n) \bigr),
\end{align} </math>
a good approximation of the differential equation <math> y' = f(t,y) </math>? More precisely, a multistep method is ''consistent'' if the local error goes to zero as the step size ''h'' goes to zero, where the ''local error'' is defined to be the difference between the result <math>y_{n+s}</math> of the method, assuming that all the previous values <math>y_{n+s-1}, \ldots, y_n</math> are exact, and the exact solution of the equation at time <math>t_{n+s}</math>. A computation using [[Taylor series]] shows out that a linear multistep method is consistent if and only if
:<math> \sum_{k=0}^{s-1} a_k = -1 \quad\text{and}\quad \sum_{k=0}^s b_k = s + \sum_{k=0}^{s-1} ka_k. </math>
All the methods mentioned above are consistent {{harv|Hairer|Nørsett|Wanner|1993|loc=§III.2}}.
If the method is consistent, then the next question is how well the difference equation defining the numerical method approximates the differential equation. A multistep method is said to have ''order'' ''p'' if the local error is of order <math>O(h^{p+1})</math> as ''h'' goes to zero. This is equivalent to the following condition on the coefficients of the methods:
:<math> \sum_{k=0}^{s-1} a_k = -1 \quad\text{and}\quad q k^{q-1} \sum_{k=0}^s b_k = s^q + \sum_{k=0}^{s-1} k^q a_k \text{ for } q=1,\ldots,p. </math>
The ''s''-step Adams–Bashforth method has order ''s'', while the ''s''-step Adams–Moulton method has order <math>s+1</math> {{harv|Hairer|Nørsett|Wanner|1993|loc=§III.2}}.
These conditions are often formulated using the ''characteristic polynomials''
:<math> \rho(z) = z^s + \sum_{k=0}^{s-1} a_k z^k \quad\text{and}\quad \sigma(z) = \sum_{k=0}^s b_k z^k. </math>
In terms of these polynomials, the above condition for the method to have order ''p'' becomes
:<math> \rho(\mathrm{e}^h) - h\sigma(\mathrm{e}^h) = O(h^{p+1}) \quad \text{as } h\to 0. </math>
In particular, the method is consistent if it has order one, which is the case if <math>\rho(1)=0</math> and <math>\rho'(1)=\sigma(1)</math>.
If the roots of the characteristic polynomial <math>\rho</math> all have modulus less than or equal to 1 and the roots of modulus 1 are of multiplicity 1, we say that the [[root condition]] is satisfied.
The method is convergent [[if and only if]] it is consistent and the root condition is satisfied. Consequently, a consistent method is stable if and only if this condition is satisfied, and thus the method is convergent if and only if it is stable.
Furthermore, if the method is stable, the method is said to be ''strongly stable'' if <math>z=1</math> is the only root of modulus 1, otherwise, it is said to be ''weakly stable''.
=== Example ===
Consider the Adams–Bashforth three-step method
:<math>y_{n+1} = y_n + h\left( {23\over 12} f(t_n, y_n) - {16 \over 12} f(t_{n-1}, y_{n-1}) + {5\over 12}f(t_{n-2}, y_{n-2})\right).</math>
The characteristic equation is thus
:<math>q(z)=z^3-z^2=z^2(z-1)\,</math>
which has roots <math>z=0, 1</math>, and the conditions above are satisfied. As <math>z=1</math> is the only root of modulus 1, the method is strongly stable.
==First and second Dahlquist barriers==
These two results were proved by [[Germund Dahlquist]] and represent an important bound for the order of convergence and for the A-stability of a linear multistep method.
'''First Dahlquist barrier'''
A [[zero stability|zero-stable]] and linear ''q''-step multistep methods cannot attain an order of convergence greater than ''q'' + 1 if ''q'' is odd and greater than ''q'' + 2 if ''q'' is even. If the method is also explicit, then it cannot attain an order greater ''q'' {{harv|Hairer|Nørsett|Wanner|1993|loc=Thm III.3.5}}.
'''Second Dahlquist barrier'''
There are no explicit [[A-stability|A-stable]] and linear multistep methods. The implicit ones have order of convergence at most 2 {{harv|Hairer|Wanner|1996|loc=Thm V.1.4}}.
== References ==
* {{citation | first1 = Francis | last1 = Bashforth | year = 1883 | title = An Attempt to test the Theories of Capillary Action by comparing the theoretical and measured forms of drops of fluid. With an explanation of the method of integration employed in constructing the tables which give the theoretical forms of such drops, by J. C. Adams | location = Cambridge }}.
* {{Citation | last1=Butcher | first1=John C. | author1-link = John C. Butcher | title=Numerical Methods for Ordinary Differential Equations | publisher=John Wiley | isbn=978-0-471-96758-3 | year=2003}}.
* {{citation | first1 = Herman H. | last1 = Goldstine | author1-link = Herman Goldstine | year = 1977 | title = A History of Numerical Analysis from the 16th through the 19th Century | publisher = Springer-Verlag | location = New York | isbn = 978-0-387-90277-7 }}.
* {{citation | first1 = Ernst | last1 = Hairer | first2 = Syvert Paul | last2 = Nørsett | first3 = Gerhard | last3 = Wanner | year = 1993 | title = Solving ordinary differential equations I: Nonstiff problems | edition = 2nd | publisher = Springer Verlag | location = Berlin | isbn = 978-3-540-56670-0 }}.
* {{Citation | last1=Hairer | first1=Ernst | last2=Wanner | first2=Gerhard | title=Solving ordinary differential equations II: Stiff and differential-algebraic problems | publisher=[[Springer-Verlag]] | location=Berlin, New York | edition=2nd | isbn=978-3-540-60452-5 | year=1996}}.
* {{citation | first1 = Arieh | last1 = Iserles | year = 1996 | title = A First Course in the Numerical Analysis of Differential Equations | publisher = Cambridge University Press | isbn = 978-0-521-55655-2 }}.
* {{citation | first1 = W. E. | last1 = Milne | year = 1926 | title = Numerical integration of ordinary differential equations | journal = American Mathematical Monthly | volume = 33 | number = 9 | pages = 455–460 | url = http://links.jstor.org/sici?sici=0002-9890%28192611%2933%3A9%3C455%3ANIOODE%3E2.0.CO%3B2-6 }}.
* {{citation | first1 = Forest R. | last1 = Moulton | author1-link = Forest Ray Moulton | year = 1926 | title = New methods in exterior ballistics | publisher = University of Chicago Press }}.
== External links ==
* {{MathWorld | urlname= AdamsMethod| title= Adams Method }}
* [http://math.fullerton.edu/mathews/n2003/AdamsBashforthMod.html Adams-Bashforth-Moulton Method]
* [http://engineersexcel.com/Tools/2nd%20Order%20LODE%20solver/Description.htm Spreadsheet ODE Solver that uses Adams-Bashforth-Moulton Method]
[[Category:Numerical differential equations]]
[[de:Mehrschrittverfahren]]
[[ru:Метод Адамса]]