February 16, 2026
the purely secular phase-space portraits computed using Hamiltonian (9). As an example of this correspondence, Fig- ure 13 illustrates a comparison between purely secular co-orbital dynamics (left panel)and the e ∆a secular dynamics embedded within the 1:1 mean motion resonance with Planet Nine15. The resemblance between purely secular and resonant-secular phase- space portraits in Figure 13 is self-evident. This inherent separation of the degrees of freedom further qualitatively explains how KBOs can smoothly enter and exit MMRs facilitated by Planet Nine without noticeably affecting the orbital clustering of the long- period objects (see Bailey et al. 2018 for further discussion).
A complete formulation of the Planet Nine hypothesis requires both a well-developed analytical understanding of the underlying physics and a detailed numerical description of the associated dynamics. The previous section provides an analytic account of the effects of Planet Nine on the orbital structure and evolution of the distant Kuiper belt. However, a comprehensive comparison between the observations and P9-sculpted or- bital structure requires a more detailed description of long-period KBO evolution. As a result, the Planet Nine hypothesis must be explored with the aid of realistic numerical simulations. The description, execution, and analysis of such numerical experimenta- tion is the primary focus of this section.
Although the Planet Nine hypothesis (in its current incarnation) was put forth only a few years ago, a vast collection of N-body simulations has already been published. This body of work employs varying degrees of approximation in order to understand the dynamical mechanisms through which Planet Nine sculpts the distant Kuiper Belt. For example, one class of calculations (Beust, 2016; Batygin and Morbidelli, 2017; Hadden et al., 2018; Bailey et al., 2018) employs a simplified model of the solar sys- tem, where the Keplerian motion of all four canonical giant planets is smoothed over. In this approximation, the solar system interior to 30 AU is replaced with a central body possessing a quadrupole field equivalent to the phase-averaged gravitational field of Jupiter, Saturn, Uranus, and Neptune. More realistic calculations (Batygin and Brown, 2016a; Brown and Batygin, 2016; Becker et al., 2017; Millholland and Laughlin, 2017; Khain et al., 2018a) include Neptune and Planet Nine directly as active bodies, while replacing the orbits of the three interior giants with massive wires; this scheme provides a compromise between capturing short-period effects in the Kuiper belt and retaining a long simulation time-step. Finally, numerical experiments that fully resolve the dy- namics of the entire outer solar system have also been performed (Batygin and Brown, 2016b; Sheppard and Trujillo, 2016; Gomes et al., 2016; Bailey et al., 2016; Becker et al., 2018; Sheppard et al., 2018).
The successive approximations to the dynamics described above represent a trade- off between increasingly realistic orbital evolution and increased computational re- sources. The model that averages out the Keplerian motion of the canonical giant planets is clearly the fastest, but the least realistic. The intermediate simulations that
Strictly speaking, the secular phase-space portrait shown on Figure 13 should be computed by confining the value of дres to its equilibrium value, which is generally a function of both e and ∆a. In the special case of the 1:1 MMR, however, secular modulation of øres is very small, and simply assuming that øres = 0 everywhere leads to an insignificant error in the computation of the resonant-secular Hamiltonian.
directly account for the full dynamics of Neptune and Planet Nine, but smooth over the other planets, are expected to yield more precise results. This approach retains a diminished computational cost relative to the final approximation of keeping all major bodies active. To this end, note that the CPU run-time for a given integration is set by the resolution of the shortest dynamical timescale and the orbital period of Neptune is longer than that of Jupiter by a factor of 14. Although the approach of including only Neptune and Planet Nine as active perturbers is often used, the extent to which this approximation captures the relevant physics remains to be determined quantitatively.
This section presents an overview of the numerical simulations of the Planet Nine hypothesis carried out to date. In order to provide a consistent description of the nu- merical results, and to sort through the computational approximations used in previous work, we have carried out an additional ensemble of simulations of the solar system with the inclusion of Planet Nine. In the first set of calculations, we adopt the interme- diate description of the dynamics, where only Planet Nine and Neptune are considered as active bodies (section 5.1). We then carry out an analogous set of simulations that self-consistently account for the N-body dynamics of all outer planets (section 5.2). With these additions, an analysis of all of the numerical simulations is presented in sec- tion 5.3, which assesses which properties of Planet Nine (mass and orbital elements) provide the best description of the observed solar system anomalies. The particular case of TNOs on high inclination orbits is taken up in section 5.4. Finally, we discuss how interactions between Planet Nine and the known planets can affect the orientation of the orbital plane of the solar system (section 5.5). A complementary discussion that concerns the dynamical stability of the observed KBOs in the presence of Planet Nine is presented in Appendix C.
As an initial numerical test of the P9 hypothesis, we simulate the effects of includ- ing Planet Nine and Neptune on a population of TNOs. The primary goal is to see which properties of Planet Nine control the onset of observed orbital anomalies of the TNOs. Following Batygin and Brown (2016a), in our preliminary suite of simulations we set the absorbing radius of the central body to = 29 AU, and endow the sun with a J2 gravitational moment given by (Burns, 1976)
7 m a2
J = 1 X j j , (14)
where where the sum runs over contributions from Jupiter, Saturn, and Uranus. In the Figures and the following text, we refer to this simulation setup with the abbreviation J2NP9.
Adopting a time-step of ∆t =14 years ( 1/10) th of Neptune’s orbital period), we evolve an initially unstructured population of test-particles for 4 Gyr, allowing it to be sculpted by the combined gravitational influence of Planet Nine and Neptune as well as the mean-field effects of the solar system interior to . Following Khain et al. (2018a),
we initialize the Kuiper belt as a disk of N = 1000 eccentric planetesimals, spanning the semi-major axis range a e (100, 800) AU, perihelion range of q e (30, 100) AU,
and adopt a half-gaussian inclination dispersion with standard deviation of σi = 15 degrees, while drawing the other orbital angles a, Ω, from a uniform distribution. With respect to the initial q-distribution in particular, we note that although secular dy- namics driven by P9 alone can act to lift the perihelia of Neptune-attached KBOs (see section 4.2), this process is also expected to naturally manifest within the solar sys- tem’s birth cluster (Morbidelli and Levison, 2004), justifying our adoption of an broad primordial spread of perihelion distances. In light of phase-space mixing associated with the chaotic nature of P9-facilitated evolution, however, we do not expect specific choices of initial conditions to affect the final results significantly. All numerical sim- ulations reported herein were performed using the mercury6 gravitational dynamics software package (Chambers, 1999), employing the hybrid Wisdom-Holman/Bulirsch- Stoer algorithm (Wisdom and Holman, 1992; Press et al., 1992).
In total, 1134 semi-averaged N-body integrations were performed, with the simula- tion grid corresponding to Planet Nine parameters a9 e (300, 1500) AU, e9 e (0.05, 0.95), i9 e (10, 35) deg, with increments ∆a9 = 100 AU, ∆e9 = 0.1, ∆i9 = 5 deg, and masses m9 = 5, 10, and 20 M⊕. Drawing on previously published results (Brown and Batygin 2016; Hadden et al. 2018; Ca´ceres & Gomes 2018 and the references therein), sim- ulations in which Planet Nine’s perihelion distance was smaller than q9 ≤ 100 AU or greater than q9 ≥ 500 AU were immediately discarded. Moreover, Planet Nine’s starting orbital angles were chosen such that its argument of perihelion at the end of
the simulation was always close to a9 Ω9 = v9 140 deg i.e., approximately 180 degrees away from the mean argument of perihelion of the dynamically stable KBOs with a ≥ 250 AU.
Orbital Footprints of Stable KBOs. Ideally, the number of particles in any given simu- lation would be so large that the synthetic Kuiper belt created after 4 Gyr of integration could be compared in a quantitatively sound manner with the observational dataset. Unfortunately, it is not computationally feasible to carry out such detailed simulations in sufficient numbers to meaningfully explore the P9 parameter space. As a result, rather than analyzing the final orbits of simulated long-term stable16 KBOs, we will instead examine the chaotic orbital footprints generated by these particles. In this con- text, the orbital footprints are generated simply by plotting the orbital state of the test particles every 1 Myr throughout the latter half of the integration. The resulting scatter of points delineates the overall orbital characteristics of surviving TNOs, including the range in angle ∆a as a function of semi-major axis a. As outlined in section 3, the observed TNOs show a distinct pattern for sufficiently large a that must be consistent with successful Plant Nine models. We further note that owing to their fundamentally stochastic nature, these footprints trace out the confines of the available phase-space, and thus define a P9-sculpted distribution of orbital parameters which can then be com- pared against the statistical properties of the observed census of distant KBOs.
The orbital footprints of ∆a = a – a9, generated by stable particles over the 2 – 4 Gyr integration timespan are shown as a function of a in Figure 14. Simulation16For the purposes of this work, we define long-term stability as being equivalent to a dynamical lifetime in excess of 4 Gyr.
30 AU 6 q 6 100 AU i < 40 deg 100 AU < q i < 40 deg 30 AU 6 q i > 40 360
270
180
90
m9 = 5M a9 = 500 AU e9 = 0.25
i9 = 20 deg
270
180
90
360
270
180
90
m9 = 5M a9 = 500 AU e9 = 0.25
i9 = 20 deg
0
200
400
600
800
1000
0
200
400
600
800 1000
m9 = 10M a9 = 800 AU e9 = 0.45
i9 = 15 deg
270
180
90
360
270
180
90
m9 = 10M a9 = 800 AU e9 = 0.45
i9 = 15 deg
0
200
400
600
800
1000
0
200
400
600
800 1000
360
270
180
90
m9 = 20M a9 = 1100 AU e9 = 0.65
i9 = 10 deg
0
200
400
600
a (AU)
800
0
1000
200
400
600
a (AU)
800 1000
Figure 14: Apsidal confinement of simulated KBOs. Each panel depicts the orbital footprints of longitude of perihelion of stable KBOs relative to Planet Nine, ∆a, as a function of semi-major axis, a. The top, middle and bottom pairs of plots depict synthetic Kuiper belts, generated with Planet Nine masses of m9 = 5 M⊕ (top
panels), 10 M⊕ (middle panels), and 20 M⊕ (bottom panels). The corresponding values of P9 orbital elements
are labeled above each row. The panels on the left-hand-side of the figure show the results of JSUNP9
simulations, where the orbital motion of Jupiter, Saturn, and Uranus is resolved self-consistently. The right- hand-side panels show analogous results, obtained within J2NP9 simulations that only treat Neptune and Planet Nine as active perturbers and replace the other planets with an effective quadrupole moment of the sun. The enhanced fuzziness of orbital footprints shown on the left panels stems primarily from our use of heliocentric rather than barycentric orbital elements. Blue points signify orbits that have 30 AU ≤ q ≤ 100 AU and i < 40 deg, and thus satisfy our criteria for observability. Pink points mark unobservable orbits with q > 100 AU and i < 40 deg. Orange dots denote high-inclination objects with i ≥ 40 deg. The large circular points depict the orbital elements of the observed KBOs, color-coded in the same way as in Figures 6 and 8. In all cases, the value of ∆a for the simulated TNOs is randomly distributed for small semi-major axes, and becomes confined to a much more limited range for larger a. The critical (transition) value of the
semi-major axis is marked on each panel with a vertical dashed line. The top of each panel further reports the apsidal confinement fraction, fa, which is given by the ratio of blue points that fall within 90 deg of exact perihelion anti-alignment with P9 (the domain between the two horizontal dashed lines) to the total number of blue points that reside outside of the critical semi-major axis.
results corresponding to three sets of P9 parameters are depicted, with J2NP9 calcu- lations illustrated on the right-hand-side. While the panels of this Figure contain the orbital footprints generated by all long-lived particles over the last two billion years of their dynamical evolution, only a fraction of these synthetic data points would have been detectable by past and ongoing surveys. Thus, we require a proxy to distinguish between observationally detectable and undetectable objects. One possibility would be to develop a survey simulator to construct samples of observable data points from the simulations. Such a simulator, however, would require adequate characterization of the surveys which could have detected these distant objects. The number of objects
detected with a & 250 AU is relatively small and they are found by several different surveys, each with differing detection efficiencies and biases. Unfortunately, a defini- tive assessment of the cumulative biases and other survey characteristics has not yet been carried out, although simulators tailored to specific surveys exist, and additional
simulators are under construction (see, e.g., Shankman et al. 2017; Hamilton et al., in prep). Here we instead make simple orbital element cuts to approximate observational biases.
Well known biases exist against the detection of both high perihelion and high in- clination objects. Maintaining consistency with our discussion of the observations, we approximate these biases by classifying the simulation data into three categories, as follows. Objects deemed discoverable by conventional ecliptic surveys are character-
ized by 30 AU ≤ q < 100 AU and i < 40 deg, and are shown with blue points in Figure
14. The undiscoverable counterpart to this population has q ≥ 100 AU and i < 40 deg, and is shown with pink points. Finally, the high-inclination component of the synthetic Kuiper belt is taken to be composed of objects with q ≥ 30 AU and i ≥ 40 deg, and is graphically depicted with orange points. Akin to Figure 8, the aggregate of observed a ≥ 250 AU KBOs is over-plotted on the panels of Figure 14, maintaining the same color-scheme as that employed in Figure 6.
The Critical Semi-Major Axis. In each simulation depicted in Figure 14, the observa- tionally discoverable subset of orbital footprints (blue points) exhibits a well-defined boundary, where low-a objects experience apsidal circulation and large-a particles clus- ter around ∆a 180 deg. Accordingly, we can characterize the suite of J2NP9 cal- culations by delineating the critical semi-major axis beyond which orbital clustering emerges, ac, as a function of a9, e9, i9, and m9.
In order to map the dependence of ac on P9 parameters, we examined the (obser- vationally discoverable) distribution of ∆a as a function of a in each simulation. We then computed the apsidal confinement fraction, ζa, defined as the fraction of parti- cle footprints that lie within 90 deg of exact anti-alignment with Planet Nine (shown with dashed horizontal lines in Figure 14) inside semi-major axis bins δa = 10 AU. As semi-major axis is increased from 150 AU outward, ζa shifts from an initial value very close to 1/2 (implying a fully randomized a distribution), to some higher value (typi- cally close to unity), signifying a clustered distribution. To parameterize this transition in the apsidal confinement fraction, we fit an an error function to the resulting (a, ζa) sequence, and interpret the crossover point in the fit as ac. In each panel depicted in Figure 14, this transition point is marked with a vertical dashed line, and marked on top of each plot.
No comments:
Post a Comment