© The Authors 2026

1 Introduction

Asteroids preserve some of the earliest solid materials to have formed in the Solar System, yet their current location in the relatively narrow region of the asteroid belt does not reflect the large chemical, compositional, and isotopic diversity of meteorites, defying any simple in situ formation scenario. Dynamical models strongly suggest that these bodies formed across a broad range of heliocentric distances, from the terrestrial planet region (∼ 1 au) outward to the Kuiper Belt (∼30 au), and were subsequently implanted into their present-day orbits by the gravitational influence and migration of the giant planets (Morbidelli et al. 2005; Bottke et al. 2006; Levison et al. 2009; Walsh et al. 2011; Raymond & Izidoro 2017b). Recent meteorite analyses have further supported this view, highlighting an isotopic dichotomy between non-carbonaceous (NCs) and carbonaceous chondrites (CCs), implying distinct formation reservoirs (Warren 2011; Budde et al. 2016; Kruijer et al. 2017). More subtle isotopic distinctions within the carbonaceous chondrite family, such as those between CI chondrites and other CC subtypes (Hopp et al. 2022; Hellmann et al. 2023; Spitzer et al. 2024), also suggest multiple reservoirs or formation epochs. Determining exactly where and when these bodies originated is thus critical for reconstructing the timing and mechanisms of giant planet formation, migration, and the delivery processes that shaped today’s Solar System architecture. Asteroids, therefore, hold the key to unlocking the chronology and dynamical evolution of planet formation.

Spectroscopic observations have linked specific CC groups to large (≳100km) asteroids in the main belt (e.g., Vernazza et al. 2016; DeMeo et al. 2022), indicating that CM-like (Mighei-like) planetesimals are Ch/Cgh-type asteroids, whereas CI-like (Ivuna-like) bodies comprise B, C, Cb, and Cg-type asteroids (DeMeo & Carry 2013; DeMeo et al. 2022; Broz et al. 2024; Anderson et al. 2025). CI chondrites are particularly interesting in their lack of chondrules and CAIs, along with their specific Fe and Ti isotopic anomalies, leading some to suggest a trichotomy between the NC-CC-CI types (Hopp et al. 2022; Hellmann et al. 2023; Spitzer et al. 2024).

Recently, it has been shown that the observed radial distributions of CM- and CI-like asteroids in the belt are distinct: CM-like objects exhibit an approximately symmetric distribution and are more abundant in the middle belt (Rivkin 2012; Anderson et al. 2025), whereas CI/Interplanetary dust particle (IDP)-like objects show a notably asymmetric distribution and dominate in the outer belt (see Fig. 1 or Anderson et al. 2025). The latter is similar to the distribution of comet-like (P-type) asteroids (DeMeo & Carry 2013), expected from latestage, outer-disk delivery (Nesvorný 2018). CI- and CM-like objects are also present in the belt at comparable overall levels (1.1:1 CI/CM). In our previous work (Anderson et al. 2025), we demonstrated that this dichotomy directly reflects two separate implantation events, reinforcing that CM and CI groups entered the asteroid belt from clearly distinct formation zones. We proposed that CM-like parent bodies formed in the Jupiter-Saturn region, in a pressure bump on the outer edge of the gap carved by Jupiter (<10au; Hellmann et al. 2023). According to this scenario, Jupiter’s early partial gap concentrated dust and pebbles at its edge, creating ideal conditions for local planetesimal formation. We demonstrated that a substantial fraction of these CM-like bodies were subsequently implanted into the asteroid belt, specifically during the rapid growth phase of Saturn, compatible with Raymond & Izidoro (2017b). By contrast, CI chondrites likely originated farther out, in the Uranus-Neptune zone (>15au, Hopp et al. 2022), and were implanted into the belt at a distinctly later stage when these ice giants underwent their final accretion and migration phases (Izidoro et al. 2015; de Sousa et al. 2020; Nesvorný et al. 2024; Anderson et al. 2025). This interpretation is consistent with their lack of chondrules, indicating that CIs directly accreted from primordial nebular dust that never experienced the intense heating events required to form chondrules.

Independent constraints also link CI-like material to a more distant and likely later-forming reservoir: CIs share buildingblock properties with cometary materials such as IDPs (100% matrix, ∼0% chondrules, ~0% CAIs; Nakamura et al. 2023; Nakashima et al. 2023; Zega et al. 2025; Bradley 1994), and anhydrous areas in CI chondrites exhibit comet-like spectral behavior (Brunetto et al. 2023). The chronology is consistent with this picture, with CM ages clustering around ~3-4 Myr after CAIs (Tanaka et al. 2024), whereas IDPs, comet-like bodies, and trans-Neptunian objects are generally associated with formation after ≳4-5 Myr (Neveu & Vernazza 2019; Davidsson et al. 2016; Bierson & Nimmo 2019). Finally, the spectral similarity among P/D-type asteroids, comets, Jupiter Trojans, and Centaurs (Vernazza et al. 2025), together with models favoring delayed Uranus-Neptune formation (∼4-6Myr after CAIs; Dodson-Robinson & Bodenheimer 2010; Helled & Bodenheimer 2014), naturally places CI-like delivery in a later evolutionary phase than CM-like implantation. Taken together, these constraints motivate treating the CM and CI source regions as radially separated reservoirs with different delivery windows, rather than neighboring subpopulations of a single carbonaceous reservoir.

While the exact mechanism of chondrule formation remains uncertain, they require intense heating and high-pressure conditions to form, implying that their production was confined to regions closer to the young Sun and not in the distant ice giant zones (Ali-Dib & Walton 2025). The relative lack of intermediate classes (CR, CK, etc.) in the belt raises the question of whether the CI formation region was geographically isolated: In short, to what extent the CI reservoir was truly “quarantined” from other CC groups, and whether CM-like material could have been scattered and implanted outward. To our knowledge, this specific phase of outward planetesimals transfer scattered by Jupiter-Saturn growth before Uranus and Neptune began significant accretion and migration has not been explicitly quantified in previous dynamical studies.

Here, we explore the potential outward transport of >100-km-sized CM-like planetesimals initially formed on the outskirts of Jupiter’s gap to the Uranus-Neptune region. By quantifying the fraction of bodies that were scattered outward during Saturn’s growth, we constrain the efficiency of outward transfer of CM-like planetesimals and assess whether any of this material can be retained on stable orbits. This provides a dynamical upper bound on the availability of chondrule-rich inner carbonaceous material to the Uranus-Neptune region and, by extension, to any populations sourced from that region. Section 2 describes the details of our dynamical simulations, including the respective disk model, planetesimal-disk interactions, and gravitational scattering by proto-planets. We present the main results in Section 3, focusing on the outward-transport efficiencies and timescales, followed by a discussion of their implications for the compositional diversity of the outer Solar System in Sect. 4. Finally, we summarize our conclusions in Section 5.

Fig. 1 Cumulative distributions of the semimajor axis N(≤a) of >100 km asteroids classified as CM-like, CI/IDP-like (including B types with the geometric albedo of pv < 0.1), and S-type. Vertical dotted lines mark major mean-motion resonances with Jupiter. CM-like objects exhibit an approximately symmetric distribution, whereas CIlike objects show a notably asymmetric distribution. Adapted from Anderson et al. (2025).

2 Methods

