Hanno Rein

REBOUND: an open-source multi-purpose N-body code for collisional dynamics

REBOUND is a new multi-purpose N-body code which is freely available under an open-source license. It was designed for collisional dynamics such as planetary rings but can also solve the classical N-body problem. It is highly modular and can be customized easily to work on a wide variety of different problems in astrophysics and beyond. REBOUND comes with three symplectic integrators: leap-frog, the symplectic epicycle integrator (SEI) and a Wisdom-Holman mapping (WH). It supports open, periodic and shearing-sheet boundary conditions. REBOUND can use a Barnes-Hut tree to calculate both self-gravity and collisions. These modules are fully parallelized with MPI as well as OpenMP. The former makes use of a static domain decomposition and a distributed essential tree. Two new collision detection modules based on a plane-sweep algorithm are also implemented. The performance of the plane-sweep algorithm is superior to a tree code for simulations in which one dimension is much longer than the other two and in simulations which are quasi-two dimensional with less than one million particles. In this work, we discuss the different algorithms implemented in REBOUND, the philosophy behind the codes structure as well as implementation specific details of the different modules. We present results of accuracy and scaling tests which show that the code can run efficiently on both desktop machines and large computing clusters.

ias15: a fast, adaptive, high-order integrator for gravitational dynamics, accurate to machine precision over a billion orbits

We present IAS15, a 15th-order integrator to simulate gravitational dynamics. The integrator is based on a Gaus-Radau quadrature and can handle conservative as well as non-conservative forces. We develop a step-size control that can automatically choose an optimal timestep. The algorithm can handle close encounters and high-eccentricity orbits. The systematic errors are kept well below machine precision and long-term orbit integrations over 10 9 orbits show that IAS15 is optimal in the sense that it follows Brouwer’s law, i.e. the energy error behaves like a random walk. Our tests show that IAS15 is superior to a mixed-variable symplectic integrator (MVS) and other high-order integrators in both speed and accuracy. In fact, IAS15 preserves the symplecticity of Hamiltonian systems better than the commonly-used nominally symplectic integrators to which we compared it. We provide an open-source implementation of IAS15. The package comes with several easy-to-extend examples involving resonant planetary systems, Kozai-Lidov cycles, close encounters, radiation pressure, quadrupole moment, and generic damping functions that can, among other things, be used to simulate planet-disc interactions. Other non-conservative forces can be added easily.

On the evolution of mean motion resonances through stochastic forcing: fast and slow libration modes and the origin of HD 128311

Aims. We clarify the response of extrasolar planetary systems in a 2:1 mean motion commensurability with masses ranging from the super Jovian range to the terrestrial range to stochastic forcing that could result from protoplanetary disk turbulence. The behaviour of the different libration modes for a wide range of system parameters and stochastic forcing magnitudes is investigated. The growth of libration amplitudes is parameterized as a function of the relevant physical parameters. The results are applied to provide an explanation of the configuration of the HD 128311 system. Methods. We first develop an analytic model from first principles without making the assumption that both eccentricities are small. We also perform numerical N-body simulations with additional stochastic forcing terms to represent the effects of putative disk turbulence. Results. We isolate two distinct libration modes for the resonant angles. These react to stochastic forcing in a different way and become coupled when the libration amplitudes are large. Systems are quickly destabilized by large magnitudes of stochastic forcing but some stability is imparted should systems undergo a net orbital migration. The slow mode, which mostly corresponds to motion of the angle between the apsidal lines of the two planets, is converted to circulation more readily than the fast mode which is associated with oscillations of the semi-major axes. This mode is also vulnerable to the attainment of small eccentricities which causes oscillations between periods of libration and circulation. Conclusions. Stochastic forcing due to disk turbulence may have played a role in shaping the configurations of observed systems in mean motion resonance. It naturally provides a mechanism for accounting for the HD 128311 system for which the fast mode librates and the slow mode is apparently near the borderline between libration and circulation.

