Details of Resonant Structures within a Nice Model Kuiper Belt: Predictions for High-perihelion TNO Detections

The eccentricity distribution of the resonances show "fingers" of descending eccentricity from the main population of outer classicals and detached TNOs (Figure 1). These populations are likely created as a result of the Kozai mechanism (Gomes et al. 2005), which also contributes to the hotter inclination distributions. This effect is stronger for the n:1 resonances, but is also apparent for the n:2 resonances. (See Section 4 for a discussion of the extensive effects of Kozai.) Simply due to large distances, TNOs in these distant resonances are sparsely sampled by current observations, particularly at low e (e.g., Gladman et al. 2012), and theoretical models of their distributions would provide useful constraints for population models derived from surveys.

Figures 2 and 3 show cumulative distributions in eccentricity and inclination for the n:1 and n:2 resonances, respectively. The resonances exterior to the 2:1 resonance do not contain very low e members, and the minimum e roughly increases with semimajor axis. The maximum e also increases with semimajor axis; this is an expected result because at larger a, a larger e is necessary for the particle to scatter off Neptune. The inclination distribution of objects in the distant resonances is visibly broader and hotter than for the 3:2 and 2:1. This different morphology between the resonances likely results from the additional Neptune encounters needed to emplace an object at large a.

Zoom In Zoom Out Reset image size

The range of eccentricities and inclinations displayed by the distant n:1 and n:2 resonances are particularly interesting because they affect population estimates and are weakly constrained by surveys (e.g., Volk et al. 2016). Survey population estimates suggest the distant resonances should have relatively large populations, albeit with large uncertainties (Gladman et al. 2012; Pike et al. 2015; Alexandersen et al. 2016). Because the highest eccentricity members of a resonance will have the closest pericenter distances and thus the brightest apparent fluxes, the highest eccentricity members of any resonance are the objects most likely to be detected in a magnitude-limited survey. In the more distant resonances, TNOs on low-eccentricity orbits are not ever bright enough to observe in OSSOS and other current surveys, so the total population of these distant resonances is impossible to constrain unless an eccentricity distribution is assumed. The B&M model presented here suggests a possible eccentricity distribution that is tied to a dynamical history for the solar system.

Several conclusions are evident from the distributions of inclination and eccentricities for the n:1 and n:2 populations in the B&M model. The B&M model suggests a hot inclination distribution is likely for all n:1 resonances beyond the 3:1, but increased inclinations are not obvious for the n:2 resonators until the 13:2, at a much larger semimajor axis. (The 2:1 resonance has an oddly biased inclination distribution, with an exaggerated high-inclination component due to cloning: see further discussion in Section 3.3.) The average eccentricities of n:1 and n:2 resonant populations increase with semimajor axis; the more distant n:1 resonances all have less than 5% of particles with $e\lt 0.35$ (see Figure 2), and the distant n:2 resonances have less than 5% of particles with $e\lt 0.5$ (see Figure 3). The bulk of the population in these resonances has nearly the maximum stable eccentricity. If this distribution represents reality, population estimates that exclude a low-e component are not significantly underestimating the population (e.g., Volk et al. 2016), while a uniform e distribution would overestimate the population by including far too many TNOs at low e.

Zoom In Zoom Out Reset image size

The n:1 resonances tend to have higher inclinations and eccentricities at a given semimajor axis than nearby n:2 resonances. We suspect this may be because capture into the n:1 resonances requires more scattering events than that for the n:2, although previous simulations (Lykawka & Mukai 2007b) find that capture probability is highest for the lowest-order resonances and also for resonances at larger a. Whether this effect is due to capture probabilities or to the timescale of diffusion to higher e and i within the resonances, and whether or not Kozai oscillations play a significant role, will require further modeling.

The libration amplitude, or range of resonant angles a particle explores (the amplitude of Equation (1) in P17), is influenced by the specific capture mechanisms and migration timescales (e.g., Chiang & Jordan 2002). The libration amplitudes of the stable B&M n:1 and n:2 resonators are shown in Figure 4. The n:2 resonators occupy a single resonant island around 180° and have libration amplitudes from 6° (for the lowest libration amplitude 5:2 resonator) to just less than 180° (for the largest amplitude 5:2 resonator). The n:1 test particles include both large-amplitude symmetric and small-amplitude asymmetric librators.

Zoom In Zoom Out Reset image size