Raymond & Izidoro (2017a) demonstrated that, in a growing giant-planet system, planetesimals can be scattered to high eccentricities along curves of constant Tisserand parameter. Under sufficient gas drag, smaller bodies (1-10 km) have their eccentricities rapidly damped and can survive on stable orbits in the outer disk. However, the gas density beyond Saturn’s orbit is too low to circularize large (<100km) planetesimals before they are ejected. Consequently, they are far more likely to be lost from the system. Building on their approach, we adopt a similar conceptual framework for outward scattering but focus on two new aspects: (i) varying the physical conditions of the protoplane-tary disk, (ii) exploring a range of planet formation timescales in which Saturn and additional ice giant cores grow from ~1 M_⊕ to their final masses, while allowing the growing cores to migrate inward, and (iii) employing different gas profiles to examine the impact of the density of the environment plays in the circularization of planetesimals’ orbits. Moreover, this model is anchored in the idea that Saturn forms at the outer edge of Jupiter’s gap, consistent with pressure-bump formation scenarios (e.g., Hellmann et al. 2023; Lau et al. 2024; Xu & Wang 2024).

2.1 Initial state of the system

As Jupiter’s core grew to about 20 M_⊕, it began carving a gap in the protoplanetary disk, thereby inducing a local pressure bump that prevented pebbles from drifting inward (Lambrechts et al. 2014). This gap grew deeper as Jupiter accumulated mass, further isolating its vicinity from the surrounding disk (Lin & Papaloizou 1986; Lee & Peale 2002; Crida & Bitsch 2017). Eventually, Jupiter settled into a Type-II migration regime, slowing its radial motion and allowing pebbles to steadily accumulate at the outer edge of its gap (Gonzalez et al. 2015). This region experienced an elevated dust-to-gas ratio, creating favorable conditions for efficient planetesimal formation via cohesion or streaming instability. Consequently, a robust planetesimal reservoir is expected to have formed beyond Jupiter’s gap, while local material close to the planet’s orbit was largely depleted.

Previous works proposed that Saturn may have formed at or near this outer boundary (Kobayashi et al. 2012; Lambrechts et al. 2014; Gonzalez et al. 2015; Lau et al. 2024), which is consistent with evidence that Saturn once orbited closer to Jupiter than its present-day location (Deienno et al. 2017). As Saturn’s core grew, its increasing gravitational reach scattered these newly formed planetesimals inward, supplying solids to Jupiter’s circumplanetary disk (CPD), potentially contributing to the formation of the Galilean satellites (Ronnet et al. 2018). Many planetesimals were delivered toward the inner Solar System, either captured on stable orbits within the asteroid belt or flung farther inward, possibly aiding water delivery to the terrestrial planet region (Meech & Raymond 2020). The relative importance of each of these pathways depends strongly on gas-drag strength (Raymond & Izidoro 2017b), which is itself influenced by both the disk density - thus its temporal evolution - and the sizes of the planetesimals, with smaller bodies experiencing more pronounced aerodynamic forces (Raymond & Nesvorný 2022). In general, outward scattering is difficult, since gas drag tends to shrink the orbit toward pericenter.

In the simulations presented here, we adopt the same disk model and growth timescales as in Anderson et al. (2025), originally based on Ronnet et al. (2018). While Jupiter is assumed to be fully formed, we grow all other planets from 1 M_⊕ to their final mass through a simplified gas accretion prescription (Sect. 2.2). The dynamical evolution of planetesimals is traced under the influence of gas drag, gravitational interactions with proto-planets, and close encounters that can lead to scattering or capture (Sect. 2.3).

2.2 Giant planet growth and dynamics

We carried out orbital integrations using the open-source REBOUND library (Rein & Liu 2012) and its high-accuracy, adaptive-timestep IAS15 integrator. These simulations include two to three giant planets and a swarm of 10 000 massless test particles that serve as proxies for planetesimals. The simulations are evolved with a timestep of 10⁻²/(2π) yr. Test particles that collide with a planet or acquire hyperbolic trajectories are removed automatically to improve computational efficiency. Terrestrial planets are omitted from the model due to their comparatively small masses and more compact orbits, which minimally affect the dynamical evolution of the outer proto-planetary disk while requiring more integration steps and longer calculations.

Initial conditions. We begin with Jupiter placed at a heliocentric distance of 5.43 au, consistent with post-disk evolution scenarios for the gas giants. Saturn’s core is positioned at 7.37 au, near the outer edge of Jupiter’s gap and in a 3:2 resonance with Jupiter. In the three-planet scenario, an additional “ice giant” core is initialized at 9.67 au to form a 3:2, 3:2 multi-resonant chain based on Deienno et al. (2017). These orbital configurations align with established constraints from the current planetary configuration, the Jupiter Trojans or the irregular satellites, and the secular tilt resonance within the cold Edgeworth-Kuiper belt observed today (Nesvorný 2018; Baguet et al. 2019). While Ronnet et al. (2018) also tested scenarios in which Saturn is formed further out and migrates inward, we do not explore that alternative here, as it has been shown to be less consistent with Saturn’s final growth pattern, as well as with the effective delivery of solids into Jupiter’s circumplanetary disk and the asteroid belt. Instead, we focus on Saturn’s in situ formation near the gap edge, which is consistent with a growing body of work arguing that giant planets commonly emerge in local pressure maxima (e.g., Hellmann et al. 2023; Lau et al. 2024; Xu & Wang 2024). While these planetary cores are likely to migrate during their formation while gas remains before assuming the positions that will lead to the Nice model configuration Morbidelli et al. (2005); Nesvorný (2018); Deienno et al. (2017), it is the most stable point from which to start.

Growth. At the outset of the simulations, Jupiter is taken to have its present-day mass (≈300 M_⊕). The other planetary cores begin at 1 M_⊕ and increase in mass over a prescribed timescale τ_growth, with Saturn reaching approximately 100M_⊕ and the additional ice giant core growing to roughly 15 M_⊕, comparable to the present-day mass of Uranus. This simplified growth prescription is adapted from Ronnet et al. (2018) and models the planet’s mass evolution as described by $$M_{\text{planet}}(t)=M_i+\mathrm{\Delta }M(1-{\mathrm{e}}^{-t/{\tau }_{\text{growth}}}),$$(1)

where M_i is the initial core mass and ∆M is the difference between each planet’s initial core mass and the final mass. This approach approximates standard core-accretion models, in which a gaseous envelope gradually contracts, transitions to runaway accretion, and then tapers off as the planet begins opening a gap (Pollack et al. 1996). We explore the impact of the formation timescale with τ_growth set to 1 × 10⁵, 5 × 10⁵, and 1 × 10⁶. Although more detailed models are needed to capture all aspects of giant planet growth, a precise reconstruction is beyond the scope of this study.

Although the exact timing for the emergence of an ice giant embryo relative to Saturn’s formation is not fully constrained, we initialize both cores together to explore how an additional growing planet might affect the outward scattering of CM-like bodies. This choice is partly motivated by computational efficiency, as a second core emerging later would require restarting the simulation at a different epoch, and by the recognized possibility of multiple protoplanetary embryos coexisting in resonance chains (Guilera et al. 2011; Lau et al. 2024). In practice, setting the ice giant core to grow in parallel with Saturn captures a range of plausible scenarios, whether this additional planet evolves into Uranus, Neptune, or the hypothetical fifth giant planet (so-called “Ice 1”) postulated by the Nice model (Nesvorný 2018; Deienno et al. 2017). Although simplified, this approach allows us to bracket the influence of a simultaneous second core on scattering dynamics while avoiding uncertainties about the precise timing or radial location of a late-emerging ice giant.

