Unlocking the secrets of chaotic planetary systems

Title: Fundamental limits from chaos on instability time predictions in compact planetary systems
Authors: Naireen Hussain, Daniel Tamayo
First Author’s Institution: Department of Astronomy and Astrophysics, University of Toronto, Toronto, ON
Status: Accepted for Publication in MNRAS, preprint on arxiv

It shows up in nearly every field of study – from weather forecasting, to physics, to economics – even sociology – and of course, astronomy. Chaos theory is the study of systems whose seemingly random behavior is the result of an extreme sensitivity to initial conditions. (For an excellent, more in-depth explanation of chaos, check out this astrobite). Chaos is a subject that commonly comes up when trying to understand the long-term stability of planetary systems.

It turns out that certain arrangements of planets are inherently unstable – that is – if you place them in a certain configuration and let them orbit their star for long enough, the gravitational interactions between the planets will fling some (or sometimes all) of the bodies clear out of the system. Unfortunately, determining how and when this will happen is not possible to work out on paper. Or at least, no one has been clever enough to figure it out yet.

Fortunately, computers make this problem somewhat tractable. By gradually evolving a collection of massive bodies over many tiny time steps, it is possible to get an incredibly accurate estimate of where and how these bodies will be moving sometime in the future (or the past, for that matter). Given enough computing power, you can simply take a planetary system and evolve it forward in time and see what happens. Does it stay stable? Do any planets get ejected? Using this technique, astronomers can try placing extra bodies in known planetary systems and see if things remain stable. If not, this sometimes can rule out the presence of additional, undetected planets.

Searching for chaos

As mentioned above, these types of systems are sometimes chaotic. If so, this means, by definition, that the outcome of whether the system is stable not, and how long it takes to become unstable, is highly sensitive to the initial conditions. The authors of today’s paper want to examine how reliable the estimates of instability timescales from these simulations actually are. If the initial conditions are tweaked just slightly, does this timescale change? And if so, is there an underlying pattern?

For this study, the authors ran a large suite of N-body simulations of compact, three planet systems (loosely inspired by the well-known TRAPPIST-1). For each simulation, the orbits and sizes of the planets were varied and the planetary configuration was then evolved for a billion or so orbits to test whether or not the system became unstable. If it did, the timescale for instability was recorded and the same configuration was run again with the initial conditions tweaked ever so slightly to probe the underlying chaotic behavior. If chaos was indeed influencing the outcome, each slight modification to the initial conditions should result in a measurably different instability time.

Figure 1: The two types of instability timescale distributions found by tweaking the initial conditions of the simulations. In a small number of cases the distribution is sharply peaked (left), while in most cases the instability times follow a log-normal shape.

After doing this, the authors found that the distribution of instability times for a given configuration of planets fell into two broad categories. This is shown in Figure 1. In some cases, the distribution was very sharply peaked around a single value. Otherwise, the distribution had a log-normal shape.

The road to dynamical instability

The difference between these two types of results can be explained by comparing the instability timescale to the Lyapunov timescale, which is how long it takes for chaotic behavior to emerge in a given system. For the sharply peaked distribution, the planets become unstable before chaos sets in. This results in an instability time that is not sensitive to slight changes in the initial conditions. For the broadly peaked distributions, chaos occurs well before the instability. Two example sets of simulations are shown in Figure 2, which demonstrate this difference.

Figure 2: The Lyapunov time as a function of the instability time for a number of simulations (top). The symbols indicate whether the instability timescales follow a peaked or a log-normal distribution. A representative example of the change in the orbital properties of a planet in each of the two cases is shown at the bottom for the peaked (left) and log-normal (right) results. The ‘shadow trajectories’ indicate how this quantity changes after the initial conditions are slightly tweaked.

Most interestingly, the log-normal distributions all have a width of ~0.4 dex, regardless of the differences in the initial conditions. The instability time distribution has the same size and shape if the instability time is short or long, or whether the planets are arranged randomly or placed in mean motion resonances. Mean motion resonances occur when the orbital periods are integer multiples of each other. This can act to substantially stabilize or destabilize an orbital configuration, which makes it even more surprising that the instability time distribution shape is not sensitive to this. The remarkable similarity between these distributions across a wide range of configurations is shown in Figure 3. The only requirements here are that the planets start off in a compact configuration and the Lyapunov time be shorter than the instability time.

Figure 3: The width of the log-normal distribution of instability times for a wide number of simulations. The top panel shows little difference between systems whose planets are arranged randomly and those which begin in mean motion resonances. The bottom panel demonstrates that there is little difference between systems that become unstable quickly and systems which take many millions of orbits for instability to occur.

It is not terribly surprising that a seemingly random process can give rise to a reproducible pattern. To quote the authors of the paper, “While individual steps in a drunkard’s random walk might be unpredictable, the cumulative effect of many steps approaches a well-defined statistical distribution.” In addition to indicating that the mean instability time in a chaotic system like this can be estimated by running only a small fraction of the simulations required to fill out the entire probability distribution, it hints at a fundamental underlying truth connecting these results produced by an extremely complicated process. The authors do not attempt to speculate on what this similarity tells us about chaotic planetary systems, but it provides a tantalizing clue about the underlying mechanisms that drive this rather abstruse process.

About Spencer Wallace

I'm a member of the UW Astronomy N-body shop working with Thomas Quinn to study simulations of planet formation. In particular, I'm interested in how this process plays out around M stars, which put out huge amounts of radiation during the pre main-sequence phase and are known to host extremely short period planets. When I'm not thinking about planet formation, I'm an avid hiker/backpacker and play bass for the band Night Lunch.

Discover more from astrobites

Subscribe to get the latest posts to your email.


  1. The obvious question is…has anyone run this test on our own solar system??

    • The key difference with the solar system is that secular resonances become important. If you take a peek in the original paper at section 4.3, the authors do try measuring the instability time distribution for something more solar system-like. To make the simulations more computationally reasonable, the eccentricities the SS planets are artificially enhanced, which greatly shortens the instability timescale. Shrinking the magnitude of this enhancement appears to make the ITD move away from a log-normal shape at some point (see figure 7), although its not immediately clear if this would still be the case with realistic eccentricities.

      This underlying pattern also appears to be broken if the planets are made to be too coplanar, so I definitely don’t think its the case that this is a global phenomenon.


Leave a Reply