way forward toward separating out the two-body Keplerian Hamiltonians is to transform to heliocentric coordinates involving the centre-of-mass and the $\mathbf{r}_i-\mathbf{r}_0$. However, rewriting the Cartesian momenta in terms of heliocentric momenta (which have an additional component along the centre-of-mass momentum), leads to several cross-terms. Alternatively, Jacobi worked out a coordinate system in which the kinetic terms are particularly clean, and the kinetic energy remains a sum of squares. For readers that may not be familiar, and because our improved accuracy is largely due to modifications of the manner in which we transform between Cartesian and Jacobi coordinates, we briefly review them [see also @Plummer1918; @SussmanWisdom2001; @solarsystemdynamics].

Jacobi Coordinates
------------------

Rather than reference planet positions to the central star, a planet’s Jacobi coordinates are measured relative to the centre-of-mass of all bodies with lower indices. For concreteness, consider a system of $N$ particles with masses $m_i$, $i=0,\ldots, N-1$. Let $\mathbf{r}_i$ be the position vector of the $i$-th particle with respect to an arbitrary origin that is fixed in an inertial frame. Here we assume that the particles are ordered such that $i=0$ corresponds to the central object, $i=1$ to the innermost object orbiting the central object and so on. The existence 