Migration and damping. We account for the migration and damping of the growing planets due to gas-disk interactions, adopting fictitious accelerations as described by Ronnet et al. (2018), based on Cresswell & Nelson (2008). Specifically, we employed the following acceleration terms for migration, eccentricity, and inclination damping: $${\mathbf{a}}_{\text{mig}}=-\frac{\mathbf{v}}{{\tau }_{\text{mig}}},$$(2) $${\mathbf{a}}_{\text{e}}=-2\frac{(\mathbf{v}\cdot \mathbf{r})\mathbf{r}}{r^2{\tau }_{\text{e}}},$$(3) $${\mathbf{a}}_{\text{i}}=-\frac{v_z}{{\tau }_i}\mathbf{k},$$(4)

where v is the velocity vector of the planet, r is its position vector, k is the unit vector in the z direction, and r is the heliocentric distance.

The prescription of the migration timescale is given by $$\begin{array}{rl}{\tau }_{\text{mig}}=& \frac{2t_{\text{wave}}}{2.7+1.1\beta }(\frac{H}{r}{)}^{-2}[P(e)+\frac{P(e)}{|P(e)|}\\ & \times (0.070(\frac{i}{h})+0/085(\frac{i}{h}{)}^4-0.080(\frac{e}{h})(\frac{i}{h}{)}^2)],\end{array}$$(5) $$P(e)=\frac{1+(\frac{e}{2.25h}{)}^{1.2}+(\frac{e}{2.84h}{)}^6}{1-(\frac{e}{2.02h}{)}^4},$$(6)

where i is the orbital inclination, e the eccentricity, h = H/r is the disk’s aspect ratio, M_* denotes the stellar mass, m_p the planet’s mass, a_p its semimajor axis, and Ω_p the Keplerian frequency. Although the sign dependence of P(e) allows torque reversal at sufficiently high eccentricities, our planets never enter that regime, and no outward migration arises. We set β to 1.5 to be consistent across all gas profiles.

The damping timescale is given by $$t_{\text{wave}}=\frac{M_{\ast }}{m_{\text{p}}}\frac{M_{\ast }}{{\mathrm{\Sigma }}_{\text{p}}{a}_{\text{p}}^2}{h}^4{\mathrm{\Omega }}_{\text{p}}^{-1}.$$(7)

The eccentricity damping time, τ_e, and inclination damping time, τ_i, are then calculated via $${\tau }_e=\frac{t_{\text{wave}}}{0.780},[1-0.14(\frac{e}{h}{)}^2+0.06(\frac{i}{h}{)}^3+0.18(\frac{e}{h}{)}^2(\frac{i}{h})],$$(8) $${\tau }_i=\frac{t_{\text{wave}}}{0.544}\times [1-0.33(\frac{i}{h}{)}^2+0.24(\frac{i}{h}{)}^3+0.14(\frac{e}{h}{)}^2(\frac{i}{h})].$$(9)

We neglect higher-order torques (e.g., heating torques from rapidly accreting protoplanets Benítez-Llambay et al. 2015; eccentric forcing from nonzero-eccentricity gaps Fendyke & Nelson 2014; or hot-trail effects Eklund & Masset 2017; Chrenko et al. 2017; Cornejo et al. 2023) since our primary aim is to track the implantation of small bodies rather than to optimize an exact migration prescription (Paardekooper et al. 2010, 2011).

This simplified disk model and set of growth prescriptions allowed us to study the implantation of planetesimals in the inner Solar System (e.g., Ronnet et al. 2018; Anderson et al. 2025), and will now allow us to focus on whether and how Saturn’s formation efficiently scatters CM-like planetesimals outward, potentially delivering them to the Uranus-Neptune zone.

2.3 Small body dynamics

Observations and dynamical studies suggest that many small bodies populating the outer Solar System originated between Jupiter’s orbit and approximately 30 au (Tsiganis et al. 2005; Gomes 2003; Levison et al. 2008; Kaib & Quinn 2008). Here, we focus on the region immediately beyond Jupiter, where the CM-like bodies in this scenario are hypothesized to form. We randomly initialized 10 000 massless test particles representing planetesimals uniformly in semimajor axis in the 7<a<10 au range (two-planet simulation) and the 7<a<11 au range (three-planet simulation). This simplified approach does not attempt to reproduce any particular radial distribution of solids, which is poorly constrained and likely consisted of multiple pressure bumps or ring-like substructures (Izidoro et al. 2022). Instead, it ensures that each radial slice in that zone is adequately sampled so we can identify how far outward these bodies might travel in response to planetary perturbations. Varying the precise radial or size distribution would not drastically change the likelihood of their scattering since planetesimal-planetesimal interactions are negligible in comparison to the gravitational and aerodynamic forces exerted by the planets and the gaseous disk (Raymond & Izidoro 2017a).

We restricted the test particles to bodies of 100km in diameter. As shown by Raymond & Izidoro (2017a), smaller (∼ 1-10 km) planetesimals are more likely to be captured onto stable orbits by aerodynamic drag once scattered outward, but they generally do not travel as far from their initial orbits (∼5 au outward). In contrast, larger bodies can, in principle, be scattered to greater heliocentric distances, although the weaker drag force they experience makes it harder for them to circularize and remain in the system. By concentrating on these 100 km objects, we effectively probe a “worst-case” limit for outward migration: if these 100 km objects fail to be implanted in the Uranus-Neptune zone, we can infer that CM-like material of similar or larger sizes is unlikely to have contributed significantly to this region.

Aerodynamic drag. We investigated how the aerodynamic drag influences each planetesimal using the methods described by Ronnet et al. (2017) by applying the following acceleration term $${\mathbf{a}}_{\text{drag}}=-\frac{1}{t_{\text{s}}}(\mathbf{v}-{\mathbf{v}}_{\text{g}}),$$(10)

where v is the planetesimal’s velocity, v_g is the local gas velocity, and t_s is the stopping time. Following Perets & Murray-Clay (2011) and Guillot et al. (2014), we calculated t_s for 100km-sized bodies of density ρ_s = 1 g cm⁻³ as $$t_{\text{s}}={\left(\frac{{\rho }_{\text{g}}{v}_{\text{th}}}{{\rho }_{\text{s}}{R}_{\text{s}}}\text{min}[1,\frac{3}{8}\frac{v_{\text{rel}}}{v_{\text{th}}}{C}_{\text{D}}]\right)}^{-1},$$(11)

where R_s is the size of the solids, here set to 100km, and ρ_s = 1 g cm⁻³ their density. The gas volumetric density, ρ_g, is obtained from the surface density, Σ_g, by assuming hydrostatic equilibrium in the vertical direction, with a disk aspect ratio of h ≡ H/r = 0.05 (Cresswell & Nelson 2008). The relative velocity between the gas and the planetesimal is v_rel, while the gas thermal velocity is given by $v_{\text{th}}=\sqrt{8/\pi }{c}_{\text{s}}$, where c_s is the isothermal sound speed. In this context, the isothermal sound speed is defined as c_s = ΩH, Ω being the Keplerian orbital frequency, expressed as $\mathrm{\Omega }=\sqrt{G{M}_{\odot }/r^3}$. The dimensionless drag coefficient, C_D, was computed as a function of the Reynolds number Re of the flow around the planetesimal (Perets & Murray-Clay 2011) $$C_{\text{D}}=\frac{24}{\mathrm{R}\mathrm{e}}(1+0.27\mathrm{R}\mathrm{e}{)}^{0.43}+0.47(1-{\mathrm{e}}^{-0.04{\mathrm{R}\mathrm{e}}^{0.38}}),$$(12) $$\mathrm{R}\mathrm{e}\equiv \frac{4R_{\text{s}}{v}_{\text{rel}}}{c_{\text{g}}{l}_{\text{g}}}.$$(13)

