Adams-Moulton methods

Implicit multistep methods that include the new node in the interpolation: full derivation of AM2 (implicit trapezoid), AM4, why they require solving a nonlinear equation and what they gain in exchange.

Implicit: using the new node

Adams-Moulton interpolates ff including the new point (tk+1,yk+1)(t_{k+1},\,y_{k+1}), which we do not yet know. With the same nodes, the interpolation is one degree richer and the method gains one order over the corresponding Adams-Bashforth; in exchange, yk+1y_{k+1} appears inside ff and an equation must be solved at every step.

Derivation of AM2

tₖtₖ₊₁fₖfₖ₊₁interpolationnew slopedepends on yₖ₊₁integrate herep₁(t)
AM2 interpolates ff at tkt_k and at the new node tk+1t_{k+1}. The area is trapezoidal, but fk+1=f(tk+1,yk+1)f_{k+1}=f(t_{k+1},y_{k+1}) depends on the unknown.
Expand diagram

AM2 interpolates ff at tkt_k and at the new node tk+1t_{k+1}. The area is trapezoidal, but fk+1=f(tk+1,yk+1)f_{k+1}=f(t_{k+1},y_{k+1}) depends on the unknown.

DerivationDerivation: one-step Adams-Moulton (AM2)View as its own page →

Step 1: integral form

tₖtₖ₊₁fₖfₖ₊₁interpolationnew slopedepends on yₖ₊₁integrate herep₁(t)
AM2 uses the whole new interval: it interpolates the slope between fkf_k and fk+1f_{k+1}, so it is the implicit trapezoidal rule.
Expand diagram

AM2 uses the whole new interval: it interpolates the slope between fkf_k and fk+1f_{k+1}, so it is the implicit trapezoidal rule.

  1. Start from the IVP y(t)=f(t,y(t))y'(t)=f(t,y(t)). Integrate both sides between tkt_k and tk+1t_{k+1}:

    tktk+1y(τ)dτ=tktk+1f(τ,y(τ))dτ\int_{t_k}^{t_{k+1}} y'(\tau)\,d\tau=\int_{t_k}^{t_{k+1}} f\bigl(\tau,y(\tau)\bigr)\,d\tau
  2. By the Fundamental Theorem of Calculus, the left-hand side is exact. The only part to approximate is the integral of ff:

    y(tk+1)=y(tk)+tktk+1f(τ,y(τ))dτy(t_{k+1})=y(t_k)+\int_{t_k}^{t_{k+1}} f\bigl(\tau,y(\tau)\bigr)\,d\tau
  3. In the numerical formula write yky(tk)y_k\approx y(t_k) and yk+1y(tk+1)y_{k+1}\approx y(t_{k+1}):

    yk+1=yk+tktk+1f(τ,y(τ))dτy_{k+1}=y_k+\int_{t_k}^{t_{k+1}} f\bigl(\tau,y(\tau)\bigr)\,d\tau

Step 2: interpolate f on the new interval

  1. Unlike AB2, AM2 interpolates at tkt_k and tk+1t_{k+1}. Thus it uses fk=f(tk,yk)f_k=f(t_k,y_k) and fk+1=f(tk+1,yk+1)f_{k+1}=f(t_{k+1},y_{k+1}):

    p1(τ)=fkLk(τ)+fk+1Lk+1(τ)p_1(\tau)=f_k L_k(\tau)+f_{k+1}L_{k+1}(\tau)
  2. The Lagrange basis functions are:

    Lk(τ)=τtk+1tktk+1=tk+1τh,Lk+1(τ)=τtktk+1tk=τtkhL_k(\tau)=\frac{\tau-t_{k+1}}{t_k-t_{k+1}}=\frac{t_{k+1}-\tau}{h},\qquad L_{k+1}(\tau)=\frac{\tau-t_k}{t_{k+1}-t_k}=\frac{\tau-t_k}{h}
  3. Substitution gives the line that approximates the slope inside the subinterval:

    p1(τ)=fktk+1τh+fk+1τtkh=fk+1τtkhfkτtk+1hp_1(\tau)=f_k\,\frac{t_{k+1}-\tau}{h}+f_{k+1}\,\frac{\tau-t_k}{h}=f_{k+1}\frac{\tau-t_k}{h}-f_k\frac{\tau-t_{k+1}}{h}