whfast: a fast and unbiased implementation of a symplectic Wisdom–Holman integrator for long-term gravitational simulations

We present WHFast, a fast and accurate implementation of a Wisdom-Holman symplectic integrator for long-term orbit integrations of planetary systems. WHFast is significantly faster and conserves energy better than all other Wisdom-Holman integrators tested. We achieve this by significantly improving the Kepler-solver and ensuring numerical stability of coordinate transformations to and from Jacobi coordinates. These refinements allow us to remove the linear secular trend in the energy error that is present in other implementations. For small enough timesteps we achieve Brouwers law, i.e. the energy error is dominated by an unbiased random walk due to floating-point round-off errors. We implement symplectic correctors up to order eleven that significantly reduce the energy error. We also implement a symplectic tangent map for the variational equations. This allows us to efficiently calculate two widely used chaos indicators the Lyapunov characteristic number (LCN) and the Mean Exponential Growth factor of Nearby Orbits (MEGNO). WHFast is freely available as a flexible C package, as a shared library, and as an easy-to-use python module.

Period ratios in multiplanetary systems discovered by Kepler are consistent with planet migration

The Kepler planet candidates are an interesting testbed for planet formation scenarios. We present results from N-body simulations of multi-planetary systems that resemble those observed by Kepler. We add both smooth (Type I/II) and stochastic migration forces. The observed period ratio distribution is inconsistent with either of those two scenarios on its own. However, applying both stochastic and smooth migration forces to the planets simultaneously results in a period ratio distribution that is similar to the observed one. This is a natural scenario if planets form in a turbulent proto-planetary disk where these forces are always present. We show how the observed period ratio and eccentricity distribution can constrain the relative strength of these forces, a parameter which has been notoriously hard to predict for decades. We make the source code of our simulations and the initial conditions freely available to enable the community to expand this study and include effect other than planetary migration.

SMACK: A NEW ALGORITHM FOR MODELING COLLISIONS AND DYNAMICS OF PLANETESIMALS IN DEBRIS DISKS

We present the Superparticle-Method/Algorithm for Collisions in Kuiper belts and debris disks (SMACK), a new method for simultaneously modeling, in three dimensions, the collisional and dynamical evolution of planetesimals in a debris disk with planets. SMACK can simulate azimuthal asymmetries and how these asymmetries evolve over time. We show that SMACK is stable to numerical viscosity and numerical heating over 107 yr and that it can reproduce analytic models of disk evolution. We use SMACK to model the evolution of a debris ring containing a planet on an eccentric orbit. Differential precession creates a spiral structure as the ring evolves, but collisions subsequently break up the spiral, leaving a narrower eccentric ring.

Traditional formation scenarios fail to explain 4:3 mean motion resonances

At least two multi-planetary systems in a 4:3 mean motion resonance have been found by radial velocity surveysy. These planets are gas giants and the systems are only stable when protected by a resonance. Additionally the Kepler mission has detected at least 4 strong candidate planetary systems with a period ratio close to 4:3. This paper investigates traditional dynamical scenarios for the formation of these systems. We systematically study migration scenarios with both N-body and hydrodynamic simulations. We investigate scenarios involving the in-situ formation of two planets in resonance. We look at the results from nely tuned planet-planet scattering simulations with gas disk damping. Finally, we investigate a formation scenario involving isolation-mass embryos. Although the combined planet-planet scattering and damping scenario seems promising, none of the above scenarios is successful in forming enough systems in 4:3 resonance with planetary masses similar to the observed ones. This is a negative result but it has important implications for planet formation. Previous studies were successful in forming 2:1 and 3:2 resonances. This is generally believed to be evidence of planet migration. We highlight the main dierences between those studies and our failure in forming a 4:3 resonance. We also speculate on more exotic and complicated ideas. These results will guide future investigators toward exploring the above scenarios and alternative mechanisms in a more general framework.