The B&M n:2 resonators have a large range of libration amplitudes in a single libration island. The 3:2 and 5:2 resonances have roughly uniform distributions from 20° to 140° and 160°. The 7:2 and 9:2 resonances have a lower median libration amplitude, and the 11:2 and 13:2 resonators have a higher mean libration amplitude than the 3:2 and 5:2 resonators. In all of these resonances, libration amplitudes below 20° are either nonexistent or an insignificant fraction of the population. This reflects a combination of reduced phase space at low amplitude and the inefficiency of deep capture because these populations are implanted through scattering capture and not through sweep-up. Volk et al. (2016) searched for low libration amplitude Plutinos (3:2 resonators) in the OSSOS survey. Their minimum libration amplitude detected TNO has an amplitude of 10°${}_{-4}^{+8}$, which may require resonance sweeping capture. This low libration amplitude Plutino discovery suggests that the Plutinos include a sweep-up capture component; the presence or absence of low libration amplitude 5:2 resonators in the real Kuiper Belt could provide useful insight into the extent of the primordial planetesimal disk. If a significant low libration amplitude component exists in the 5:2, the primordial disk would have had to extend to large enough distances that sweep-up could operate effectively in the 5:2, suggesting a much wider primordial disk than often theorized. A small, swept-up, low libration amplitude component resulting from a lower-density extended disk (not included in B&M) would be consistent with the Neptune's migration evolution in this simulation, but quantifying this component would require additional simulations of the instability.

The dynamical behavior of the n:1 resonators is more complex than that of the n:2 test particles because the n:1 resonators can occupy three libration islands: the symmetric island centered at 180° as well as the leading (centered at small angles <180°) and trailing (centered at large angles >180°) asymmetric libration islands. Because of the segregated phase space, symmetric librators all have large libration amplitudes ($\gtrsim 120^\circ$), and asymmetric librators have smaller amplitudes. As a result, the relative fractions of symmetric and asymmetric resonators are apparent in Figure 4, as a large gap separates the libration amplitude of the symmetric and asymmetric islands.

The fraction of TNOs in the different islands, particularly in the 2:1, has been theorized to be a diagnostic of the timescale of Neptune's migration (Section 5.1; Chiang & Jordan 2002). Due to the large number of test particles in our simulation, we are able to measure this fraction for six of the n:1 resonances (Table 1). In all of these resonances, the leading asymmetric islands contain more particles than the trailing islands. To assess the significance of these values, we compute the probability of randomly finding at least the measured number of leading librators when an equal number of particles are expected in the leading and trailing islands ("Likelihood" column in Table 1). Using this, we show that most of the n:1 resonances have a leading island occupation that is much too high to be explained by random fluctuations. The leading islands are more populated in each of the n:1 resonances due to emplacement during this simulation.

3.3.1. Leading Island Enhancement Is Not an Artifact of Integration Methods

3.3.2. Leading Island Enhancement Is Not Caused by Neptune's Orbital Circularization

We explored the effect of the circularization of Neptune's orbit on the populations in the leading and trailing 2:1 resonant islands with an additional simulation. This simulation used an initial Neptune eccentricity of 0.3, damped to $e\sim 0.01$ over 10 Myr. The other giant planets were included with their current orbital configuration. A total of 26,450 test particles were given initial conditions covering the phase space of the 2:1 resonance. After Neptune's orbit was circularized, the end state of the simulation was integrated without migration for an additional 10 Myr for classification, as we did for the B&M simulation. Thus 901 particles were identified as librating in the 2:1 resonance, with 218 in the leading and 275 in the trailing libration island, a 26% enhancement of the trailing island (Figure 5), the opposite enhancement from that seen in the B&M simulation. If equal populations are expected in both islands, this enhancement has a statistical likelihood of only 0.5%. This migration simulation involving only Neptune circularization results in a statistically significant enhancement of the trailing libration island.

Zoom In Zoom Out Reset image size

The test particle initial conditions were selected to initially populate the 2:1 resonance and lose resonators only as a result of Neptune's circularization, with no outward migration. This evolution did not result in the enhancement of the leading island seen in the B&M simulation; as shown in Figure 5, it resulted in the opposite enhancement. The full emplacement scenario is clearly important for creating the B&M asymmetry; it cannot be attributed to any isolated portion of the migration and emplacement simulation. The particle initial conditions may also play a role in their resonant behavior, or perhaps the scattering process in B&M combined with the large e of Neptune and the particles is necessary to enhance the leading island. The circularization of the large-e Neptune clearly can play a role in creating asymmetries in resonances, but that alone does not reproduce the enhancement in the leading island of the 2:1 resonance. Further simulations must be undertaken to understand the enhancement effect we observe here.

