Run the simulation

Once the JSON file is created, we are ready to start our N-body simulation. Posidonius implements the symplectic integrator WHFAST (Rein & Tamayo, 2015) without correctors (since tidal forces are velocity dependent). This algorithm, based on the Wisdom-Holman integrator but carefully designed with finite double floating-point precision problems in mind, has been proven to be fast and accurate for long-term orbit integrations of planetary systems. It assumes a single gravitationally dominating body (e.g., exoplanetary system with a single star). Its symplectic nature (i.e., simulations obeys Hamilton's equations to a high degree of accuracy and the energy is conserved) is formally lost when velocity dependent forces are included (such as tidal forces). It is unbiased (i.e., the errors are random and uncorrelated), and has a very slow error growth. For sufficiently small time steps, it achieves Brouwer's law (i.e., the energy error behaves like a random walk).

The integration with WHFAST can take place in three different coordinates (Hernandez and Dehnen, 2017):

  • Jacobi coordinates: it leads to a better precision if orbits are well separated and do not cross each other.
  • Democratic-Heliocentric coordinates (a.k.a. canonical heliocentric, Poincare or mixed-variables coordinates): better precision if there are close orbits and/or encounters.
  • WHDS: it splits the Hamiltonian into three parts and solves the two body problem exactly, contrary to the democratic heliocentric one.

The default coordinates can be changed from the Python script by explicitly indicating the desired reference system:

    output_filename = "target/custom.json"
    whfast_alternative_coordinates="DemocraticHeliocentric"
    whfast_alternative_coordinates="WHDS"
    whfast_alternative_coordinates="Jacobi"
    universe.write(output_filename, integrator="WHFast", whfast_alternative_coordinates=whfast_alternative_coordinates)

Posidonius also implements the IAS15 integrator (Rein & Spiegel, 2015), which can handle conservative and non-conservative forces. Its time step size is automatically adapted to the most optimal one, it handles close encounters and high-eccentricity orbits. It does not assume that there is a single massive body that dominates the evolution (e.g., exoplanetary system with a single star), it preserves the symplecticity of Hamiltonian systems even better than symplectic integrators but it is significantly slower. It also follows Brouwer's law. To select his integrator:

    output_filename = "target/custom.json"
    universe.write(output_filename, integrator="IAS15")

Posidonius WHFAST and IAS15 implementations are equivalent to the REBOUND code (Rein & Liu, 2011), in WHFAST case it is equivalent to the safe mode (required by the tidal effects) and without correction (i.e. 2nd order integrator, comparable to mercury symplectic part of the hybrid integrator).

In a consistent way, the spin evolution is computed for all the bodies by using a midpoint integrator (Stoer & Bulirsch, 2013).

Start the simulation of a JSON case

The simulation can be started using the JSON file (which contain a description of the initial conditions). When starting a simulation, the recovery and historic snapshot file names should be specified. The former will contain the information needed to resume interrupted simulations, while he latter stores the evolution of the simulation over the years.

posidonius start target/case3.json target/case3.bin target/case3_history.bin
posidonius start target/case4.json target/case4.bin target/case4_history.bin
posidonius start target/case7.json target/case7.bin target/case7_history.bin
posidonius start target/example.json target/example.bin target/example_history.bin

The recovery snapshot will be backed up every 12 hours, thus for example files named like case3.20171126TAM.bin and case3.20171126TPM.bin will appear in the same directory where case3.bin is saved.

The flag --silent can be added to avoid printing the current year of the simulation. An execution time limit can also be specified with the flag --limit, this can be useful for supercomputers that only allow processes to last a given amount of real time (not simulation time).

Resume an interrupted simulation

Interrupted simulations can be restored using the recovery snapshot file. The historic snapshot filename has to be specified also to continue storing the history of the simulation.

posidonius resume target/case3.bin target/case3_history.bin
posidonius resume target/case4.bin target/case4_history.bin
posidonius resume target/case7.bin target/case7_history.bin
posidonius resume target/example.bin target/example_history.bin

The flag --silent can be added to avoid printing the current year of the simulation. n execution time limit can also be specified with the flag --limit, this can be useful for supercomputers that only allow processes to last a given amount of real time (not simulation time). In case the user wants to change historic or recovery snapshot periods when resuming a simulation, it can be done with the flags --historic-snapshot-period and --recovery-snapshot-period plus the new period (in days) after each one. The flag --time-limit can be used to change the simulation time limit (e.g., to increase it for a previous short simulation that looks promising).