Step 3: change of variable and integrals

  1. To integrate, use the same change as in your notes: s=τtks=\tau-t_k. Then τ=tk+s\tau=t_k+s, dτ=dsd\tau=ds, and the limits are s=0s=0 and s=hs=h.

    s=τtk,τ=tk+s,τtk+1=shs=\tau-t_k,\qquad \tau=t_k+s,\qquad \tau-t_{k+1}=s-h
  2. With this change, the interpolant becomes:

    p1(tk+s)=fk+1shfkshh=fk+1sh+fkhshp_1(t_k+s)=f_{k+1}\frac{s}{h}-f_k\frac{s-h}{h}=f_{k+1}\frac{s}{h}+f_k\frac{h-s}{h}
  3. Integrate the fk+1f_{k+1} term:

    0hfk+1shds=fk+11h[s22]0h=h2fk+1\int_0^h f_{k+1}\frac{s}{h}\,ds=f_{k+1}\frac{1}{h}\left[\frac{s^2}{2}\right]_0^h=\frac{h}{2}f_{k+1}
  4. Integrate the fkf_k term:

    0hfkhshds=fk1h[hss22]0h=h2fk\int_0^h f_k\frac{h-s}{h}\,ds=f_k\frac{1}{h}\left[hs-\frac{s^2}{2}\right]_0^h=\frac{h}{2}f_k
  5. Therefore the slope integral is approximated by the sum of those two weights:

    tktk+1f(τ,y(τ))dτtktk+1p1(τ)dτ=h2(fk+1+fk)\int_{t_k}^{t_{k+1}} f\bigl(\tau,y(\tau)\bigr)\,d\tau\approx \int_{t_k}^{t_{k+1}}p_1(\tau)\,d\tau=\frac{h}{2}\bigl(f_{k+1}+f_k\bigr)
  6. Substituting into the integral form gives AM2:

    yk+1=yk+h2(fk+1+fk)y_{k+1}=y_k+\frac{h}{2}\bigl(f_{k+1}+f_k\bigr)

Step 4: why it is implicit

  1. The key detail is that fk+1f_{k+1} has not been computed yet:

    fk+1=f(tk+1,yk+1)f_{k+1}=f(t_{k+1},y_{k+1})
  2. So the formula is really an equation for the unknown yk+1y_{k+1}:

    yk+1=yk+h2(f(tk+1,yk+1)+f(tk,yk))y_{k+1}=y_k+\frac{h}{2}\Bigl(f(t_{k+1},y_{k+1})+f(t_k,y_k)\Bigr)
  3. If solving with Newton, define a residual R(z)R(z), keeping ff only for the ODE function:

    R(z)=zykh2(f(tk+1,z)+fk)R(z)=z-y_k-\frac{h}{2}\Bigl(f(t_{k+1},z)+f_k\Bigr)
  4. Newton updates the approximation until R(z)R(z) is essentially zero:

    z(m+1)=z(m)R(z(m))R(z(m)),R(z)=1h2fy(tk+1,z)z^{(m+1)}=z^{(m)}-\frac{R(z^{(m)})}{R'(z^{(m)})},\qquad R'(z)=1-\frac{h}{2}\,f_y(t_{k+1},z)
  5. In a predictor-corrector pair, one usually predicts yk+1y_{k+1} with AB2 and inserts that prediction into AM2:

    yk+1(p)=yk+h2(3fkfk1),yk+1(c)=yk+h2(f(tk+1,yk+1(p))+fk)y_{k+1}^{(p)}=y_k+\frac{h}{2}(3f_k-f_{k-1}),\qquad y_{k+1}^{(c)}=y_k+\frac{h}{2}\Bigl(f(t_{k+1},y_{k+1}^{(p)})+f_k\Bigr)

Step 5: error and link with the trapezoid

The formula obtained is exactly the trapezoidal rule applied to the slope integral y(t)y'(t). Since the trapezoid has local error proportional to the third derivative,

ek+1=112h3y(ξ)=O(h3)e_{k+1}=-\frac{1}{12}h^3y'''(\xi)=\mathcal{O}(h^3)
Local error of AM2; after accumulating steps, the global error is order 2.

Thus AM2 has the same global order as AB2, but a smaller error constant. The improvement is paid for by solving the implicit equation or using a corrector.

yk+1=yk+h24(9fk+1+19fk5fk1+fk2)y_{k+1}=y_k+\frac{h}{24}\bigl(9f_{k+1}+19f_k-5f_{k-1}+f_{k-2}\bigr)
AM4, order 4: cubic interpolation including the new node.

Solving the implicit equation

Since yk+1y_{k+1} appears on both sides, each step requires solving a nonlinear equation, often with the Newton-Raphson method:

R(yk+1)=yk+1ykh2(f(tk+1,yk+1)+f(tk,yk))=0R(y_{k+1})=y_{k+1}-y_k-\frac{h}{2}\bigl(f(t_{k+1},y_{k+1})+f(t_k,y_k)\bigr)=0