Artificial precession and instability in solar system and planetary simulations: analytic and numerical results

Abstract

Wisdom-Holman (WH) methods are algorithms used as a basis for a wide range of codes used to solve problems in solar system and planetary dynamics. The problems range from the growth and migration of planets to the stability of the solar system. In many cases, these codes work with Democratic Heliocentric Coordinates (DHC) which offer some advantages. However, it has been noted these coordinates affect the dynamics of solar system bodies in simulations, in particular Mercury's, and introduce artificial precession which affects solar system stability. In this work, we analytically derive the two-body artificial precession induced by DHC. We show the effect is small for solar system bodies, but the artificial effect on Jupiter is $242$ times larger than on Mercury. In a two-body Mercury-Sun system with general relativity (GR), artificial precession is negligible compared to GR precession, even with extreme timesteps that amplify the numerical effects. In a full solar system model, numerical effects are amplified further.

Figures (4)

• Figure 1: Numerical verification of the dependence on eccentricity of the two-body precession rate, Eq. \ref{['eq:omfin']}. In units of solar mass, au, and yr, we let $m_0 = 1$, $m_1 = 0.001$, $a = 1$, $h = P/100$, $P = 1.00$, and $t_{\mathrm{max}} = 10^4 P$, where $t_{\mathrm{max}}$ is the runtime. The expected slope of $\log_{10} (1-e^2)$ versus $\log_{10} \left \langle \dot{\varpi}_{\mathrm{num}} \right \rangle$ is $s \in [-3.25, -3]$ as $e$ is increased to $1$. At large $e$, $e > 0.74$, our precession derivation breaks down. A least squares $s = -3.21$ for $e \le 0.74$ is plotted, and is consistent with theory.
• Figure 2: The relative error between the predicted ($\left \langle \dot{\varpi} \right \rangle_{\mathrm{an}}$) and numerical $(\left \langle \dot{\varpi} \right \rangle_{\mathrm{num}})$ precession rate as a function of $e$. The system and parameters are the same as those for Fig. \ref{['fig:escale']}. While Fig. \ref{['fig:escale']} tests the scaling with $e$, here we test the full Eq. \ref{['eq:omfin']}. The agreement is good for $e \le 0.74$, and agreement at higher $e$ is reestablished with smaller timesteps.
• Figure 3: Error in numerical precession as compared to Eq. \ref{['eq:omgr']}. A two-body problem consisting of Mercury and Sun, with $e = 0.6$, and a GR potential, is run. Even at large timesteps, artificial precession makes no significant impact.
• Figure 4: The artificial precession measured in a simplified two-planet Mercury--Jupiter system, as a function of timestep. Mercury's eccentricity has been increased to $e = 0.6$, and there is no GR. Strong prograde artificial precession for $h \ge 3.9$ days, surpassing the GR contribution (except at the last point), is found.