Numerical Analysis

The methods of this period rested on replacing a hard function by a polynomial. Newton and Lagrange gave formulas for the polynomial through given points, and from it came rules for integrating and differentiating tabulated data. Intuition said that more points would give a better fit. In 1901 Carl Runge showed that it can give a much worse one. For $1/(1 + 25x^2)$ on $[-1, 1]$, polynomials through equally spaced points oscillate more and more wildly near the ends as the number of points grows. Placing the points more densely near the ends cures it.

Runge's example set the tone for the subject. A method can be exact in principle and useless in practice. Numerical analysis is as much about when methods fail as about the methods themselves.

Electronic computers arrived in the 1940s and changed the problem. A computer stores each number with a fixed number of binary digits, so almost every operation rounds. A single rounding is tiny, but a large computation does millions of them unseen. Early estimates added up the worst possible error at each step and concluded that long computations would be meaningless. James Wilkinson showed a better way to look at it. The computed answer is usually the exact answer to a slightly perturbed problem. If the perturbation is as small as the rounding itself, the method is as good as the data allows. His book of 1963 made this backward error analysis the standard test of a method.

Arithmetic itself was still chaotic. Each manufacturer rounded in its own way, and some machines gave different answers to the same program. William Kahan led the design of a standard for floating point, adopted as IEEE 754 in 1985. It fixed the formats, required every basic operation to be rounded correctly, and defined what happens on overflow or division by zero.

Forgotten rounding can still kill. The Patriot missile battery at Dhahran counted time in tenths of a second in a 24-bit register. Chopped to that length, 0.1 is too small by about $9.5 \times 10^{-8}$. After 100 hours, or 3.6 million ticks, the clock was off by 0.34 seconds, in which a Scud missile travels more than half a kilometre. On 25 February 1991 the battery failed to track an incoming Scud, and 28 soldiers died.

The number $\sqrt{2}$ is the positive root of $f(x) = x^2 - 2$. Newton's method replaces $f$ by its tangent line at the current guess $x_n$ and takes the tangent's root as the next guess:

This is the old Babylonian rule: average the guess with 2 divided by the guess. Start at $x_0 = 1$:

The exponent of the error roughly doubles at each step: $10^{-3}$, $10^{-6}$, $10^{-12}$, $10^{-24}$. This is quadratic convergence. Writing $e_n = x_n - \sqrt{2}$, the formula gives exactly

so each error is about the square of the one before, divided by $2\sqrt{2} \approx 2.8$. Five steps give 24 correct decimal places, more than double precision can store. A method that merely halved the error each step would need about 80 steps for the same accuracy.

The catch is the start. Quadratic convergence holds only close to the root. From a poor guess, Newton's method can wander, cycle or jump to a different root, and for polynomials of degree three and higher the regions of starting points that lead to each root have fractal boundaries. Arthur Cayley first asked about those regions in 1879, and his question is one of the roots of complex dynamics.

By the 1960s numerical analysis had become a mathematical subject with its own journals and a central question: how much can a computed answer be trusted? Its answers split into branches. Solving large systems of linear equations became numerical linear algebra. Solving the equations of physics on grids became numerical methods for partial differential equations. Computing with random samples became Monte Carlo methods, and finding the lowest point of a function became continuous optimisation. All four rest on the same two questions Newton's method raised: how fast does the method converge, and how much does rounding cost along the way?