3.3.3. Leading Island Enhancement Is Not an Artifact of Cloning

Kozai resonators are identified within several MMRs in the B&M simulation. Figure 6 shows the eccentricity–inclination distribution for three mean-motion resonances with Kozai oscillators: the 3:2, 5:2, and 2:1 resonances. Kozai oscillators in each of the resonances occupy a significant fraction of the highest e and i space.

Zoom In Zoom Out Reset image size

Over the course of a Kozai cycle, as e and i cycle to higher and lower values, the quantity

$$\begin{eqnarray}&&{L}_{z}\propto \cos i\sqrt{1-{e}^{2}}\end{eqnarray} \tag{ 1 }$$

is conserved. Lines of constant L_z are shown in Figure 6. Some particles in the simulation that were cloned while in a Kozai cycle are distributed along curves of constant L_z. This is particularly obvious in the 5:2 and 2:1 distributions (see cloning discussion in Section 3.3).

Particularly for the Kozai Plutinos (3:2 resonators), larger i correlates with a constrained range of e values (see Figure 6). This dependence is not usually included in population models, which is a major weakness of typical resonant population modeling. We find that 21% of the 1640 B&M Plutinos are in Kozai, which agrees well with the highest-likelihood Kozai fraction of 20% for the OSSOS survey found by Volk et al. (2016). Though the Kozai fraction matches observations, the inclination distribution within the B&M Plutinos does not include particles sufficiently excited to match known Plutinos. A slower migration rate for which Kozai is more effective may be key to a sufficiently hot inclination distribution for the Plutinos (as in Kaib & Sheppard 2016).

The inclination distribution of the B&M model 5:2 resonators appears to have multiple discrete components as a result of Kozai resonance and particle cloning during the simulation. The $e\lt 0.3$ particles all have large i; two components for both the eccentricity and inclination are apparent. We find that 29% of the 337 5:2 resonators are in Kozai. The inclination distribution of the non-Kozai fraction of the 5:2 is less dynamically hot than the 5:2 population as a whole; this distinction has important implications for population models.

Even though Kozai is not operating for all of the resonant particles at the end state of the simulation, its effect is clearly seen in the distribution. The interdependence of inclination and eccentricity as well as the presence of high-inclination particles is indicative of Kozai evolution. As discussed below, Kozai resonance dropout may be just as important for shaping the orbital element distribution within resonances as resonant dropout is for shaping the distribution of detached objects (Kaib & Sheppard 2016; Nesvorný et al. 2016).

The effects of temporary evolution in Kozai resonance are apparent in the regions beyond the main classical belt. In the B&M simulation, all particles in this region have been emplaced as a result of the planetary migration; the initial planetesimal disk was truncated at 34 au. The classification of these particles is based on Gladman et al. (2008). As shown in Figure 1, this region is populated by outer classical objects (at $e\lt 0.24$, assumed to be an extension of the main classical belt in this classification scheme), objects in resonance with Neptune, detached objects, and scattering objects, which have semimajor axes that are changing on short timescales as a result of Neptune's influence. In Figure 1, the outer classical objects and some detached objects can be seen located near the resonance fingers. The q = 40 au boundary is shown as a black arc starting at a = 40 au. These outer classical particles all have $q\gt 40\,\mathrm{au}$ and were likely emplaced by the same mechanism that created the detached objects in the same q range; these particles cannot be primordial. The high-a ($a\gt 62.5$, outside the 3:1 MMR) particles with $q\gt 40\,\mathrm{au}$, which we refer to as "high-q" objects, seem to be a more informative dynamical subdivision and include test particles from the detached and outer classical populations. These subpopulations have different dynamical characteristics as a result of their emplacement histories. Figure 7 shows the inclination and eccentricity distributions of particles in the various subpopulations in the distant Kuiper Belt.

Zoom In Zoom Out Reset image size

The Gladman et al. (2008) dynamical classification system beyond the main classical belt is not optimal for understanding the B&M simulation results. This classification separates test particles that share a common origin into different subcomponents, and it also places into the same classification some particles with very different evolutionary histories. The "outer classical" objects in the B&M simulation are clearly just the low-e extension of the high-q population, which is much more dynamically distinct. With the initial conditions of the B&M model, objects in the outer classical region cannot be primordial. There is a shared dynamical origin for the outer classical and some detached objects, so a pericenter cut for dynamical classification is more informative than the Gladman et al. (2008) classification in the distant solar system.

