Multi-Point Hermite Methods for the N-Body Problem

TL;DR

The paper tackles the challenge of accurate and efficient time integration for self‑gravitating $N$‑body systems. It introduces a family of multi‑point, multi‑derivative Hermite integrators that unify Adler–Bashforth–Moulton and Hermite methods, including a $3$‑point $6$th‑order, $9$th‑order, and $12$th‑order set with variable timesteps and a broader class of time‑symmetric collocation schemes up to $18$th order. In cluster‑like tests, the $3$‑point $6$th‑order scheme matches or exceeds the performance of the standard $2$‑point $4$th‑order Hermite at similar $O(N)$ cost, and higher‑order collocation methods show potential for long‑term accuracy and parallelization. Overall, the approach promises substantial efficiency gains for large‑$N$ collisional simulations and improved long‑term accuracy for few‑body integrations, especially when leveraging time symmetry and higher‑order derivative information.

Abstract

Numerical integration methods are central to the study of self-gravitating systems, particularly those comprised of many bodies or otherwise beyond the reach of analytical methods. Predictor-corrector schemes, both multi-step methods and those based on 2-point Hermite interpolation, have found great success in the simulation of star clusters and other collisional systems. Higher-order methods, such as those based on Gaussian quadratures and Richardson extrapolation, have also proven popular for high-accuracy integrations of few-body systems, particularly those that may undergo close encounters. This work presents a family of high-order schemes based on multi-point Hermite interpolation. When applied as a multi-step multi-derivative schemes, these can be seen as generalizing both Adams-Bashforth-Moulton methods and 2-point Hermite methods; I present results for the 6th-, 9th-, and 12th-order 3-point schemes applied in this manner using variable time steps. In a cluster-like test problem, the 3-point 6th-order predictor-corrector scheme matches or outperforms the standard 2-point 4th-order Hermite scheme at negligible O(N) cost. I also present a number of high-order time-symmetric schemes up to 18th order, which have the potential to improve the accuracy and efficiency of long-duration simulations.

Figures (9)

• Figure 1: A schematic diagram illustrating the points in time used by the 3-point variable-timestep Hermite schemes and the notation used in this section. To simplify notation, I define $\zeta\equiv\Delta t_0/\Delta t_1$.
• Figure 2: The convergence of the 3-point 6th-order variable-timestep scheme with respect to the dimensionless timestep control $\eta$ and the number of force evaluations during a 100-orbit eccentric 2-body test problem. Each scheme converges at the expected rate, or perhaps slightly faster when taking larger timesteps. The different timestep criteria result in minor differences in performance, though the generalized Aarseth criterion is slightly less efficient at lower accuracies.
• Figure 3: The convergence of the 3-point 9th-order variable-timestep scheme with respect to the dimensionless timestep control $\eta$ and the number of force evaluations in a 100-orbit eccentric 2-body test problem. In this case the PRS and Aarseth timestep criteria typically perform better than the generalized Aarseth criterion, resulting in fewer force evaluations to achieve a particular degree of accuracy.
• Figure 4: The energy errors incurred over the course of short-term (until $t=10$ in Hénon units) simulations of a 1024-object Plummer sphere using 2-point and 3-point Hermite integrators. The higher-order scheme is unequivocally more efficient, although the improvements it yields are more substantial at error tolerances below $\sim10^{-6}$ or $\sim10^{-4}$ depending on the timestep criterion employed.
• Figure 5: A schematic diagram illustrating the points in time used by the $n=2$, $r=1$ time-symmetric Hermite schemes and the notation used in this section. Equation (\ref{['eq:3x2mid']}) integrates from $t_0$ to $t_{1/2}$, while Equation (\ref{['eq:sym3x2']}) integrates from $t_0$ to $t_1$, though both use the same underlying interpolant.