Multistep methods

Multistep methods are an approach to exploit the history of solutions \(y_{\le i}\) in the calculation of the next step \(y_{i+1}\). Note that since they exploit a history, they require bootstrapping in order to get started.

In general we know that

\[\begin{split} \begin{align} \frac{dy}{dt} &= f(t,y) \\ \int dy &= \int f(t,y) dt\\ y_{i+1} &= y_i + \int f(t,y) dt \end{align} \end{split}\]

but remember we have formulae for integration which amounts to weighted sums of the function at different points! We can then expand the right hand side and write a General Linear Multistep formula:

\[ \sum_{j=0}^{s} \alpha_j y_{i+j} = h \sum_{j=0}^{s} \beta_j f(t_{i+j}, y_{i+j}) \]

where the \(\alpha_j\) and \(\beta_j\) coefficients are chosen for accuracy and balance of computation. Note that this indexing is a matter of stype; in general we are interested in the final solution \(y_{i+s}\), and \(\alpha_s=1\) as a matter of normalization and convenience.

Explicit linear multistep methods

When \(\beta_s = 0\), the final solution, \(y_i+s\) depends only on previous solutions and is therefore explicit.

This scheme is often used to generate a predictor based solely on what has happened before. While a good guess, it must be corrected using some kind of scheme like Heun’s method in order to capture anything that happened during the interval.

This scheme does not handle unevenly spaced steps (a significant shortfall!), and doesn’t self-start. Moreover, explicit single-step methods generally outperform these, and they are seldom used.

Backward Difference Formulae: Implicit linear multistep methods

If we say \(\beta_{j\ne s}=0\), we arrive at the implicit linear multistep method, better known as the Backward Difference Formulae which is the default for many modern computational tools.

Starting from: \(\frac{dy}{dx} = f(x,y)\), the \(\alpha_j\) coefficients are found from the derivative of a Lagrange interpolation polynomial fit to the distory: \(<x_n, y_n> ... <x_{n+s}, y_{n+s}>\).

The first 5 orders are:

\[\begin{split} \begin{align} y_{n+1}^{(1)} &= y_n + h f(x_{n+1}, y_{n+1})\\ y_{n+1}^{(2)} &= \frac{4}{3} y_n - \frac{1}{3} y_{n-1} + \frac{2}{3} h f(x_{n+1}, y_{n+1})\\ y_{n+1}^{(3)} &= \frac{18}{11} y_n - \frac{9}{11} y_{n-1} + \frac{2}{11} y_{n-2} + \frac{6}{11} h f(x_{n+1}, y_{n+1})\\ y_{n+1}^{(4)} &= \frac{48}{25} y_n - \frac{36}{25} y_{n-1} + \frac{16}{25} y_{n-2} - \frac{3}{25} y_{n-3} + \frac{12}{25} h f(x_{n+1}, y_{n+1})\\ y_{n+1}^{(5)} &= \frac{300}{137} y_n - \frac{300}{137} y_{n-1} + \frac{120}{137} y_{n-2} - \frac{25}{137} y_{n-3} + \frac{12}{137} y_{n-4} + \frac{60}{137} h f(x_{n+1}, y_{n+1}) \end{align} \end{split}\]

Note that the Backward Euler method is BDF1.

The benefits of BDFs are:

• There is no requirement for a constant step size

• It is implicit and therefore good for stiff equations

• One can dynamically change between orders to self-start and restart if the physics change.