Search

We present trace , an almost time-reversible hybrid integrator for the planetary N -body problem. Like hybrid symplectic integrators, trace can resolve close encounters between particles while retaining many of the accuracy and speed advantages of a fixed time-step symplectic method such the Wisdom–Holman map. trace switches methods time-reversibly during close encounters following the prescription of Hernandez & Dehnen. In this paper we describe the derivation and implementation of trace and study its performance for a variety of astrophysical systems. In all our test cases, trace is at least as accurate and fast as the hybrid symplectic integrator mercurius. In many cases, trace 's performance is vastly superior to that of mercurius. In test cases with planet–planet close encounters, trace is as accurate as mecurius with a 12× speed-up. If close encounters with the central star are considered, trace achieves good error performance while mercurius fails to give qualitatively correct results. In ensemble tests of violent scattering systems, trace matches the high-accuracy IAS15 while providing a 15× speed-up. In large N systems simulating lunar accretion, trace qualitatively gives the same results as ias15 but at a 41× speed-up. We also discuss some cases such as von Zeipel–Lidov–Kozai cycles where hybrid integrators perform poorly and provide some guidance on which integrator to use for which system. trace is freely available within the rebound package. [ABSTRACT FROM AUTHOR]

We present new 'almost' time-reversible integrators for solution of planetary systems consisting of 'planets' and a dominant mass ('star'). The algorithms can be considered adaptive generalizations of the Wisdom–Holman method, in which all pairs of planets can be assigned time-steps. These time-steps, along with the global time-step, can be adapted time-reversibly, often at no appreciable additional compute cost, without sacrificing any of the long-term error benefits of the Wisdom–Holman method. The method can also be considered a simpler and more flexible version of the symba symplectic code. We perform tests on several challenging problems with close encounters and find the reversible algorithms are up to 2.6 times faster than a code based on symba. The codes presented here are available on Github. We also find adapting a global time-step reversibly and discretely must be done in block-synchronized manner or similar. [ABSTRACT FROM AUTHOR]

We present a simple algorithm to switch between N -body time integrators in a reversible way. We apply it to planetary systems undergoing arbitrarily close encounters and highly eccentric orbits, but the potential applications are broader. Upgrading an ordinary non-reversible switching integrator to a reversible one is straightforward and introduces no appreciable computational burden in our tests. Our method checks whether the integrator during the time-step violates a time-symmetric selection condition and redoes the step if necessary. In our experiments, a few per cent of steps would have violated the condition without our corrections. By eliminating them, the algorithm avoids long-term error accumulation, of several orders of magnitude in some cases. [ABSTRACT FROM AUTHOR]

Numerical integrations of the Solar system have been carried out for decades. Their results have been used, for example to determine whether the Solar system is chaotic, whether Mercury's orbit is stable, or to help discern Earth's climate history. We argue that all of the past studies we consider in this work are affected by numerical chaos to different degrees, affecting the possible orbits and instability probability of Mercury, sometimes significantly. We show how to eliminate the effects of numerical chaos by resolving Mercury's pericentre passage. We also show that several higher order symplectic maps do not exhibit significant differences in resolving pericentre passage of Mercury (at fixed time-step), making their advantages suspect for calculating long-term orbits. Resolving pericentre passage affects a wide array of orbital numerical studies, like exoplanet studies, studies of the galactic centre, and other N -body problems. [ABSTRACT FROM AUTHOR]

When fitting N -body models to astronomical data – such as transit times, radial velocity, and astrometric positions at observed times – the derivatives of the model outputs with respect to the initial conditions can help with model optimization and posterior sampling. Here, we describe a general purpose symplectic integrator for arbitrary orbital architectures, including those with close encounters, which we have recast to maintain numerical stability and precision for small step sizes. We compute the derivatives of the N -body coordinates and velocities as a function of time with respect to the initial conditions and masses by propagating the Jacobian along with the N -body integration. For the first time, we obtain the derivatives of the transit times with respect to the initial conditions and masses using the chain rule, which is quicker and more accurate than using finite differences or automatic differentiation. We implement this algorithm in an open source package, NbodyGradient.jl , written in the Julia language, which has been used in the optimization and error analysis of transit-timing variations in the TRAPPIST-1 system. We present tests of the accuracy and precision of the code, and show that it compares favourably in speed to other integrators that are written in C. [ABSTRACT FROM AUTHOR]

