Finite difference derivatives

Approximating derivatives with finite differences

Whenever the derivative $f'(x)$ of a function is not available in closed form (or is expensive to derive/implement), it can be approximated numerically from function values alone. The idea is to replace the derivative — the slope of the tangent — by the slope of a nearby secant through points that are a small distance $h$ apart. The left plot below shows the function together with the exact tangent and the secants that correspond to each formula. The right plot shows how the approximation error behaves as $h \to 0$.

Forward, backward and central differences

The three basic one-dimensional finite difference formulas are:

$$ \begin{align*} f'(x) &\approx \frac{f(x+h)-f(x)}{h} &&\text{(forward difference)}\\[4pt] f'(x) &\approx \frac{f(x)-f(x-h)}{h} &&\text{(backward difference)}\\[4pt] f'(x) &\approx \frac{f(x+h)-f(x-h)}{2h} &&\text{(central difference)} \end{align*} $$

Truncation error: why central differences are more accurate

The accuracy of each formula follows directly from a Taylor expansion of $f$ around $x$:

$$ \begin{align*} f(x+h) &= f(x) + h f'(x) + \tfrac{h^2}{2} f''(x) + \tfrac{h^3}{6} f'''(x) + O(h^4)\\ f(x-h) &= f(x) - h f'(x) + \tfrac{h^2}{2} f''(x) - \tfrac{h^3}{6} f'''(x) + O(h^4) \end{align*} $$

Inserting the first expansion into the forward difference formula leaves an error term proportional to $h$:

$$\frac{f(x+h)-f(x)}{h} - f'(x) = \tfrac{h}{2} f''(x) + O(h^2) \;\Rightarrow\; \text{error} = O(h).$$

The odd-order terms cancel when the two expansions are subtracted for the central difference, leaving an error term proportional to $h^2$ — one order more accurate for (almost) the same cost:

$$\frac{f(x+h)-f(x-h)}{2h} - f'(x) = \tfrac{h^2}{6} f'''(x) + O(h^4) \;\Rightarrow\; \text{error} = O(h^2).$$

Rounding error: why $h$ cannot be made arbitrarily small

The formulas above suggest choosing $h$ as small as possible. In floating point arithmetic this is not true: $f(x+h)$ and $f(x)$ are nearly equal for small $h$, so their subtraction suffers from catastrophic cancellation. With machine epsilon $\varepsilon \approx 2.2\times10^{-16}$ (IEEE double precision), the rounding error of a difference quotient grows like $\varepsilon/h$, so the total error is roughly

$$\text{error}_\text{forward}(h) \approx \tfrac{h}{2}|f''(x)| + \tfrac{2\varepsilon}{h}|f(x)|, \qquad \text{error}_\text{central}(h) \approx \tfrac{h^2}{6}|f'''(x)| + \tfrac{\varepsilon}{h}|f(x)|.$$

Both expressions have a minimum at an optimal step size: $h^\ast_\text{forward} \approx \sqrt{\varepsilon} \approx 1.5\times10^{-8}$ and $h^\ast_\text{central} \approx \varepsilon^{1/3} \approx 6\times10^{-6}$ (for values of $f$, $f''$, $f'''$ of order 1). This is exactly the characteristic V-shaped curve visible in the right plot: the error first drops with slope $1$ (forward) or $2$ (central) on the log-log plot, then turns around and grows again once rounding error takes over. Try dragging the step size slider all the way to the left to see the curve pass through this minimum.

Some functions are harder than others

The four test functions illustrate different failure modes:

• Smooth polynomial $f(x)=x^3-2x^2+x+1$: a well-behaved baseline. Truncation error follows the theoretical $O(h)$ / $O(h^2)$ slopes closely over many orders of magnitude before rounding error takes over.
• Oscillatory function $f(x)=\sin(15x)$: higher derivatives are large ($|f'''| \le 15^3$), which inflates the constant in front of $h$ resp. $h^2$ in the truncation error. The whole error curve is shifted upward, so a much smaller $h$ is needed to reach the same accuracy as for the polynomial — and if $h$ is chosen too large relative to the oscillation period, the difference quotient can even get the sign of the slope wrong.
• Steep transition $f(x)=\tanh(25x)$: away from $x=0$ the function is almost flat and easy, but inside the narrow transition region the second and third derivatives are huge, so $h$ must be small relative to the width of the transition, not just small in absolute terms. Move $x_0$ towards $0$ and watch the error curve shift upward.
• Non-smooth function (kink) $f(x)=|x|$: at $x_0=0$ the derivative simply does not exist. No choice of $h$ helps here — forward difference always returns $+1$, backward difference always returns $-1$, and central difference always returns $0$, independent of $h$. This is a fundamentally different failure than truncation or rounding error: finite differences silently produce a number even where no derivative exists, which is an important pitfall to know about.

References

• [LeV07] Randall J. LeVeque. Finite Difference Methods for Ordinary and Partial Differential Equations. SIAM, 2007.