Ruination by rounding

The following is an extract from
Seven Deadly Sins of Numerical Computation, by Brian J. McCartin

Interacting with a digital computer can be an unsettling affair. Despite all our school years spent studying the real number continuum, we are rudely awakened to the fact that we must, in practice, content ourselves with some finite set of numbers. This consequence of finite-precision arithmetic introduces roundoff error into numerical computations (even upon input!). Let us explore, via example, how this affects traditional mathematical arguments based upon properties of the real number system.

All of the commonly used numerical differentiation rules possess error estimates of the form

$$D(f) = f^\prime(a) = D_h + ch^r + o(h^r)$$ (1)

where $D_h$ is our finite difference approximation, $c$ is the aymptotic error constant, $h$ is the mesh width, and $r$ is the order of approximation. However, such estimates assume infinite precision arithmetic, which is not available on a fixed-word-length computer.

The portion of the error account for in (1) is called truncation or discretization error. The additional error due to finite word length is referred to as the roundoff error. The difficulty lies in the fact that, while truncation error goes to zero with $h$, the roundoff error becomes unbounded. Consequently, error terms such as that appearing in (1) present the illusion that unlimited accuracy is achievable simply by refining the mesh.

For the sake of definiteness, consider the central difference formula

$$D_h = \frac{ f(a + h) - f(a-h) }{2h}$$ (2)

with

$$c = -f^{\prime\prime\prime}(a)/6, r=2$$

If we account only for the roundoff error, $E_\pm$, incurred by function evaluations then we obtain a computed value of

$$f^\prime_c = \frac{ f(a+h)+E_+ - f(a-h) - E_-}{2h} = D_h + \frac{ E_+ - E_-}{2h}$$

$$= D(f) + \frac{ E_+ -E_-}{2h} + ch^2 + o(h^2)$$ (3)

Thus, we determine that the actual error is composed of two quite dissimilar pieces; a parabolic truncation error and a hyperbolic roundoff error. As seen in figure, this implies the existence of an optimum mesh width below which the total error increases (but the precise determination of this optimum is impossible in practice).

In practice we probably dont have to worry about going below the optimum mesh size since the precision of modern computers is quite high. If you feel that this issue is important in some situation then it is better to use double precision for your computations.