N -body integrations are used to model a wide range of astrophysical dynamics, but they suffer from errors which make their orbits diverge exponentially in time from the correct orbits. Over long time-scales, their reliability needs to be established. We address this reliability by running a three-body planetary system over about 200 e-folding times. Using nearby initial conditions, we can construct statistics of the long-term phase-space structure and compare to rough estimates of resonant widths of the system. We compared statistics for a wide range of numerical methods, including a Runge–Kutta method, Wisdom–Holman method, symplectic corrector methods, and a method by Laskar and Robutel. 'Improving' an integrator did not increase the phase-space accuracy, but simply increasing the number of initial conditions did. In fact, the statistics of a higher order symplectic corrector method were inconsistent with the other methods in one test. [ABSTRACT FROM AUTHOR]

Symplectic integrators are a foundation to the study of dynamical N -body phenomena, at scales ranging from planetary to cosmological. These integrators preserve the Poincaré invariants of Hamiltonian dynamics. The N -body Hamiltonian has another, perhaps overlooked, symmetry: it is smooth, or, in other words, it has infinite differentiability class order (DCO) for particle separations greater than 0. Popular symplectic integrators, such as hybrid methods or block adaptive stepping methods do not come from smooth Hamiltonians and it is perhaps unclear whether they should. We investigate the importance of this symmetry by considering hybrid integrators, whose DCO can be tuned easily. Hybrid methods are smooth, except at a finite number of phase space points. We study chaotic planetary orbits in a test considered by Wisdom. We find that increasing smoothness, at negligible extra computational cost in particular tests, improves the Jacobi constant error of the orbits by about 5 orders of magnitude in long-term simulations. The results from this work suggest that smoothness of the N -body Hamiltonian is a property worth preserving in simulations. [ABSTRACT FROM AUTHOR]

Hybrid symplectic integrators such as MERCURY  are widely used to simulate complex dynamical phenomena in planetary dynamics that could otherwise not be investigated. A hybrid integrator achieves high accuracy during close encounters by using a high-order integration scheme for the duration of the encounter while otherwise using a standard second-order Wisdom–Holman scheme, thereby optimizing both speed and accuracy. In this paper we reassess the criteria for choosing the switching function that determines which parts of the Hamiltonian are integrated with the high-order integrator. We show that the original motivation for choosing a polynomial switching function in MERCURY is not correct. We explain the nevertheless excellent performance of the MERCURY integrator and then explore a wide range of different switching functions including an infinitely differentiable function and a Heaviside function. We find that using a Heaviside function leads to a significantly simpler scheme compared to MERCURY  , while maintaining the same accuracy in short-term simulations. [ABSTRACT FROM AUTHOR]

Calculating the long-term solution of ordinary differential equations, such as those of the N-body problem, is central to understanding a wide range of dynamics in astrophysics, from galaxy formation to planetary chaos. Because generally no analytic solution exists to these equations, researchers rely on numerical methods that are prone to various errors. In an effort to mitigate these errors, powerful symplectic integrators have been employed. But symplectic integrators can be severely limited because they are not compatible with adaptive stepping and thus they have difficulty in accommodating changing time and length scales. A promising alternative is time-reversible integration, which can handle adaptive time-stepping, but the errors due to time-reversible integration in astrophysics are less understood. The goal of this work is to study analytically and numerically the errors caused by time-reversible integration, with and without adaptive stepping. We derive the modified differential equations of these integrators to perform the error analysis. As an example, we consider the trapezoidal rule, a reversible non-symplectic integrator, and show that it gives secular energy error increase for a pendulum problem and for a H'enon-Heiles orbit. We conclude that using reversible integration does not guarantee good energy conservation and that, when possible, use of symplectic integrators is favoured. We also show that time-symmetry and time-reversibility are properties that are distinct for an integrator. [ABSTRACT FROM AUTHOR]

The symplectic Wisdom–Holman map revolutionized long-term integrations of planetary systems. There is freedom in such methods of how to split the Hamiltonian and which coordinate system to employ, and several options have been proposed in the literature. These choices lead to different integration errors, which we study analytically and numerically. The Wisdom–Holman method in Jacobi coordinates and the method of Hernandez, H16, compare favourably and avoid problems of some of the other maps, such as incorrect centre-of-mass position or truncation errors even in the one-planet case. We use H16 to compute the evolution of Pluto's orbital elements over 500 million years in a new calculation. [ABSTRACT FROM AUTHOR]