Using the mean free path of the gas, l_g, from Supulver & Lin (2000) and the size of 100 km, we find a value of C_D = 0.44, which is then fixed as a constant to save computational time.

Fig. 2 Comparison of the gas profiles used in the various simulations. The initial orbits of Jupiter, Saturn’s core, and the additional ice giant core are indicated by vertical gray lines. Beyond 10 au, the profiles can be seen as Flat (Black, Raymond & Izidoro 2017b), Sloped (Blue, Cresswell & Nelson 2008), and with a Sharp Drop-off (Orange, Desch et al. 2018).

2.4 Gas profiles

Accurately defining the disk’s gas distribution is notoriously challenging due to the complex interplay of processes such as radiative transfer, dust condensation, viscous heating, and planetdisk interactions (e.g., Nelson et al. 2000; Bitsch et al. 2015; Johansen et al. 2015). As a result, no single disk model has emerged as definitive. In our previous study Anderson et al. (2025), we showed that subtle differences in adopted gas profiles can have a profound impact on how, when, and where planetesimals are implanted into the inner Solar System. Here, we extend that analysis to the outer region beyond Saturn, where the physical conditions are even less certain. We employ the same three different disk models in order to capture the range of outcomes for planetesimal scattering and planetary embryo migration. A comparison of all these gas surface density profiles is shown in Fig. 2. Since the gas profile of the outer disk is even less constrained than the inner system, we can describe them as follows:

Sloped profile. A commonly used 1D model following a powerlaw decrease in gas surface density, typically with Σ_g ∝ r^−0.5 normalized to Σ_g = 300 g cm⁻² at 1 au (Ronnet et al. 2018) along with the aspect ratio h = 0.05 (Cresswell & Nelson 2008). Jupiter’s partial gap is introduced by reducing the surface density by a factor of 10⁻¹ near the planet’s orbit (see also Bitsch et al. 2015), mimicking the deepening of the gap as Jupiter is already close to its final mass.

Flat profile. Raymond & Izidoro (2017b) assumed the surface density profile from Morbidelli & Crida (2007), starting after Jupiter and Saturn were both formed. The starting disk profile has a Σ_g of 200 g cm⁻² at 1 au, a gap opened by Jupiter at 5.4 au, and a partial gap opened around Saturn at 7 au. While they decreased the surface density uniformly in radius, on a 2 × 10⁵ yr exponential timescale, until it was entirely removed by 2 × 10⁶ yr, we chose to examine their starting disk profile, as the gap opened by Jupiter and Saturn was already much wider than in the other profiles. The model is practically flat beyond 10 au, and allows us to test how a higher residual gas density may allow some larger solids to damp their eccentricities enough to survive outward scattering, although it also influences the migration speed of the growing embryos.

Sharp drop-off profile. Desch et al. (2018) considered a selfsimilar disk model (Hartmann et al. 1998) that evolves in time through viscous diffusion. It features a higher initial surface density at small radii and a steep slope beyond 8 au. While they examined a gap location of ~3 au for Jupiter, we apply a similar prescription to a gap at 5.4 au, consistent with our Saturn-growth scenario and the notion that Jupiter opened its gap closer to its current orbit. In this model, the precipitous decline in Σ_g beyond Saturn reduces the efficacy of aerodynamic damping for large planetesimals, potentially limiting their ability to remain in stable orbits at high eccentricities. We select their 3 Myr gas profile, as it is an intermediary profile, and should approximate the time after Jupiter has formed.

The disk viscosity is parameterized by the classical α prescription (Shakura & Sunyaev 1973), with a typical α ≈ 2 × 10⁻³. We do not investigate this parameter as we found it has a modest effect on our previous simulations. More details can be found in our Anderson et al. (2025) paper as well as the respective papers that originated the individual profiles.

By comparing these three scenarios, we aim to bracket the potential range of outcomes for how far CM-like planetesimals might be transported outward and whether they can be implanted in the Uranus-Neptune range. Given the uncertainties inherent to protoplanetary disk evolution at large radii, these models should provide a robust test of how the outer disk profile shapes the fate of planetesimals initially located near Jupiter’s gap edge.

3 Results

3.1 Two planet scenario

In all of the two-planet runs, type-I migration drives a modest inward migration of Saturn, though this inward shift amounts to only a few tenths of an astronomical unit by the end of the simulations. However, under the Desch et al. (2018)-inspired disk profile, where the surface density remains higher < 10 au, both Jupiter and Saturn migrate more substantially, essentially because the deeper reservoir of gas exerts a stronger net torque on them.

As expected, the gravitational influences of Jupiter and the growing Saturn excite the eccentricities of nearby planetesimals, leading to scattering both inward (toward the asteroid belt and terrestrial region) and outward (toward semimajor axes beyond 10 au). Right away, Saturn’s core scatters many planetesimals that happen to lie within roughly ±1 au of its orbit, depending on initial configuration (Fig. 3). These objects follow characteristic “wings” in (a, e)-space (Fig. 4), corresponding to curves of constant Tisserand parameter (Levison & Duncan 1997).

While we ran the simulations out to 1 × 10⁶ years, we note that the most significant scattering and planetesimal removal occurred within the first 2.5 × 10⁵ years. By this point, the bulk of strongly perturbed bodies have either collided with Saturn, become trapped within Jupiter’s circumplanetary region, or been ejected from the system entirely. Subsequent evolution, beyond a few hundred thousand years, largely results in the slow inward migration of Saturn without any appreciable new scattering events. Since many simulations also lost all of their planetesimals before reaching the 1 Myr mark, we focus on outputs after 2.5 × 10⁵ years, where the fate of individual planetesimals is already well established and further integration offers diminishing returns. A direct comparison of the final distributions across the different gas profiles at t = 300 kyr is provided in Table 1.

Initial position. All planetesimals that have been completely removed from the system originated in the 7.2-8.6 au zone (Fig. B.1). These represent a quarter to a third of all initialized planetesimals. Their fates vary: some are accreted onto Saturn, some are driven onto orbits that deliver them to Jupiter’s circumplanetary disk, and others are ejected from the Solar System. Any planetesimals formed in this range that remain after this initial scattering event typically settle onto stable orbits either in the inner Solar System, becoming implanted in the asteroid belt or even the terrestrial formation zone, or occupy more distant, excited orbits, in the 10-20 au range.

Planetesimals having formed beyond ∼8.6 au generally retain semimajor axes close to their initial locations but develop higher eccentricities. Although Saturn stirs these outer objects, it does not scatter them strongly inward or outward on these timescales. Meanwhile, any bodies initiated inside ∼7.5 au are quickly caught co-orbiting Saturn (Fig. 4) on trajectories that would ultimately deliver them into Saturn’s Hill sphere and feed the circumplanetary region. In these simulations, which do not include an explicit circumplanetary disk, such bodies are treated as ‘absorbed’ once they collide with the planetary surface. In principle, a more detailed prescription that includes a circumplanetary disk (as performed in Ronnet et al. 2018) could track whether some fraction remains in bound orbits around Saturn, potentially contributing to satellite formation.

