Richardson Extrapolation

If we know the form of a method's error, we can cancel its leading term by combining two results. For a central difference, $$D(h) = f'(x) + c_2h^{2} + c_4h^{4} + \cdots$$ The constants $c_2, c_4$ are unknown, but they are the same for every $h$, and that is enough.

Key idea

If $A(h) = A + ch^{p} + O(h^{q})$ with $q > p$, then $$A \approx \frac{2^{p}A(h/2) - A(h)}{2^{p} - 1},$$ with error $O(h^{q})$. For an $h^{2}$ error: $\frac{4A(h/2) - A(h)}{3}$. For an $h$ error: $2A(h/2) - A(h)$. As a bonus, $\frac{A(h/2) - A(h)}{2^{p} - 1}$ estimates the error of $A(h/2)$.

Why it works

Write the two results: $A(h) = A + ch^{p} + \cdots$ and $A(h/2) = A + \frac{c}{2^{p}}h^{p} + \cdots$. Multiply the second by $2^{p}$ and subtract the first: the $ch^{p}$ terms cancel, leaving $(2^{p} - 1)A$ plus higher-order terms. Divide by $2^{p} - 1$. For the error estimate, subtract instead: $A(h/2) - A(h) \approx -\frac{2^{p} - 1}{2^{p}}ch^{p} = (2^{p} - 1)\big(A(h/2) - A\big)$.

Blue: central differences for eˣ at 0 (slope 2). Purple: Richardson-extrapolated values (slope 4), many digits better for the same h.

Worked example: Extrapolating a central difference

For $f(x) = e^{x}$ at $x = 0$ (exact derivative $1$): $D(0.2) = 1.0066800$ and $D(0.1) = 1.0016675$. Extrapolate.

The error is $O(h^{2})$, so use $\frac{4D(0.1) - D(0.2)}{3} = \frac{4.0066700 - 1.0066800}{3}$.
$= 0.9999967$, an error of about $3.3 \times 10^{-6}$.
$D(0.1)$ alone was off by $1.7 \times 10^{-3}$: one subtraction bought about three more digits.

Answer: $0.9999967$

Worked example: Estimating the error for free

From the same two numbers, estimate the error of $D(0.1)$.

$\frac{D(0.1) - D(0.2)}{3} = \frac{1.0016675 - 1.0066800}{3} \approx -0.0016708$.
So the true value is about $1.0016675 - 0.0016708 \approx 0.99999$; indeed the actual error of $D(0.1)$ is $0.0016675$.

Answer: about $-0.00167$

Worked example: An order-1 method

Forward differences for the derivative of $e^{x}$ at $x = 1$ give $A(0.1) = 2.8588$ and $A(0.05) = 2.7874$. Their error is $O(h)$. Extrapolate.

For $p = 1$: $2A(h/2) - A(h) = 2(2.7874) - 2.8588 = 2.7160$.
The exact value is $e \approx 2.7183$: the error fell from $0.069$ to $0.002$.
Using $\frac{4A(h/2) - A(h)}{3}$ here would cancel a term that is not there and leave the $h$ term only partly removed.

Answer: $2.7160$

Extrapolation is cheap: it needs no new function values, only arithmetic on results already computed. It is also a safety check. If the extrapolated value and $A(h/2)$ differ by much more than the estimated error, the assumed form of the error is wrong, perhaps because $h$ is not yet small enough.

The idea goes further: with three results one can cancel two error terms, and repeating it on the trapezoid rule gives Romberg integration. The one requirement is that the error really behaves like $ch^{p}$; for a function that is not smooth enough, or with $h$ so small that round-off matters, extrapolation does not help.

Common mistake

Using the weights for the wrong order. $\frac{4A(h/2) - A(h)}{3}$ cancels an $h^{2}$ term; $2A(h/2) - A(h)$ cancels an $h$ term; $\frac{16A(h/2) - A(h)}{15}$ cancels $h^{4}$. Check the order of your method first (or measure it from three results), then use $2^{p}$.

Richardson Extrapolation

Practice

Builds on

Leads to

Sources and further reading