# Backward differentiation formulas

Bill Gear (2007), Scholarpedia, 2(8):3162. | doi:10.4249/scholarpedia.3162 | revision #91024 [link to/cite this article] |

## Contents |

## Backward Differentiation Methods

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

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

\[ y(t) \approx y(t_n) + (t - t_n)\frac{y_n - y_{n-1}}{t_n - t_{n-1}} \]

so the approximation to the derivative is the familiar

\[ y^\prime_n \approx \frac{y_n - y_{n-1}}{t_n - t_{n-1}} \]

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

\[\tag{1} y^\prime = f(y,t) \]

by replacing the derivative on the left hand side of equation
(1), one obtains the Backward Euler method

\[\tag{2} y_n = y_{n-1} + (t_n - t_{n-1})f(y_n,t_n) \]

If \(y_{n-1}\) is known, then equation (2)
is *implicit* in \(y_n\) --- it occurs on both sides of
the equation. (Implicitness is essential for arbitrarily
Stiff Systems.) Because equation (2) is based on a linear
approximation to \(y(t)\ ,\) it is a first-order method.

BDFs can also be used directly to solve some DAEs of the form \[\tag{3} F(y^\prime,y,t) = 0 \]

by replacing the derivative term by an expression such as in equation (2) to get \[\tag{4} F(\frac{y_n - y_{n-1}}{t_n - t_{n-1}},y_n,t_n) = 0 \]

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
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, \(h_n = t_n - t_{n-1}\ ,\) and they
are constant if \(h_n = h\) independent of
\(n\ .\) The \(k\)-th degree polynomial
interpolating \(y_n, y_{n-1}, \cdots, y_{n-k}\) can be
written using backward differences as

\(\tag{5} 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\)

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

\[\tag{6} hy^\prime_n = \sum_{j = 0}^k \alpha_{kj} y_{n-j}\]

where the coefficients \(\{\alpha_{kj}\}\) are

Order k |
\(\alpha_{k0}\ !\)! \(\alpha_{k1}\ !\)
\(\alpha_{k2}\ !\)! \(\alpha_{k3}\ !\) \(\alpha_{k4}\ !\)! \(\alpha_{k5}\ !\) \(\alpha_{k6}\) | ||||||
---|---|---|---|---|---|---|---|

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 |

While equation (6) 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 \(y\) at
\(t_n\) 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 \(h\lambda\) such that when the method is applied to the *test* equation
\(y^\prime = \lambda y\) with complex \(\lambda\ ,\) the numerical
solution is non-increasing in magnitude. The stability regions for the methods of order 1 through 6 are shown inFigure 1 and Figure 2. 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)

**Internal references**

- 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.