Growth timescale. An additional parameter we tested is the growth timescale, τ_growth, for Saturn. Unsurprisingly, faster core accretion means that Saturn reaches a higher mass sooner, which tends to enhance its migration and early scattering of nearby planetesimals. In simulations with a faster τ_growth, a third of the test particles are removed, while a slower τ_growth allows us to retain three quarters. Yet these results show that the final distribution of surviving objects is almost indistinguishable across different τ_growth values. The system simply evolves more quickly when Saturn is growing faster but ultimately settles into a similar end state, with the same fraction of planetesimals stranded in the inner Solar System or scattered outward. In particular, for growth timescales ranging from τ_growth = 1 × 10⁵, 5 × 10⁵, and 1 × 10⁶ years, the overall clearing zone and the fraction of removed bodies remain roughly the same. This implies that how quickly Saturn reaches its final mass has little bearing on which planetesimals survive, aside from minor variations in their removal timescales.

Saturn’s orbital evolution remained modest, consistent with the findings of Raymond & Izidoro (2017b), thanks to its resonance lock with Jupiter and the partial gap carved by the two planets. Planetesimals within roughly 1 au on either side of Saturn’s orbit were rapidly scattered, following constant-Tisserand pathways and eventually reaching either the inner Solar System or high-eccentricity orbits in the 10-20 au region. However, since neither the surface density nor the gas drag was reduced over time in our simulations, the abundant gas supply acted to circularize many orbits rather than eject them outright, allowing a fraction of larger bodies to survive beyond 20 au. Interestingly, we also find that scenarios with lower gas density beyond ∼ 10 au can allow scattered planetesimals to reach even larger semimajor axes, yet these same runs yield fewer objects with e < 0.4 on stable orbits. Strong aerodynamic damping is needed to retain them at such distances: a higher disk gas density increases the fraction of bodies retained on orbits beyond 10 au, since aerodynamic drag can circularize orbits before they are fully ejected.

Because gas drag damps eccentricity at roughly fixed q, scattered bodies tend to fall back toward their pericenters rather than circularize at large a. In the sharp drop-off case, no particles attain q > 10 au ; in the sloped (canonical) disk, <1% reaches q ≃ 12 au ; and in the flat disk, <1% achieve q ≃ 15 au. Even when a > 20, the small q means the orbits remain vulnerable to loss. This behavior is also consistent with the two-planet, gas-rich disk experiments of Ribeiro de Sousa et al. (2024), who examined scattering during Jupiter-Saturn growth in an Raymond & Izidoro (2017b)-like disk context. Although their study emphasizes asteroid-belt implantation rather than retention in the ice giant zone, their reported outer-disk outcomes similarly show that essentially none are circularized onto stable orbits beyond 12 au over million-year timescales.

Moreover, suppose the outer disk contained partial gaps or ring-like density drops beyond Saturn, such as those that might arise at ice lines. In that case, these results indicate there would be little chance for substantial outward transfer of CM-rich material. In such lower-density regions, aerodynamic drag and scattering events would be too weak to bridge the ∼8.6 au boundary effectively, as evidenced by the gas-poor, sharp drop-off profile, where 0% of the initial CM population can reach stable e < 0.4 orbits in the 15-20 au zone. Hence, the presence of any pronounced gas deficit outside Saturn’s orbit would further limit CM-like objects from migrating into the Uranus-Neptune zone, reinforcing the idea that distinct ring structures in the disk could isolate compositional reservoirs until another core-forming event reshaped the local dynamics.

Overall, even in the most gas-rich scenario, only a modest fraction of planetesimals survive outward scattering with sufficiently low eccentricity to remain long-term (<2%). Most remain on highly eccentric paths or else are lost from the system, implying that a truly robust implantation of large CM-like bodies into the Uranus-Neptune region requires a combination of high disk density and additional dynamical events (e.g., interactions with later-forming embryos). Moreover, the very gas levels required to help particles remain would also enhance Type-I migration of embryos, further pulling pericenters inward and working against durable implantation. In sum, only exceedingly high outer-disk gas would have produced any outward retention at all, and even then at a percent-level at best; the contribution of CM-like material to the Uranus-Neptune formation zone and the present-day TNO population would therefore be negligible.

Taken together, these findings confirm that the two-planet configuration allows Saturn to scatter a significant fraction of CM-like bodies both inward to the asteroid belt and also outward into the 10-20 au region, though whether or not their orbits will circularize depends on the amount of gas remaining in the Uranus-Neptune formation zone. Neither the precise growth rate of Saturn nor the final integration time alters the distribution of scattered planetesimals. Once the scattering window closes around 2.5 × 10⁵ years, the fate of each particle is, for all practical purposes, sealed. We estimate planetesimalplanetesimal collisions in Appendix A.

Fig. 3 Fraction (%) of planetesimals from each initial formation zone that end up within 1-au-wide semimajor-axis bins at t = 300, kyr for the two-planet model (log scale). We consider implanted objects as those with e < 0.4. In all cases, most bodies formed between 9 and 10 au remain near their birth locations, reflecting a dynamically quiet zone largely unaffected by Saturn’s growth. The prominent blue feature near −7 au corresponds to particles trapped in Saturn’s co-orbital region.

Fig. 4 Eccentricity of remaining planetesimals as a function of semimajor axis at t = 300 kyr for each gas profile and Saturn growth timescale for our two-planet model. The dotted line represents Jupiter-crossing orbits, the dashed line represents Saturn-crossing orbits. The color represents the planetesimals’ formation location, with cooler colors representing earlier loss and yellow for those that remain. The more gas in the simulation, the more easily the planetesimals are circularized (lower final e). The cluster of purple points near ~7 au represents planetesimals trapped near Saturn’s orbit; these bodies are dynamically influenced by the planet and are likely to be absorbed into its circumplanetary disk over time.

Fig. 5 Same as Fig. 3, but for the three-planet model. We consider implanted objects as those with e < 0.4. Faster planetary growth clears the local planetesimal population. Objects that formed beyond 10 au are less likely to interact with the forming planetesimals.

3.2 Three planet scenario

In the three-planet scenario, we extend the two-planet model by adding an ice giant core at 9.6 au, which grows from 1, M_⊕ to an approximately Uranus-mass (∼15M_⊕) over the prescribed τ_growth. Whether this core represents a hypothetical “Ice 1” planet postulated by recent Nice Model variations (Deienno et al. 2017) or an early Uranus-analog would require migration to match the current orbital configuration of the ice giants. We initialize planetesimals uniformly between 7 au and 11 au to capture the broader region where multiple planetary embryos could form.

The additional core’s dynamical impact is qualitatively similar to Saturn’s. It excites eccentricities for planetesimals in the vicinity of its orbit, scattering them along constant-Tisserand pathways. In the flat disk model (where Jupiter and Saturn carve wide gaps), the ice giant’s net inward drift remains moderate, and its scattering zone closely mirrors that of Saturn in the two-planet case. By contrast, since the sharp drop-off model sees intense inward migration due to the high amount of gas in the inner system, this core even manages to push Jupiter inward by nearly 1 au, demonstrating how strongly the gas reservoir can affect overall giant-planet architecture.

We observe a somewhat stronger dependence on the growth timescale when the additional core is present, as is demonstrated in Fig. 5. Rapidly forming ice giant cores reach higher masses earlier, thereby scattering more distant planetesimals and promoting greater orbital reshuffling at higher heliocentric distances. In slower-growth runs, the ice giant’s milder migration exerts comparatively less perturbation on bodies beyond ∼10.4 au. Nevertheless, once scattered planetesimals do enter the ice giant’s sphere of influence, they show a similar pattern of partial circularization to that seen in the two-planet case (Fig. 6). Notably, objects initially beyond Saturn’s scattering limit (∼8.6 au) can now be pushed deeper into the Uranus-Neptune zone and beyond (>30 au) with eccentricities e < 0.4 in the 15-20 au range. However, their qs remain low, none exceeding q > 20 au for any gas profile.

