Burden & Faires, 9th Ed. - §4.1–4.2
Chapter 8: Numerical Differentiation & Richardson Extrapolation
Finite-difference formulas · The truncation/round-off tradeoff · Optimal step size · Richardson extrapolation
Numerical differentiation looks like the simplest possible idea - approximate $f'(x_0)$ by the slope between two nearby points - and it is also the clearest place in this entire course to watch Chapter 1's two error sources fight each other directly. Shrinking the step $h$ always shrinks truncation error, but past a certain point it makes round-off error explode, because differentiation subtracts two nearly-equal function values and divides by an even smaller number. This chapter derives the formulas, characterizes exactly where that crossover happens, and then shows a trick - Richardson extrapolation - that gets high accuracy without needing $h$ to be dangerously small.
Section 4.1
Finite-Difference Formulas
Two-Point Formulas
Forward, Backward & Central Difference
Taylor-expanding $f(x_0+h)$ about $x_0$ (Chapter 1, Theorem 1.1.4) and solving for $f'(x_0)$ gives the forward-difference formula:
$$f'(x_0) = \frac{f(x_0+h)-f(x_0)}{h} - \frac{h}{2}f''(\xi), \qquad \xi\in(x_0,x_0+h)$$
so the forward-difference approximation has truncation error $O(h)$. Its mirror image (expand $f(x_0-h)$) gives the backward-difference formula, also $O(h)$. Subtracting the two expansions instead of using just one cancels the $O(h)$ term entirely, giving the central-difference formula:
$$f'(x_0) = \frac{f(x_0+h)-f(x_0-h)}{2h} - \frac{h^2}{6}f'''(\xi), \qquad \xi\in(x_0-h,x_0+h)$$
Central difference is $O(h^2)$ - for small $h$, halving $h$ quarters the truncation error, twice as good an exchange rate as forward or backward difference.
📘 Example 8.1 - Three Formulas, Same Problem
Approximate $f'(1)$ for $f(x)=\sin x$ (true value $\cos 1 = 0.5403023$) using $h=0.1$.
Solution
Central difference is roughly $47\times$ more accurate than either one-sided formula at the identical step size $h=0.1$ and identical cost (two function evaluations) - there is essentially never a reason to prefer a one-sided formula except at a boundary node where the other side's data doesn't exist.
| Formula | Approximation | Error |
|---|---|---|
| Forward: $[f(1.1)-f(1)]/0.1$ | 0.4973638 | $4.29\times10^{-2}$ |
| Backward: $[f(1)-f(0.9)]/0.1$ | 0.5814408 | $4.11\times10^{-2}$ |
| Central: $[f(1.1)-f(0.9)]/0.2$ | 0.5394023 | $9.00\times10^{-4}$ |
Three- and Five-Point Formulas
Higher-Order Formulas from More Points
Differentiating the Lagrange interpolating polynomial itself (Chapter 5) through more nodes produces higher-order finite-difference formulas. The five-point central-difference formula:
$$f'(x_0) = \frac{f(x_0-2h)-8f(x_0-h)+8f(x_0+h)-f(x_0+2h)}{12h} + \frac{h^4}{30}f^{(5)}(\xi)$$
is $O(h^4)$ - dramatically more accurate for the same $h$, at the cost of two extra function evaluations per point. (Three-point endpoint and midpoint formulas, $O(h^2)$, fill the gap between the two-point and five-point formulas and are the standard choice at domain boundaries where five points aren't available.)
📘 Example 8.2 - Five-Point Formula on the Same Problem
Repeat Example 8.1's estimate of $\cos(1)$ using the five-point central formula with $h=0.1$.
Solution
$$f'(1) \approx \frac{\sin(0.8)-8\sin(0.9)+8\sin(1.1)-\sin(1.2)}{1.2} = 0.5403005$$
Error $\approx1.80\times10^{-6}$ - about $500\times$ better than central difference's $9.00\times10^{-4}$ at the exact same $h$, confirming the jump from $O(h^2)$ to $O(h^4)$.
Second-Derivative Formula
Central-Difference Formula for $f''(x_0)$
$$f''(x_0) = \frac{f(x_0-h)-2f(x_0)+f(x_0+h)}{h^2} - \frac{h^2}{12}f^{(4)}(\xi)$$
$O(h^2)$ accurate, and the natural building block for the finite-difference solution of boundary-value ODEs and PDEs (beyond this course's scope, but the same three-point stencil reappears constantly in applied numerical work).
Section 4.2
The Round-off/Truncation Tradeoff
Why "Smaller $h$" Is Not Always "Better"
Every formula above bounds truncation error, which shrinks as $h\to0$. But every formula also divides by $h$ (or $h^2$), and the function values themselves carry round-off error of size $\sim\varepsilon_{\text{mach}}|f(x_0)|$ (Chapter 1). For central difference, the total error behaves like
$$E(h) \;\approx\; \underbrace{\frac{h^2}{6}|f'''(\xi)|}_{\text{truncation, shrinks with }h} \;+\; \underbrace{\frac{\varepsilon_{\text{mach}}|f(x_0)|}{h}}_{\text{round-off, grows as }h\to0}$$
These pull in opposite directions - there is an optimal $h^*$ that minimizes the sum, and pushing $h$ below it makes the answer worse, not better.
Finding the Optimal Step Size
Minimizing $E(h) = \dfrac{M h^2}{6} + \dfrac{\varepsilon_{\text{mach}}|f(x_0)|}{h}$ over $h$ (differentiate w.r.t. $h$, set to zero):
$$\frac{dE}{dh} = \frac{Mh}{3} - \frac{\varepsilon_{\text{mach}}|f(x_0)|}{h^2} = 0 \;\Rightarrow\; h^* = \left(\frac{3\varepsilon_{\text{mach}}|f(x_0)|}{M}\right)^{1/3}$$
For double precision ($\varepsilon_{\text{mach}}\approx2.2\times10^{-16}$) and $|f|,M=O(1)$, this gives $h^*\sim10^{-5}$ to $10^{-6}$ - far larger than intuition ("as small as the computer allows") would suggest.
Live IEEE double-precision error of the central-difference estimate of $\cos(1)$ from $f(x)=\sin x$, computed in real floating-point arithmetic as $h$ shrinks. Drag left past the minimum (the "V" bottom, around $h\sim10^{-5}$–$10^{-6}$) and watch the error climb back up - pure round-off, exactly as the formula above predicts.
Richardson Extrapolation
Richardson Extrapolation
The Idea - Cancel the Leading Error Term Algebraically
Suppose an approximation $N(h)$ satisfies $N(h) = M + K_1h^2 + K_2h^4 + \cdots$ (true for central difference - only even powers of $h$ appear, by symmetry). Then $N(h/2) = M + K_1h^2/4 + K_2h^4/16+\cdots$. Eliminate the $h^2$ term algebraically:
$$M = \frac{4N(h/2) - N(h)}{3} + O(h^4)$$
This new combination is accurate to $O(h^4)$ - without needing a smaller, round-off-dangerous $h$, and without needing any new formula - just two evaluations of the same $O(h^2)$ formula. Repeating the trick on the new sequence cancels the $h^4$ term too, reaching $O(h^6)$, and so on - precisely the same triangular-table idea as Neville's Method (Chapter 5), applied to step size instead of interpolation nodes.
📘 Example 8.3 - A Richardson Table for $f'(1)$, $f(x)=e^x$
Using central differences at $h=0.1,\,0.05,\,0.025,\,0.0125$ to approximate $f'(1)=e\approx2.7182818285$, build a Richardson extrapolation table.
Solution
The raw $h=0.1$ central-difference estimate has error $4.5\times10^{-3}$. After just one Richardson step ($N_2$ column), error drops to $\sim10^{-8}$ - better than what the five-point formula (Example 8.2) achieves, using only data the central-difference formula already produced, and without ever using an $h$ small enough to risk the round-off blow-up of §8.2.
| $h$ | $N_1$ (central, $O(h^2)$) | $N_2$ ($O(h^4)$) | $N_3$ ($O(h^6)$) | $N_4$ ($O(h^8)$) |
|---|---|---|---|---|
| 0.1 | 2.72281456 | |||
| 0.05 | 2.71941459 | 2.71828126 | ||
| 0.025 | 2.71856499 | 2.71828179 | 2.71828183 | |
| 0.0125 | 2.71835262 | 2.71828183 | 2.71828183 | 2.71828183 |
Practice
Practice Problem
Try this problem yourself before expanding the solution below.
✏️ Practice Problem 8.1 - Richardson Extrapolation for $f'(2)$, $f(x)=\ln x$
Use the centered-difference formula to estimate $f'(2)$ for $f(x)=\ln x$ with $h=0.1$ and again with $h=0.05$. Then apply Richardson extrapolation to combine the two into an $O(h^4)$ estimate, and compare all three to the exact derivative.
Solution
The exact derivative is $f'(2)=1/2=0.500000$. The centered formula $D(h)=\dfrac{f(x_0+h)-f(x_0-h)}{2h}$ gives
$$D(0.1) = 0.500417, \qquad D(0.05) = 0.500104.$$
Both are $O(h^2)$, and the error indeed shrinks by about a factor of 4 when $h$ is halved ($4.17\times10^{-4}\to1.04\times10^{-4}$). Richardson extrapolation eliminates the leading $O(h^2)$ term:
$$R = \frac{4D(0.05)-D(0.1)}{3} = 0.499999834.$$
$D(0.1)=0.500417$ (error $4.2\times10^{-4}$), $D(0.05)=0.500104$ (error $1.0\times10^{-4}$), Richardson $R=0.499999834$ (error $1.7\times10^{-7}$) - two cheap $O(h^2)$ estimates combine into an $O(h^4)$ estimate almost 1000 times more accurate.
Connections Across the Course
- Every finite-difference formula here comes from differentiating an interpolating polynomial from Chapter 5 - the error terms are direct descendants of the interpolation error formula (Theorem 5.2).
- Richardson extrapolation reappears essentially unchanged as Romberg integration in Chapter 10 - same table, same algebra, applied to the trapezoidal rule instead of central differences.
- The optimal-step-size analysis is the same truncation-vs-round-off reasoning first introduced in Chapter 1, now made fully quantitative.