The wrong orbit beside a right one

The error that doubles on a schedule ends with a promise it does not keep. Two double pendulums released a hundred-millionth of a degree apart follow one curve for ten seconds and then have nothing to do with each other, and it notes, in passing, that the same is true of every computer simulation of a chaotic system — and that such a simulation is nevertheless “usually close to some true trajectory of the same system, just not the one it was asked for”. That sentence carries a great deal of weight. Weather forecasts, simulations of molecules in a liquid, integrations of the solar system over billions of years and every numerical study of turbulence are computations of chaotic systems, run for far longer than their computed trajectories stay near the true ones. If the sentence were false, all of them would be pictures of nothing.

It is a theorem, for a large and important class of systems, and it is called the shadowing lemma. Its proof is a construction, and the construction is simple enough to draw.

One equation, written two ways #

The logistic map $x \to rx(1-x)$ is the one the dripping tap turns out to be, and at $r = 4$ it is fully chaotic: nearby points separate by a factor of two per step on average, a Lyapunov exponent of $\ln 2$. The experiment in the drawing does not even perturb the starting point. It computes the same function twice, written in two algebraically identical forms, and the only difference between the two computations is where the rounding falls. The rounding is a part in $10^{16}$. By step 52 it is everything.

The arithmetic of that is worth doing once. An error of $10^{-16}$ that doubles every step reaches order one after $\log_2 10^{16} \approx 53$ steps. Computing in quadruple precision, with errors of $10^{-34}$, would stretch that to about 113 steps; each extra decimal digit of precision buys about three and a third more steps. No attainable precision buys a thousand. After a few dozen steps, a computed orbit of this map is not an approximation to the orbit from its starting point. It is an unrelated sequence of numbers.

The computer’s doubling map ends at zero #

The doubling map, $x \to 2x \bmod 1$, is the simplest chaotic map there is — the logistic map at $r = 4$ is the doubling map seen through the substitution $x = \sin^2\pi\theta$ — and in floating-point arithmetic it does something that looks worse than anything in the first drawing. From every starting point, it goes to zero and stays there.

This time there is no rounding at all. Multiplying a binary number by two is exact: it shifts the bits one place left. Dropping the integer part is exact: it discards the bit that crossed the binary point. A double-precision number carries 53 significant bits, and after 53 doublings none are left. The computed orbit is not an approximation to anything; it is exactly the orbit of the number the computer stored. The trouble is that every number a computer stores is a fraction whose denominator is a power of two, and those are exactly the starting points whose orbits end at zero. They are a set of measure zero — a starting point picked at random from the real interval lands on one with probability zero — and they are the only starting points a computer can represent.

The lesson is sharper than “rounding errors grow”. An exactly correct computation of a chaotic map can be completely unrepresentative, because the computer can only ever reach an atypical corner of the space of initial conditions. Whatever makes chaotic simulation trustworthy, it cannot be the accuracy of each step.

The noisy orbit a true orbit stays beside #

Here the noise is deliberately large, a thousandth at every step, so that it can be seen. The red curve compares the noisy computation with the exact orbit from the same starting point, and it is what the first drawing already showed: the two part within a handful of steps and never come back.

The blue curve is a different exact orbit — the shadow — and it stays beside the computation for all four hundred steps, never further away than the noise added at a single step. It is found by a construction that uses the chaos rather than fighting it. Run the doubling map forwards and every error is doubled at every step. Run it backwards — halve, and add back the integer the forward step removed — and every error is halved. So start at the far end of the noisy computation, take its last point as the shadow’s last point, and step backwards: at each step the shadow is pulled towards the computed point by the halving, and the difference between them can never exceed the error made at a single step. At the start of the computation the shadow is at a slightly different place from where the computation began — by the time the backwards run gets there, the halvings have fixed where it has to be.

That is the whole of the shadowing lemma for this map. The noisy computation is not an approximation to the orbit it started on, but it is uniformly close to an exact orbit that started nearby. Every point it produces is a faithful picture of where a real orbit of the map goes, and every statistical question about real orbits — how often they visit a region, how fast they separate, what their averages are — can be answered from the computation. The only question it cannot answer is the one a forecaster most wants: where this orbit, from this starting point, goes.

Stretching errors undone backwards, squeezing ones forwards #