Overall, the presence of this third core modestly widens the region of active scattering compared to the two-planet scenario, particularly under high-gas conditions and rapid core growth (Table 2). At the same time, once again, objects forming even further out typically remain out of reach unless additional planetary embryos emerge or the local gas density remains sufficiently high for longer, allowing repeated scattering and aerodynamic damping to take hold.

Fig. 6 Same as Fig. 4, but for the three-planet model. The dotted line represents Jupiter-crossing orbits, and the dashed line represents Saturn-crossing orbits. The color represents the planetesimals’ formation location, with cooler colors representing earlier loss and yellow for those that remain. The more gas in the simulation, the more easily the planetesimals are circularized.

3.3 Reservoir contamination

Even if only a small fraction of CM-like planetesimals can be scattered outward and circularized, it is useful to compare that contribution to the scale of the reservoir it would enter. In our working picture, the CM source region is comparatively narrow (7-11 au), whereas the CI-like source region is typically placed farther out and spans a broader annulus (15-25 au; Hopp et al. 2022), possibly already containing a native mass of order M_CI ~ 20 M_⊕. Adopting the CM mass budget inferred in Anderson et al. (2025) for loose planetesimal material between Saturn and Uranus (M_CM,tot ~ 1 M_⊕), our measured upper limit of ≲2% (two-planet) to ≲4% (three-planet) of 100-km bodies attaining long-lived, low-e orbits beyond Saturn implies at most M_CM ≲ 0.04 M_⊕ reaching the 15-25 au region, i.e., a mass fraction ≲2 × 10⁻³ once diluted into the native ring. Such a contribution is too small to measurably affect the observed CM-CI dichotomy in the asteroid belt; any CM material that briefly resided in the outer ring would be both rare and, if it experienced outer-disk processing, potentially difficult to recognize unambiguously in the present-day, implanted population.

3.4 Comparison with the five planet scenario

In our previous Anderson et al. (2025) study, we explored a model in which five planetary cores formed simultaneously.

Additional cores were placed for Uranus at 15.4 au and Neptune at 20.3 au in the “Wide” configuration (3:2, 3:2, 2:1, 3:2 resonance chain), and at 12.6 au and 16.2 au, respectively, in the ‘Tight’ variant (3:2, 3:2, 3:2, 3:2 resonance chain) Deienno et al. (2017). Apart from these added cores, the underlying gas profiles (flat, sloped, or sharp) and τ_growth were analogous to those presented in this work. A natural consequence of having gas present in the disk is that all planetary embryos migrated inward and drastically reshaped the distribution of planetesimals. Whereas the previous paper focused on the delivery of planetesimals to within 5 au, we can also re-examine these scenarios from the perspective of potential CM-like bodies that survived beyond 10 au. For consistency with our current work, we will only focus on planetesimals with the initial starting semimajor axis, a_init, between 7 and 11 au. These results are summarized in Table 3. Notably, none of these runs produce long-lived, low-eccentricity objects on stable orbits beyond 20 au.

What stands out in the five-planet scenario is that if Uranus’s and Neptune’s cores are allowed to migrate (Fig. 7) while Saturn is still forming, while many planetesimals are ejected, a non-negligible fraction are quickly circularized into stable orbits between 15-20 au within 300 kyr. This is a clear departure from the two- and three-planet cases. In a flat gas profile, more than 10% of the original 7-11 au population survives on relatively stable orbits in the 15-20 au region, with Uranus and Neptune helping to capture and damp their orbits.

To isolate the impact of migration, we ran a control “Wide, no-migration” simulation, in which the ice giant cores’ semimajor axes were held fixed while they accreted (akin to Nesvorný et al. 2024). In this case, less than <1% of the 7-11 au cohort is circularized between 15-20 au by 300 kyr. The contrast arises because migrating cores sweep across the 10-20 au band, providing repeated, gentle encounters that lengthen planetesimals’ residence time at large a, so gas drag can damp e. Since more gas leads to faster core migration and stronger drag, this also leads to a larger stabilized fraction (up to 10%), whereas suppressed migration yields minimal implantation (<1%).

4 Discussion

Confronting simulations with observations. Cycle 1–4 JWST/NIRSpec observations of ∼100 Centaurs and TNOs, primarily from the DiSCo-TNOs program (e.g., Pinilla-Alonso et al. 2025; Licandro et al. 2025), but also from complementary GO and GTO programs (e.g., Emery et al. 2024; Holler et al. 2025), have revolutionized our understanding of these populations by revealing remarkable compositional diversity in volatiles and organics across the outer Solar System. Cluster analyses consistently identify three distinct compositional groups among TNOs – H₂O-type Bowl-TNOs, CO₂-type Double-dip-TNOs, and organics-type Cliff-TNOs (Holler et al. 2025; Pinilla-Alonso et al. 2025) – that correlate strongly with molecular abundances. In particular, no phyllosilicates were identified at the surfaces of these bodies, whose refractory component appears totally anhydrous based on JWST/MIRI observations (Vernazza et al. 2025).

If CM bodies are present among TNOs, they must be hidden in their interior. However, this scenario is not supported by near- and mid-infrared observations of Centaurs and TNOs. Midinfrared observations show the systematic presence of anhydrous dust at the object surfaces, which appears similar to CP IDPs (Vernazza et al. 2025). Yet, CP IDPs differ in texture from CM chondrites, being chondrule-poor. A surface consisting of the anhydrous precursors of CM chondrites would somewhat resemble spectrally CO and CV chondrites, which show a mixture of matrix and chondrules; that is, an olivine band would be clearly present in the near-infrared range. Such a band is, however, not seen in the aforementioned surveys. Note that the presence of CI-type material is not seen either among current surveys of Centaurs and TNOs. However, CI-like material is present at the surfaces of some irregular satellites of Jupiter, which are thought to be, similarly to the Jupiter Trojans, former TNOs (Sharkey et al. 2025). Overall, there is little to no evidence of CM-type material in the Kuiper Belt based on current surveys, fully in agreement with these simulations strongly suggesting inefficient mixing.

Less common carbonaceous chondrites. These simulations confirm that large (∼100km) planetesimals formed beyond ∼8.6 au experience minimal dynamical disturbance if Saturn is the only currently forming core. This potentially leaves ample time for undisturbed planetesimals to undergo compositional changes or to be later scattered by an additional embryo that emerges in the outer disk, such as our three-planet model, in which case the undisturbed population exists >10.4 au. Such intermediary populations might correspond to the less common carbonaceous chondrite subgroups (CK, CR, CH, CB) or even ungrouped carbonaceous materials, bridging the gap between the CM reservoir and a more distant CI-like region. Although this remains tentative, given the limited number of known large asteroids associated with these types. Current observations do not provide strong support for or against this possibility.

Sequential formation of giant planets. The five-planet case highlights why coeval formation is difficult to reconcile with the observational constraints. When multiple ice giant cores grow in a gas-rich disk, their Type-I migration and repeated encounters can raise the fraction of CM-like bodies temporarily “parked” at 15-20 au to —5-10%. However, this parking is inseparable from the same gas-rich, migratory phase that also promotes contemporaneous delivery from trans-Saturnian reservoirs. In that limit, CM-like and CI-like material are redistributed into the asteroid belt under similar disk conditions and at similar epochs, which would tend to erase the distinct implantation signatures inferred from the present-day radial distributions of CM- and CI-like asteroids.

