Intriguing title, no? These are the first eleven words of Neal Stephenson’s novel Seveneves, which set up the remaining 600 pages as an extended treatise on the future of humanity as it copes with certain annihilation. I thoroughly recommend it, as long as you can deal with hundreds of pages of orbital mechanics. In this post I will numerically explore this post-lunar age, to verify for myself if it would be as deadly as described.

In the novel, one day the moon breaks up into 7 roughly equal-sized pieces. These pieces continue peacefully orbiting the Earth for a while, and eventually two pieces collide. This collision causes a piece to fragment, making future collisions more likely. The process repeats, at what Stephenson says is an exponential rate, until the Earth is under near-constant bombardment from meteorites, wiping out (nearly) all life on Earth.

How likely is this? Let’s simulate the process numerically.

Simulation design

At the beginning of the simulation, there are 7 masses orbiting Earth in near circular orbits. In this model I include

A static mass at the origin of a 2D co-ordinate system – this is Earth.

Gravitational interactions between all bodies. I use an extraordinarily naive direct summation approach, which scales as and is horrible. A huge improvement would be to approximate forces from distant bodies with some kind of multipole expansion, and store masses in a balanced tree – this is known as a Barnes-Hut simulation.

and is horrible. A huge improvement would be to approximate forces from distant bodies with some kind of multipole expansion, and store masses in a balanced tree – this is known as a Barnes-Hut simulation. Pairwise collisions between bodies – most of the collisions are approximated as elastic, for which the post-collision velocities are easily calculated.

Possible fragmentation of bodies – in this case one of the masses breaks into 2 equally-sized pieces, which fly off into different directions. Energy is conserved by subtracting kinetic energy from the other body in the collision. Additionally, I set a kinetic energy threshold for the fragmentation, such that smaller bodies have to collide at larger relative velocities for a fragmentation to occur.

for a fragmentation to occur. Energy-conserving integration – the temporal integration is performed with a Verlet (leapfrog) technique. This is important as it is the simplest of a class of numerical integrators called sympletic integrators, which have inherently bounded errors. (This is because they preserve volume in the phase space of the particle dynamics, and so are not susceptible to damping or growth of energy). Note that this is not true for simple integrators like the Euler or Runge-Kutta methods.

Simulation results

Let’s have a look at what this setup looks like – the case after a couple of orbits is plotted below. The Earth is at the origin, and the moon fragments start at a position (1,0), orbiting in a circular orbit in the anti-clockwise direction.

I’ve set the ‘radius’ of Earth to be one fifth of the initial moon distance. This is much larger than reality, but it helped keep the simulation short. With my terrible and slow code, I needed to help everything evolve a bit more quickly!

You can see the fragments have already started to interact, spreading out in the transverse and forward/rearward directions. There have also been a couple of fragmentations.

A video of the entire simulation is below. On the left is a close-up view of the Earth, on the right is a running tally of meteorite impacts. An impact is defined as when a mass collides with the circle at the centre, and the mass is removed from the simulation. Every impact, the outline of Earth flashes red.

As you can see, things start out gently but suddenly the number of impacts skyrockets, until eventually tailing off. This cumulative plot is expanded below.

One of the events depicted in the book is the ‘white sky’, where the number of impacts grows exponentially until they are near-constant. Plotted on a logarithmic scale, there does appear to be an initial exponential phase in the simulation, lasting a relatively short time. Soon after, the impact rate slows but not to zero, with impacts occurring every so often for the remainder of the simulation. This is another scenario predicted in the book, where the planet takes many thousands of years to cool down after the initial heavy bombardment.

However, as time goes on more and more of these collisions should be by lighter fragments, due to the increased fragmentation. The distribution of masses in orbit evolves as follows.

Because collisions must happen at higher and higher relative velocities for smaller fragments, the mass distribution bottoms out here at around a ten-thousandth of the initial fragment mass. (This was again a deliberate decision in order to reduce the total number of fragments in the simulation.) After the initial burst of fragmentation, the population remains relatively stable. Continuing collisions don’t often cause fragmentation, but do direct fragments towards the Earth for a long time.

As expected then, the mass of impacts drops with time in a very similar fashion:

We can also track where the mass in the simulation ends up. As plotted below, over 25% ends up on the surface of the Earth! This is an enormous amount, the energy in gravitational potential alone is > J, or the equivalent of 10 billion Tsar Bombas.

The rest of the mass, liberated by the spent gravitational potential energy of the rest of the fragments, is sent spinning out into space. We can observe this by plotting the distribution of distances of the fragments from Earth:

The fragments hitting Earth eventually reach the dotted line, but you can also see the bulk of the population is rapidly speeding away. (In fact, in the simulation any fragment which reaches a large enough distance above the escape velocity there is removed from the simulation to save on computational time.)

Next steps

Writing this simulation and producing a nice video was a decent challenge, but there is much to improve. First is optimising the simulation to handle more than a thousand or so fragments, either by switching to a more efficient algorithm or parallelising with CUDA or similar. Second is to use a more realistic collision and fragmentation scenario, allowing for inelastic collisions and using a physically relevant model for the masses of daughter fragments.

However, I think this simple model is enough to add some realism to the scenario Neal Stephenson wrote about. It also makes me glad the moon is still in one piece. Last I checked.