Figure 7 shows that the high-q subpopulation has different characteristics than the detached and scattering populations. Because of their selection criteria, these high-q particles have a much colder eccentricity distribution, reaching down to e = 0.16. However, this population also has a hotter inclination distribution, with a median of 28°. Based on the anticorrelated e and i distributions, it appears that these particles were emplaced through resonance dropout or diffusion, after the Kozai mechanism reduced their eccentricities and increased their inclinations. In the e–i plot in Figure 7, the lines of constant L_z resulting from Kozai resonance are apparent. When objects evolve to low e, where the resonance is narrower, particles are more likely to exit the resonance. Some of the detached objects are not particularly close to a resonance, and these particles (as well as others near resonance) may have been emplaced through dropout as Neptune (and the resonances) migrated outward (Gomes 2003). A circularization of Neptune's orbit can also shrink the width of the resonance, particularly at low eccentricity, and cause resonance dropout (Gomes et al. 2008). If particles drop out during migration, they should be primarily located Sunward of resonances, while dropout during circularization would have a more symmetric signature around resonances (see Section 5.2). These high-q particles are consistent with formation through resonant dropout after modification by Kozai within distant mean-motion resonances.

The detached population includes the majority of these high-q objects in addition to a large $q\lt 40\,\mathrm{au}$ population. The eccentricity distribution of the detached objects is somewhat hotter and the inclination distribution is somewhat colder than the $q\gt 40\,\mathrm{au}$ subpopulation, indicating the larger sample includes particles less modified as a result of Kozai. In Figure 7, only the $q\lt 40\,\mathrm{au}$ subcomponent of the detached objects is presented. These particles still may be emplaced through resonant dropout, but the Kozai mechanism likely did not play a significant role in the $q\lt 40\,\mathrm{au}$ portion of the detached population.

The scattering objects, in contrast, have hotter eccentricity distributions and colder inclination distributions. These particles were also likely modified by the Kozai mechanism, increasing the eccentricity and decreasing the inclination during a resonance capture. The mean-motion resonances are weaker at the extreme values of e and i, so escape is more likely for these particles. The lines of constant L_z are not apparent in this population in Figure 7, because the scattering process erases this signature. After the particle diffuses or drops out of resonance, it is scattered by Neptune on relatively short timescales.

The distant particles in the B&M simulation have a range of eccentricity and inclination distributions indicative of their emplacement history. The Kozai mechanism within MMRs has resulted in populations anticorrelated in i and e, which is evident in the plot of all objects in i versus q in Figure 7. This distribution of objects in the B&M simulation in e and i is consistent with the observed Kuiper Belt, and the effect of the Kozai mechanism has been noted in previous simulations, which found a noticeable fraction of distant resonant objects evolving to low e (Gomes 2003; Pike et al. 2015; Kaib & Sheppard 2016; Nesvorný et al. 2016).

The relative fraction of leading and trailing asymmetric librators for the 2:1 resonance may be diagnostic of Neptune's mode of migration. In a different Neptune migration scenario with a low-eccentricity Neptune, Chiang & Jordan (2002) tested a range of simulation migration speeds, from 50,000 yr to 5 Myr, and found that the faster migration speeds increased the fraction of 2:1 resonators in the trailing libration island. Both Chiang & Jordan (2002) and Murray-Clay & Chiang (2005) found that the trailing 2:1 libration island should be more populated than the leading island, an effect that became more amplified with faster Neptune migration (up to ∼300%). However, the B&M model has a factor of two enhancement in the leading island over the trailing island. The preference for the leading libration island continues through all n:1 resonances that have enough test particles for a reliable ratio to be determined (though not all are statistically significant; see Table 1). The B&M simulation results confirm that migrations can produce asymmetries in libration island capture, but migration rate is not the only factor that can create an asymmetry. The enhancement in particular libration islands must depend on additional migration parameters than purely migration speed, such as the starting orbital parameters of all of the giant planets and the initial orbits of the particles.

Measuring the libration island population fractions is a goal of characterized surveys; the complete OSSOS survey will provide strong constraints on this ratio (e.g., Volk et al. 2016). Based on the first results from the OSSOS survey, Volk et al. (2016) provides the only statistically robust constraint on the fraction of 2:1 asymmetric resonators in the leading island at <90%; refining this estimate is a goal of the full OSSOS survey. Determining relative population sizes of the different resonant components in the Kuiper Belt based on survey results will provide constraints on different solar system migration models.