We study analytically and experimentally certain symplectic and time-reversible N-body integrators which employ the Kepler solver for each pairwise interaction, including the method of Hernandez & Bertschinger. Owing to the Kepler solver, these methods treat close two-body interactions correctly, while close three-body encounters contribute to the truncation error at second order and above. The second-order errors can be corrected to obtain a fourth-order scheme with little computational overhead. We generalize this map to an integrator which employs the Kepler solver only for selected interactions and yet retains fourth-order accuracy without backward steps. In this case, however, two-body encounters not treated via the Kepler solver contribute to the truncation error. [ABSTRACT FROM AUTHOR]

We apply one of the exactly symplectic integrators, which we call HB15, of Hernandez & Bertschinger, along with theKepler problem solver ofWisdom&Hernandez, to solve planetary system N-body problems. We compare the method to Wisdom-Holman (WH) methods in the MERCURY software package, the MERCURY switching integrator, and others and find HB15 to be the most efficient method or tied for the most efficient method in many cases. Unlike WH, HB15 solved N-body problems exhibiting close encounters with small, acceptable error, although frequent encounters slowed the code. Switching maps like MERCURY change between two methods and are not exactly symplectic. We carry out careful tests on their properties and suggest that they must be used with caution. We then use different integrators to solve a three-body problem consisting of a binary planet orbiting a star. For all tested tolerances and time steps, MERCURY unbinds the binary after 0 to 25 years. However, in the solutions of HB15, a time-symmetric HERMITE code, and a symplectic Yoshida method, the binary remains bound for >1000 years. The methods' solutions are qualitatively different, despite small errors in the first integrals in most cases. Several checks suggest that the qualitative binary behaviour of HB15's solution is correct. The Bulirsch-Stoer and Radau methods in the MERCURY package also unbind the binary before a time of 50 years, suggesting that this dynamical error is due to a MERCURY bug. [ABSTRACT FROM AUTHOR]

We derive and present a fast and accurate solution of the initial value problem for Keplerian motion in universal variables that does not use the Stumpff series. We find that it performs better than methods based on the Stumpff series. [ABSTRACT FROM AUTHOR]

We present a new symplectic integrator designed for collisional gravitational N-body problems which makes use of Kepler solvers. The integrator is also reversible and conserves nine integrals of motion of the N-body problem to machine precision. The integrator is second order, but the order can easily be increased by the method of Yoshida. We use fixed time step in all tests studied in this paper to ensure preservation of symplecticity. We study small N collisional problems and perform comparisons with typically used integrators. In particular, we find comparable or better performance when compared to the fourth-order Hermite method and much better performance than adaptive time step symplectic integrators introduced previously. We find better performance compared to SAKURA, a non-symplectic, non-time-reversible integrator based on a different two-body decomposition of the N-body problem. The integrator is a promising tool in collisional gravitational dynamics. [ABSTRACT FROM AUTHOR]

Symplectic methods, in particular the Wisdom–Holman map, have revolutionized our ability to model the long-term, conservative dynamics of planetary systems. However, many astrophysically important effects are dissipative. The consequences of incorporating such forces into otherwise symplectic schemes are not always clear. We show that moving to a general framework of non-commutative operators (dissipative or not) clarifies many of these questions, and that several important properties of symplectic schemes carry over to the general case. In particular, we show that explicit splitting schemes generically exploit symmetries in the applied external forces, which often strongly suppress integration errors. Furthermore, we demonstrate that so-called 'symplectic correctors' (which reduce energy errors by orders of magnitude at fixed computational cost) apply equally well to weakly dissipative systems and can thus be more generally thought of as 'weak splitting correctors'. Finally, we show that previously advocated approaches of incorporating additional forces into symplectic methods work well for dissipative forces, but give qualitatively wrong answers for conservative but velocity-dependent forces like post-Newtonian corrections. We release REBOUNDx , an open-source C library for incorporating additional effects into REBOUND N -body integrations, together with a convenient python wrapper. All effects are machine independent and we provide a binary format that interfaces with the SimulationArchive class in REBOUND to enable the sharing and reproducibility of results. Users can add effects from a list of pre-implemented astrophysical forces, or contribute new ones. [ABSTRACT FROM AUTHOR]