The doubling map only stretches. Most chaotic systems stretch in some directions and squeeze in others: a patch of initial conditions in the fold that has to be there is drawn out along one direction and flattened along another, which is how a bounded region can hold exponentially separating trajectories. For such a system the backwards construction cannot work alone, because running backwards now stretches the squeezed directions.

The cat map — a square stretched into a parallelogram and cut and reassembled back into the square, like a slice of toffee pulled and folded — is the cleanest example of both at once. The fix is to split every error into its two components and correct each one in the direction in which it shrinks. The component along the stretching direction is corrected backwards from the end of the computation, as for the doubling map. The component along the squeezing direction is corrected forwards from the start, since running forwards squeezes it. Both corrections stay small, bounded by the noise per step divided by how much the map contracts in the direction concerned, and together they define an exact orbit that stays within that bound of the computation at every step.

This is what the word hyperbolic means in dynamics, and why it matters. A system is hyperbolic if at every point the directions split cleanly into ones that are stretched and ones that are squeezed, by rates bounded away from one, and the split varies smoothly as the orbit moves. Dmitri Anosov and Rufus Bowen proved in the 1960s and 1970s that such systems can always be shadowed: every sufficiently accurate computation, however long, has an exact orbit beside it, at a distance proportional to the error per step.

Every step wrong, the histogram right #

The payoff of shadowing is statistical. The logistic map at $r = 4$ has an exactly known invariant density — the fraction of time a typical orbit spends near each point — which piles up near the ends of the interval, because the map is slow there. Two hundred thousand computed steps, not one of which after the first fifty resembles the orbit from its own starting point, reproduce that density within the scatter expected from that many samples. The computation is a faithful sample of how the system behaves on average.

The logistic map is not uniformly hyperbolic — at $x = \tfrac12$ it neither stretches nor squeezes, and orbits passing close to that point lose the protection — so the drawing is evidence rather than a consequence of the theorem. That it works anyway is typical of a wide class of chaotic systems that are “hyperbolic almost everywhere”, for which shadowing can be established for very long times, if not for ever, by checking the computed orbit itself.

This is also what gives the time averages of the equation that only runs forwards their authority. A molecular-dynamics simulation of a liquid follows a trajectory that is wrong within picoseconds and computes a diffusion coefficient, a pressure and a heat capacity that agree with experiment. What it computes is an average over the invariant measure, and shadowing, where it holds, is the guarantee that the computed trajectory samples that measure as a true trajectory would.

Finding a shadow for a computation that has already run #

For the doubling map and the cat map the shadow can be written down, because the stretching and squeezing directions are the same everywhere. For a general system they turn from point to point along the orbit, and the construction has to be done numerically — which raises the obvious worry that the shadow is being computed with the same untrustworthy arithmetic as the orbit it is supposed to vindicate.

The way round it is to compute the shadow with the errors made explicit. Take the computed orbit as given. Recompute each single step in much higher precision, which is cheap for one step, and record how far the computed next point misses the true image of the current one: those misses are the noise in the first drawing, measured rather than assumed. Then solve for a correction to every point at once, using the linearised map at each point to split the correction into stretching and squeezing parts, and correcting the first backwards and the second forwards exactly as the cat map’s shadow was built. Repeat until the corrections stop changing. What remains is an orbit whose every step obeys the map to within a known tiny error, lying within a known distance of the original computation — and with interval arithmetic, which carries rigorous bounds on every rounding, the result is a proof that a true orbit lies within that distance.

That procedure has been carried out for orbits of millions of steps of the Hénon and logistic maps and for flows like the Lorenz system. It is expensive, and it answers the question for one computation at a time, but it turns “the simulation is probably representative” into a checked statement.

A forecast and a shadow are different questions #

Weather forecasting makes the distinction concrete. The atmosphere’s errors grow with a doubling time of a day or two at the scales that matter for a forecast, so a single forecast from the best available analysis of today’s weather loses its skill after one to two weeks, however good the model. No shadowing theorem changes that: the forecast is being asked about the orbit that starts at the observed state, and that is exactly the question a shadow cannot answer.

What forecasting centres do instead is lean on the statistical side. They run an ensemble — typically fifty or so forecasts from slightly different starting states, each consistent with the observations to within their errors — and read the spread of the ensemble as a measure of how predictable the weather currently is. The procedure quietly assumes something shadowing-like: that each member of the ensemble, although wrong as a prediction, is a plausible evolution of a real atmosphere, so that the ensemble samples the range of real possibilities rather than a range of numerical artefacts. When the members agree for ten days, the atmosphere is in a predictable state; when they fan out after three, it is not, and the forecast says so.