The dynamical origin of the multi-planetary system HD 45364

The recently discovered planetary system HD 45364, which consists of a Jupiter and Saturn-mass planet, is very likely in a 3:2 mean motion resonance. The standard scenario for forming planetary commensurabilities is convergent migration of two planets embedded in a protoplanetary disc. When the planets are initially separated by a period ratio larger than two, convergent migration will most likely lead to a very stable 2:1 resonance. Rapid type III migration of the outer planet crossing the 2:1 resonance is one possible way around this problem. In this paper, we investigate this idea in detail. We present an estimate of the required convergent migration rate and confirm this with N-body and hydrodynamical simulations. If the dynamical history of the planetary system had a phase of rapid inward migration that forms a resonant configuration, we predict that the orbital parameters of the two planets will always be very similar and thus should show evidence of that. We use the orbital parameters from our simulation to calculate a radial velocity curve and compare it to observations. Our model provides a fit that is as good as the previously reported one. The eccentricities of both planets are considerably smaller and the libration pattern is different. Within a few years, it will be possible to observe the planet-planet interaction directly and thus distinguish between these different dynamical states.

Planet–disc interaction in highly inclined systems

We study the interaction of a proto-planetary disk and a planet on a highly inclined orbit in the linear regime. The evolution of the planet is dominated by dynamical friction for planet masses above several Earth-masses. Smaller planets are dominated by aerodynamic drag, especially for very high inclinations and retrograde orbits. The time-scales associated with migration and inclination damping are calculated. For certain values of the inclination, the inclination damping time-scale is longer than the migration time-scale and the disk lifetime. This result shows that highly inclined planets can not (re-)align with the proto-planetary disk. We discuss the dependence of numerical simulations on the gravitational softening parameter. We find only a logarithmic dependence, making global three dimensional simulations of this process computationally feasible. A large fraction of Hot Jupiters is on highly inclined orbits with respect to the rotation axis of the star. On the other hand small-mass planetary systems discovered by the Kepler mission have low mutual inclinations. This shows that there are two distinct formation mechanisms at work. The process that creates inclined Hot Jupiters does not operate on small mass planets because the damping timescales are so long that these systems would still be inclined today.

The validity of the super-particle approximation during planetesimal formation

The formation mechanism of planetesimals in protoplanetary discs is hotly debated. Currently, the favoured model involves the accumulation of meter-sized objects within a turbulent disc, followed by a phase of gravitational instability. At best , one can simulate a few million particles numerically as opposed to the several trillion meter-sized particles expected in a real protopl anetary disc. Therefore, single particles are often used as super-partic les to represent a distribution of many smaller particles. I t is assumed that small-scale phenomena do not play a role and particle collisions are not modelled. The super-particle approximation is not always valid when applied to planetesimal formation because the system can be marginally collisional (of order one collision per particle per orbit). The super-particle approximation can only be valid in a collisionless or strongly collisional system, althoug h, in many recent numerical simulations this is not the case. In this work, we present new results from numerical simulations of planetesimal formation via gravitational instabili ty. A scaled system is studied that does not require the use of super-part icles. This system is simplified for computational practica lity and proper identification of important processes: 1) the evolution of p articles is studied in a local shearing box; 2) the particle- particle interactions such as gravity, physical collisions, and gas drag are solve d directly assuming a constant background shear flow without any feedback from the particles. We find that the scaled particles can be us ed to model the initial phases of clumping if the properties of the scaled particles are chosen such that all important timescales in t he system are equivalent to what is expected in a real protoplanetary disc. Constraints are given for the number of particles needed in order to achieve numerical convergence. We compare this new method to the standard super-particle approach. We find that the super-particle approach produces un reliable results that depend on artifacts such as the gravitational s oftening in both the requirement for gravitational collaps e and the resulting clump statistics. Our results show that short-range intera ctions (collisions) have to be modelled properly.