More theoretical work is needed to understand the additional constraints that influence the specifics of where the asymmetries appear. This is unfortunate as many authors have hoped that measuring the ratios of populations in the various libration islands could provide a straightforward diagnostic of migration speed; this result shows the interpretation will be more complicated. Recent "grainy" migration simulations (Kaib & Sheppard 2016; Nesvorný et al. 2016) have shown promise in reproducing aspects of the observed Kuiper Belt, but these works have not presented the resulting n:1 resonant island fractions, so no direct comparisons can yet be made.

Recent simulations by Nesvorný et al. (2016) and Kaib & Sheppard (2016) (hereafter referred to as N16 and KS16, respectively) present detailed analyses of resonance population by smooth and grainy Neptune migration models. Our Figure 1 can be directly compared with Figure 2 in KS16 and Figures 1 and 2 in N16. These figures highlight the difference in the orbital distribution of near-resonant distant TNOs.

The main difference in near-resonant structure of the models is that the B&M migration model produces a significant number of resonant dropouts at low e on both sides of each resonance. Due to its relative isolation, this resonant "beard" is most easily visible by eye around the 3:1 resonance (Figure 1). This relatively symmetric distribution of low-e test particles around each resonance happens because Neptune's very high eccentricity in the B&M simulation initially produces wide, powerful resonances. As Neptune's orbit circularizes and the resonances narrow, particles drop out of resonance and produce this beard-like structure (as discussed in Gladman et al. 2012). In the smooth and grainy migration simulations in KS16 and N16, where Neptune stays below a more moderate e of ∼0.1, the resulting features around resonances are not symmetric beard-like structures; there are more particles dropped on the sunward side of each resonance.

Figure 8 highlights this near-resonant structure in the high-q ($q\gt 40$ au) population. (Compare this directly with Figures 3 and 6 in KS16.) While all of the simulations presented in N16 and KS16 have many more particles located just Sunward of the resonances (due to resonant dropout during Neptune's migration phase), the B&M simulation produces many particles on both sides of the resonances (Figure 8). This is likely due to a much more significant amount of dropout occurring during the circularization phase of Neptune's orbit evolution during the B&M simulation, which attains much higher eccentricity during the Nice model instability ($e\sim 0.3$) as compared with the N16 and KS16 simulations ($e\sim 0.1$).

Zoom In Zoom Out Reset image size

As discussed extensively in KS16, when exactly the low-e dropouts primarily happen depends on the timescale of Neptune's orbital evolution relative to the Kozai cycling timescale within MMRs. In the N16 and KS16 analyses, more dropouts occur Sunward of the resonances in the slower migration simulations. This is because the timescale for particle capture and Kozai cycling is shorter than the slow migration timescale (∼100 Myr) but longer than the fast migration timescale (∼10 Myr).

Particles are captured into Neptune's migrating resonances and then pushed to higher e by Neptune's continued outward migration, into the phase space where Kozai oscillations are stronger. Kozai cycling will eventually drive the particles to higher i and lower e, where the resonance is weaker and dropout is more likely. How Neptune's orbit evolves on the timescale of a Kozai cycle within an MMR (several megayears) will determine where low-e dropouts are deposited. In the case of faster migration (10 Myr timescale simulations in N16 and KS16), Neptune's orbit migrates outward significantly before Kozai cycling can dramatically modify particle orbits, so the particles are likely to remain in resonance. In the case of slower migration (100 Myr timescale N16 and KS16 simulations), Kozai cycling has time to modify particle orbits to lowere as the resonance slowly marches outward, resulting in preferential deposition of low-e particles on the sunward side of resonances. In the B&M simulation, Neptune's orbital circularization is more significant than in the N16 and KS16 simulations and occurs slowly enough that particles are Kozai cycled to lower e while the resonances narrow, leading to the more symmetric resonant beard structure we observe.

The low-e Neptune simulations thus far do not contain nearly as many test particles as the B&M simulation, so the details of resonance populations, in particular resonant island fractions, are not currently able to be compared in a statistically meaningful way. Since the B&M simulation predicts a large leading–trailing asymmetry among the n:1 resonances, in the future it would be very interesting to compare this to predictions of resonant island occupation from these grainy/smooth lower-e simulations.

The resonant fingers in our Figure 1 extending to low e by eye most closely match the fluffier near-resonant structure visible in the N16 and KS16 simulations where Neptune had a slower migration timescale (∼100 Myr). It appears that the most important factor in determining how well populated the near-resonant trails or beards are is the timescale of the evolution of Neptune's orbit relative to the Kozai cycling timescale for particles in Neptune's MMRs. Future detailed observations of the populations of TNOs within and immediately surrounding these distant resonances will be an important and powerful diagnostic for the exact mode and timescale of Neptune's migration.