The solar system gives a second example, on a longer clock. Integrations of the planets’ orbits by Jacques Laskar and others since 1989 find the inner solar system chaotic, with a Lyapunov time of about five million years, so no integration can say where Mercury will be in a billion years. What the integrations can do is run thousands of neighbouring solar systems forward for five billion years and count how many end with Mercury’s orbit elongated enough to cross Venus’s — its eccentricity, which a lone planet round a lone sun would keep fixed for ever, pumped up by the slow chaotic pull of the others. The answer, about one in a hundred, is a statistic of exactly the kind shadowing protects — and it is the only kind of statement about the far future of the planets that a chaotic system allows anyone to make.

The same trust underlies the most famous numerical experiment in the subject. The energy that refuses to be shared describes how Fermi, Pasta, Ulam and Tsingou’s computer, in 1955, returned almost all the energy of a chain to its first mode, when every expectation was that it would spread. Whether a computed recurrence of that kind was physics or arithmetic was a real question, and it was settled by recomputing with different arithmetic and getting the same answer — the practical form of asking whether a true orbit lies beside the computed one.

Where shadowing fails #

Non-hyperbolic points. A map that neither stretches nor squeezes somewhere — the logistic map at its maximum, the flat parts of the last curve to go in the standard map, where chaotic regions border islands of regular motion — gives the backwards construction nothing to work with there. Computed orbits passing close to such points can separate from every true orbit, and near a point where the dynamics is only marginally unstable they can linger, as they do in the calm that is the ghost of a cycle, long enough for the arithmetic’s errors to matter: in 1987 Stephen Hammel, James Yorke and Celso Grebogi showed, for the logistic map and the Hénon map, that shadows exist for very long stretches but are interrupted by “glitches” near such points, and that the length of time a computation can be shadowed depends on how accurate it is.

Changing numbers of unstable directions. Worse, in some systems the number of stretching directions itself varies from place to place along an orbit — two unstable directions here, one there. Around 1994 Stephen Dawson and colleagues showed that such systems can be unshadowable in a strong sense: no true orbit stays near a computed one for longer than a short time, however accurate the arithmetic. High-dimensional chaotic systems, including models of the atmosphere, are suspected of having this property.

Averages that are not typical. Shadowing guarantees that the computation is near some true orbit, not that the true orbit is typical. The doubling map’s march to zero is the extreme case: a shadowed, even exact, computation that samples nothing like the invariant measure. Most practical computations are protected against this by the rounding itself, which behaves like small random noise; the doubling map is dangerous exactly because its arithmetic is too exact.

What the drawings cannot show #

The drawings are one- and two-dimensional maps chosen because their shadows can be constructed exactly. What they cannot show is whether a particular large computation — a climate model, a simulation of a turbulent flow — is shadowed, which nobody can currently check. Such systems are almost certainly not uniformly hyperbolic, and the evidence that their statistics are trustworthy comes from convergence tests and comparison with observation rather than from a theorem.

Nor do the drawings show time, in the sense a forecaster needs. A shadow exists beside the whole computation; it says nothing about how long the computation stays near the orbit that began at the observed starting point. That is still governed by the Lyapunov exponent, which the error that doubles on a schedule measures, and no amount of shadowing extends a forecast past it.

Still open: which chaotic systems can be trusted #

Hyperbolic systems are shadowable, and they are rare among the systems physics cares about. Most real chaotic dynamics has regular islands mixed into chaotic seas, directions that are only sometimes unstable, and attractors whose structure varies from place to place. For these, whether long computations are faithful is known case by case: some are shadowable for times far longer than any computation, some are not shadowable at all, and there is no general criterion that says which, applicable to a system of the size of a weather model. A related open question is whether the statistics of such a system are robust — whether a small change in the equations, of the kind every numerical model makes, changes the long-run averages by a small amount or by a large one. For hyperbolic systems it is small; for many others it is not known, and it is part of the question of how far climate projections from different models should be expected to agree.

The habit worth carrying away is to separate a trajectory from a statistic. A chaotic computation is worthless as a forecast past a few Lyapunov times and can still be exact as a sample, because the errors that destroy the forecast are stretched along the very directions the system explores anyway. What must be asked of it is not “is this the right orbit” — it never is — but “is there a right orbit beside it”, and for a hyperbolic system the answer can be constructed.