Avoiding this outcome would require that the CI source region lies well beyond the dynamically accessed zone during the gas-rich phase, implying either a much more distant (≫20 au) CI reservoir or substantial intermediate reservoirs between CM and CI, neither of which is strongly supported by current constraints. These results thus favor a sequential giant-planet formation model: Saturn’s formation and scattering of CM-like bodies likely occurred before Uranus and Neptune reached significant mass, or began migrating. This sequence of events preserves the distinct implantation profiles of CM and CI groups in the asteroid belt. It is also consistent with models where giant planets form one generation at a time in pressure maxima (e.g., Lau et al. 2024), with Saturn’s gap first halting inward flux of solids and only later a new structure emerging beyond it to seed the ice giants.

Arguments against implantation of CI-like bodies. Any interpretation based on implantation efficiencies must also be consistent with simple mass budgets. If the CM source region contained of order MCM ~ 1 M_⊕ of leftover planetesimal mass, then an efficient inward delivery of —5-10% would correspond to ~0.05-0.1 M_⊕ of CM-like material implanted in the asteroid belt. By contrast, adopting a trans-Saturnian CI source region with M_CI ~ 20 M_⊕, a —2-5% delivery efficiency would imply —0.41 M_⊕ of CI-like material. This would yield a CI:CM implanted mass ratio of 10:1 rather than the 1.1:1 observed today.

This not only differs substantially from the near-unity ratio but also far exceeds the current total main-belt mass (—4 × 10⁻⁴, M_⊕ Pitjeva & Pitjev 2018). Some works (e.g., Nesvorný et al. 2024) thus assume additional mechanisms, preferentially destroying CI-like bodies, for instance, planetesimal ablation by gas. However, this mitigation requires sufficient gas to operate, and again reintroduces the problem with explaining differing distributions of CM- versus CI-like bodies.

Small (1-10 km) planetesimals. Although these simulations focus on 100 km-size planetesimals, the conclusion that CM-like material is difficult to implant in the ice giant region is not sensitive to this choice. Raymond & Izidoro (2017b) found that such smaller planetesimals can more readily circularize at moderate distances when scattered outward, but do not necessarily extend their orbits as far as larger bodies. For small (—1-10 km) bodies, their shorter gas-drag timescales allow them to lose orbital energy before they are fully ejected; consequently, many of them remain in stable orbits over ~10 au, though none reach beyond — 15 au. Meanwhile, 10 km objects behave somewhat similarly to 100 km bodies in their ultimate radial extent - often maxing out near —20 au, but experience more effective eccentricity damping than our reference case.

Such objects would be particularly susceptible to reprocessing, collisional grinding, or re-accretion onto later-forming planetary embryos. Future work incorporating explicit size distributions - including pebbles and dust - will be crucial for assessing whether such processes could help explain the compositional diversity of distant objects and the persistence of chondritic materials in orbits beyond Saturn.

Fig. 7 Same as Fig. 4, but for the five-planet models (Anderson et al. 2025). In both simulations, Uranus and Neptune are initialized in the ‘Wide’ configuration. On the left, they are allowed to migrate naturally as they grow due to Type-I torque. On the right, the semimajor axes are fixed. Planetesimals are color-coded based on their formation location. In the simulations in which the ice giant cores are allowed to migrate, surviving planetesimals are circularized quickly.

5 Conclusions

We show that large (100-km) CM-like planetesimals formed at the outer edge of Jupiter’s gap are efficiently scattered by Saturn but only weakly implanted beyond Saturn’s orbit. In two-planet runs, Saturn clears a narrow annulus around its path and sends a few percent (4-6%) inward to the asteroid belt, while bodies that formed beyond ~8.6 au largely retain semimajor axes close to their birth locations, experiencing little scattering or long-range transport on the timescales explored, and could be considered dynamically sheltered. Adding a single ice giant core at —9.6 au modestly extends Saturn’s reach, and while up to 10% are scattered inward, outward implantation of their bodies into long-lived orbits beyond 10 au remains intrinsically inefficient: once scattered to high-e, gas drag damps toward pericenter rather than circularizing at larger radii. Even in gas-rich cases, only up to <4% of 100 km objects achieve stable orbits in the Uranus-Neptune zone.

When multiple ice giant embryos co-form and migrate inward, circularization in the 15-20 au region becomes faster and more common, up to 10%. The embryos are subject to TypeI migration, which continues to reduce their orbits and to steer the planetesimals along to lower q, so even in this favorable case it remains difficult to “park” substantial CM-like mass on wide, long-lived orbits at large distances. Conversely, if extensive CM contamination of the outer disk had occurred at the same epoch as CI formation, subsequent outward migration during the instability would have imprinted similar distributions of CM and CI bodies in the asteroid belt, contrary to observations. A key observational constraint is the lack of CM-like spectra among trans-Neptunian objects in current surveys. Any CM-like population implanted into the 15-20 au region would also represent only a small contaminant (≲0.2%) of the much broader outer reservoir, often taken to span ∼15-25 au, and only a further subset (≲0.05%) would be expected to reach the main belt during later phases, with an upper limit of ≲2 × 10⁻³ M_⊕.

Taken together, these results highlight that attempts to place CM-like material far beyond Saturn’s orbit must contend with either substantial gas availability or additional planetary embryos that can shepherd planetesimals into stable orbits. In this picture, pressure bumps and gaps act as semi-permeable barriers that preserve compositional rings with limited inter-ring mixing. The absence of implantation of CM-like bodies even under gas-dense conditions, together with the distinct CM and CI distributions in the main belt and the lack of CM-like spectra among TNOs (Vernazza et al. 2025), strongly indicates that the CM and CI source reservoirs were spatially and temporally isolated.

Acknowledgements

The work of S.E.A. has been supported by the Institut Origines of Aix Marseille University, and the work of P.V. by CNES, CNRS/INSU/PNP, and the Institut Origines. The work of M.B. has been supported by the Czech Science Foundation (grant number 25-16507S).

References

Appendix A Collisional Efficiencies

A.1 Planetesimal-planetesimal collisions

Our simulations neglect planetesimal-planetesimal collisions, which could, in principle, modify implantation outcomes by fragmenting bodies into smaller sizes that couple more strongly to gas drag. This effect is particularly relevant during phases of giant-planet growth and migration: Carter & Stewart (2020) demonstrated that giant planet migration can drive a substantial number of highvelocity, potentially disruptive collisions. In their study, collisional processing can significantly reshape the population ultimately delivered to the inner Solar System and asteroid belt, even when planets are not strongly evolving.

We estimate mutual collision rates among the outwardly scattered CM-like bodies using a particle-in-a-box approximation. We first identify the subset of survivors satisfying a > 15 au and e < 0.4 at 300 kyr for each simulation s, and convert these fractions f_s into a physical number N_s of bodies by adopting the total CM reservoir mass of M_CM ≈ 1 M_⊕: $$N_s=f_s\frac{M_{\mathrm{C}\mathrm{M}}}{\frac{\pi }{6}\rho D^3}$$(A.1)

The corresponding number density n_s is defined as n_s = N_s/V, where V is the volume of the annulus sampled by these orbits, approximated as a 10 au wide annular torus centered on a = 20 au. The collisional cross section is defined as $\sigma =\frac{\pi }{4}{D}^2$. We take the characteristic relative velocity $v_{\mathrm{r}\mathrm{e}\mathrm{l},\mathrm{s}}=v_K\sqrt{e_s^2+i_s^2}$ where e_s and i_s are the mean values of e and i, respectively, of the surviving particles in each simulation. The resulting collisional timescale is: $$t_{\mathrm{c}\mathrm{o}\mathrm{l}\mathrm{l}}=\frac{n_s}{\sigma v_{\mathrm{r}\mathrm{e}\mathrm{l},\mathrm{s}}}$$(A.2)

