View source for Backward differentiation formulas

==Backward Differentiation Methods==

These are numerical integration methods based on '''Backward Differentiation Formulas''' (BDFs).  They are particularly useful for
[[Stiff Systems|stiff differential equations]] and
[[Differential-Algebraic Equations]] (DAEs).  BDFs are
formulas that give an approximation to a derivative of a variable at
a time <math>t_n</math> in terms of its function values <math>y(t)
</math>  at <math>t_n</math>  and earlier times (hence the "backward"
in the name).  They are derived by forming the <math>k</math>-th degree
interpolating polynomial approximating the function <math>y(t)
</math> using <math>y(t_n), y(t_{n-1}), \cdots, y(t_{n-k})\ ,</math>
differentiating it, and evaluating it at <math>t_n\ .</math>

For example, the linear interpolating polynomial through <math>y(t_n)
</math>  and <math>y(t_{n-1})</math>   is

:<math>
y(t) \approx y(t_n) + (t - t_n)\frac{y_n - y_{n-1}}{t_n - t_{n-1}}
</math>

so the approximation to the derivative is the familiar

:<math>
y^\prime_n \approx \frac{y_n - y_{n-1}}{t_n - t_{n-1}}
</math>

If this is used to obtain a numerical approximation to the
ordinary differential equation

:<math>\label{ode}
y^\prime = f(y,t) </math>


by replacing the derivative on the left hand side of equation
\eqref{ode}, one obtains the Backward Euler method

:<math>\label{bem}
y_n  = y_{n-1} + (t_n - t_{n-1})f(y_n,t_n) </math>


If <math>y_{n-1}</math> is known, then equation \eqref{bem}
is ''implicit'' in <math>y_n</math> --- it occurs on both sides of
the equation.  (Implicitness is essential for arbitrarily 
[[Stiff Systems]].) Because equation \eqref{bem} is based on a linear
approximation to <math>y(t)\ ,</math> it is a first-order method.

BDFs can also be used directly to solve some DAEs of the form
:<math>\label{dae}
F(y^\prime,y,t) = 0
</math>

by replacing the derivative term by an expression such as in equation
\eqref{bem} to get
:<math>\label{daesol}
F(\frac{y_n - y_{n-1}}{t_n - t_{n-1}},y_n,t_n) = 0
</math>

The reader should be aware that this method will work only if the DAE has
''index'' no greater than 1 or has other special properties (see
[[Differential-Algebraic Equations|DAEs]]).

When higher-order methods are discussed, it is easier to discuss
''constant step size'' methods.  The step size is the gap between
adjacent time points,  <math>h_n = t_n - t_{n-1}\ ,</math> and they
are constant if <math>h_n = h</math> independent of
<math>n\ .</math>  The <math>k</math>-th degree polynomial
interpolating <math>y_n, y_{n-1}, \cdots, y_{n-k}</math> can be
written using backward differences as

<math>\label{interp}
y(t) \approx y_n + \frac{1}{h}(t-t_n)\nabla y_n +
\frac{1}{2 h^2}(t-t_n)(t-t_{n-1})\nabla^2y_n + \cdots
+\frac{1}{h^k k!}(t-t_n)\cdots(t-t_{n-k+1})\nabla^k</math>


where <math>\nabla</math> is the backward difference operator
<math>\nabla y_n = y_n - y_{n-1}</math> and <math>\nabla^p y_n =
\nabla^{p-1}y_n - \nabla^{p-1}y_{n-1}\ .</math>  If equation
\eqref{interp} is differentiated and <math>t</math> is set to
<math>t_n\ ,</math> the <math>k</math>-th order BDF is obtained.  The
formula has the form

:<math>\label{bdf}
hy^\prime_n = \sum_{j = 0}^k \alpha_{kj}
y_{n-j}</math>


where the coefficients <math>\{\alpha_{kj}\}</math> are

{| border="1"
|+BDF Coefficients
!Order ''k''!! <math>\alpha_{k0}\ !</math>! <math>\alpha_{k1}\ !</math>
 <math>\alpha_{k2}\ !</math>! <math>\alpha_{k3}\ !</math>
 <math>\alpha_{k4}\ !</math>! <math>\alpha_{k5}\ !</math>
 <math>\alpha_{k6}</math>
|-
!1
|1||-1|| || || || ||
|-
!2
|3/2 ||-2 || 1/2 || || || ||
|-
!3
|11/6 ||-3 || 3/2 || -1/3 || || ||
|-
!4
| 25/12 ||-4 ||3 || -4/3 || 1/4 || ||
|-
!5
| 137/60 || -5 || 5 ||-10/3 || 5/4 || -1/5 ||
|-
!6
| 49/20 || -6 || 15/2 || -20/3 || 15/4 || -6/5 || 1/6
|}

[[Image:Stabs1-3.jpg|thumb|400px|right|F1|Stability regions for BDF orders 1  - 3]]

[[Image:Stabs4-6.jpg|thumb|400px|right|F2|Stability regions for BDF orders 4 - 6]]

While equation \eqref{bdf} provides an easy way to discuss BDF's, quality codes
implement a variable step size (and variable order) version of these methods, often using
''modified divided differences,'' which are an unequal step size version of the backwards
differences used above.  Another representation often used (with either equal or unequal
step sizes) is the ''Nordsieck array,'' in which the derivatives of <math>y</math> at
<math>t_n</math> are approximated by difference quotients.

The importance of the BDF-based methods is their stability.  The ''region of absolute'' stability
of a method is that set of <math>h\lambda</math> such that when the method is applied to the ''test'' equation
<math>y^\prime = \lambda y</math> with complex <math>\lambda\ ,</math> the numerical
solution is non-increasing in magnitude.  The stability regions for the methods of order 1 through 6 are shown in<figref>Stabs1-3.jpg</figref> and  <figref>Stabs4-6.jpg</figref>.  The stable regions are the ''exterior'' of the contours indicated
(the full contour is plotted only for orders 1 and 2).  The other contours close in the right-half plane.)
Note that these six methods are stable along the whole of the negative real axis, the fact that makes them
suitable for stiff equations.  The contour of the BDF method of order 7 crosses the negative real axis,
making it and any higher order BDF methods of no value.

The first use of BDF methods appears to date back to Curtiss
and Hirschfelder (1952), although they were not
given that name at the time.  Later, Henrici (1962), in Section 5.1-4 of his book,
discussed "methods based on differentiation."  The methods were dismissed by Henrici
because they are "less accurate that the corresponding Adams-Moulton formula."  For non-stiff
equations this is a valid point.  For stiff systems, the value of BDF methods lies in their superior stability properties which allow them to take much larger stepsizes than would be possible with explicit methods.


===References===
* C. F. Curtiss and J. O. Hirschfelder, Integration of Stiff Equations, Proc US Nat. Acad. Sci, 38, #3 pp235-243, March, 1952
*P. Henrici, Discrete variable methods in ordinary differential equations, New York, Wiley (1962)

<!-- Authors, please check this list and remove any references that are irrelevant. This list is generated automatically to reflect the links from your article to other accepted articles in Scholarpedia. -->
<b>Internal references</b>

* Lawrence F. Shampine and Skip Thompson (2007) [[Initial value problems]]. Scholarpedia, 2(3):2861.
* Philip Holmes and Eric T. Shea-Brown (2006) [[Stability]]. Scholarpedia, 1(10):1838.
* Lawrence F. Shampine and Skip Thompson (2007) [[Stiff systems]]. Scholarpedia, 2(3):2855.