For all disk profiles tested, we obtain t_coll > 6 × 10⁹ yr, longer than the age of the Solar System. This equates to a per-object collision probability of P ≃ 1 - exp(−300kyr/t_coll) < 5 × 10⁻³% over 300 kyr for each planetesimal, with the highest probability being in the flat-gas-profile case. Whether the simulations had two or three planets has less of an impact on the collisional probability than the gas profile itself, with sharp-drop off profiles having 0 - 10⁻⁵% probability, sloped profiles having 10⁻⁴% probability, and flat profiles having 10⁻³%. In the specific context of the outer disk, collisions are unlikely to alter our main conclusions: Across all of our runs, no bodies are circularized onto stable orbits at a > 20 au with e < 0.4, and only up to < 4% reach the 15-25 au region on low-eccentricity orbits in our most generous case. In other words, while collisions may strongly affect the size distribution and implantation efficiency inside ∼10 au (as emphasized by Carter & Stewart 2020), the scarcity of long-lived CM-like material at large heliocentric distances in our simulations makes collisional reprocessing an inefficient pathway for enhancing CM-CI mixing in the ice giant zone.

A.2 Collisions with the native disk

Collisions with a pre-existing “native” CI outer-disk population could be more frequent if such a population had already formed. To place an upper bound on this channel, we assumed a planetesimal disk of total mass M_CI = 20M_⊕ between 15 and 25 au (e.g., Hopp et al. 2022; Nesvorný et al. 2024). Using the same approach as above, we change the approximation of $v_{\mathrm{r}\mathrm{e}\mathrm{l}}=v_K\sqrt{e_s^2+e_{\mathrm{C}\mathrm{I}}^2+i_s^2+i_{\mathrm{C}\mathrm{I}}^2}$. We take an eccentricity of e_CI = 0.05 and inclination of i_CI = 3° for the native disk.

In this case, the resulting collisional timescales all lie within a narrow range: we estimate $t_{\mathrm{c}\mathrm{o}\mathrm{l}\mathrm{l},\mathrm{C}\mathrm{M}-\mathrm{C}\mathrm{I}}\sim (1.3-1.8)\times {10}^7$ yr, so over 300kyr, the per-object collision probability is P ≃ 1 - exp(−300kyr/t_coll) ≈ 2%. Such impacts would likely be highly disruptive at the relative velocities inferred here (v_rel = 2 − 3 km.s⁻¹), and therefore would tend to destroy rather than reprocess the material. Taken together, both the small implanted fraction and the long CM-CM collision timescales indicate that outward contamination by large CM-like bodies is dynamically and collisionally inefficient, and would be unlikely to alter the bulk character of a CI-dominated outer reservoir measurably.

We note that Carter & Stewart (2020) also found that outcomes can differ in configurations where Jupiter and Saturn form substantially farther out (e.g., beyond ∼10 au), even when adopting disk conditions comparable to those explored here. Such architectures increase the overlap between the scattering region and the outer disk, potentially raising collision rates. Since our nominal setup places Jupiter and Saturn in a more compact configuration, our results imply that collisional grinding is not a limiting factor for the outer-disk mixing problem: the bottleneck is the dynamical inefficiency of placing large CM-like bodies onto low-eccentricity, long-lived orbits beyond Saturn in the first place.

Appendix B Supplementary figures

Fig. B.1 Semimajor axis of planetesimals at t = 300 kyr as a function of their initial formation zone for each gas profile and Saturn growth timescale for our two-planet model (top) and three-planet model (bottom). The color represents the moment the test particle was lost from the simulation, with cooler colors representing earlier loss and yellow for those that remain. The extent of Saturn’s “reach” is nearly identical across all simulations. The objects outside this reach, the somewhat horizontal yellow bar, will retain their initial semimajor axis, while planetesimals formed around Saturn may “leapfrog” over them to higher heliocentric distances.

All Tables

All Figures

Fig. 1 Cumulative distributions of the semimajor axis N(≤a) of >100 km asteroids classified as CM-like, CI/IDP-like (including B types with the geometric albedo of pv < 0.1), and S-type. Vertical dotted lines mark major mean-motion resonances with Jupiter. CM-like objects exhibit an approximately symmetric distribution, whereas CIlike objects show a notably asymmetric distribution. Adapted from Anderson et al. (2025).

In the text

	Fig. 2 Comparison of the gas profiles used in the various simulations. The initial orbits of Jupiter, Saturn’s core, and the additional ice giant core are indicated by vertical gray lines. Beyond 10 au, the profiles can be seen as Flat (Black, Raymond & Izidoro 2017b), Sloped (Blue, Cresswell & Nelson 2008), and with a Sharp Drop-off (Orange, Desch et al. 2018).
In the text

Fig. 3 Fraction (%) of planetesimals from each initial formation zone that end up within 1-au-wide semimajor-axis bins at t = 300, kyr for the two-planet model (log scale). We consider implanted objects as those with e < 0.4. In all cases, most bodies formed between 9 and 10 au remain near their birth locations, reflecting a dynamically quiet zone largely unaffected by Saturn’s growth. The prominent blue feature near −7 au corresponds to particles trapped in Saturn’s co-orbital region.

In the text

Fig. 4 Eccentricity of remaining planetesimals as a function of semimajor axis at t = 300 kyr for each gas profile and Saturn growth timescale for our two-planet model. The dotted line represents Jupiter-crossing orbits, the dashed line represents Saturn-crossing orbits. The color represents the planetesimals’ formation location, with cooler colors representing earlier loss and yellow for those that remain. The more gas in the simulation, the more easily the planetesimals are circularized (lower final e). The cluster of purple points near ~7 au represents planetesimals trapped near Saturn’s orbit; these bodies are dynamically influenced by the planet and are likely to be absorbed into its circumplanetary disk over time.

In the text

	Fig. 5 Same as Fig. 3, but for the three-planet model. We consider implanted objects as those with e < 0.4. Faster planetary growth clears the local planetesimal population. Objects that formed beyond 10 au are less likely to interact with the forming planetesimals.
In the text

	Fig. 6 Same as Fig. 4, but for the three-planet model. The dotted line represents Jupiter-crossing orbits, and the dashed line represents Saturn-crossing orbits. The color represents the planetesimals’ formation location, with cooler colors representing earlier loss and yellow for those that remain. The more gas in the simulation, the more easily the planetesimals are circularized.
In the text

Fig. 7 Same as Fig. 4, but for the five-planet models (Anderson et al. 2025). In both simulations, Uranus and Neptune are initialized in the ‘Wide’ configuration. On the left, they are allowed to migrate naturally as they grow due to Type-I torque. On the right, the semimajor axes are fixed. Planetesimals are color-coded based on their formation location. In the simulations in which the ice giant cores are allowed to migrate, surviving planetesimals are circularized quickly.

In the text

Fig. B.1 Semimajor axis of planetesimals at t = 300 kyr as a function of their initial formation zone for each gas profile and Saturn growth timescale for our two-planet model (top) and three-planet model (bottom). The color represents the moment the test particle was lost from the simulation, with cooler colors representing earlier loss and yellow for those that remain. The extent of Saturn’s “reach” is nearly identical across all simulations. The objects outside this reach, the somewhat horizontal yellow bar, will retain their initial semimajor axis, while planetesimals formed around Saturn may “leapfrog” over them to higher heliocentric distances.

In the text