Global Modeling of Nebulae With Particle Growth, Drift, and Evaporation Fronts.
III. Redistribution of Refractories and Volatiles

In Figure 2 we first examine the compositional evolution for a compact particle growth model with our fiducial turbulent intensity $\alpha_{\rm{t}}=10^{-3}$, also assuming a gaussian particle collision velocity pdf (case sa3g, see Paper II). Plotted in the top set of panels are the log of the solids mass fractions $x^{\rm{d}}_{i}+x^{\rm{m}}_{i}$ of the various species $i$, at three different times, 0.1, 0.3 and 0.5 Myr. The bottom panels are the corresponding vapor fractions $x^{\rm{v}}_{i}$. Though locations of EFs can be seen in both sets of panels, they are the most easily seen in the vapor. In both sets of panels, the black dotted curves are the total fractions of solids or vapors, respectively, which for the solids also represents the locally varying metallicity $Z$ (Sec. 2.2).

Mass Loss: The first thing to note is the systematic depletion of solid material outside the snowline, that occurs over time due to inward radial drift of material. As a result the inner disk regions become enhanced with solids while the outer disk becomes depleted. The depletion is especially rapid slightly beyond the water ice EF (cyan) because the mass-dominant particles here have the largest Stokes numbers ($\sim 0.02-0.03$) and thus the fastest radial drift (see Fig. 1 panel f, and Paper II, ). Recall that for the models shown, particles and aggregates can grow much larger everywhere beyond the water snowline due to the stickiness of water ice (Sec. 2.2). The particle sizes that correspond to these Stokes numbers are $\sim 15$ cm for the compact particle case, and $\sim 27$ m for the fractal case (see Appendix A, Paper II). As a result, the local $Z$ drops by about an order of magnitude from the initial value in only $3\times 10^{5}$ years. In the inner disk, the solids enhancement appears to peak around 0.3 Myr with $Z\sim 0.03-0.04$ between $3-4$ AU, primarily due to enhancement of solids outside the organics EF (orange). However, because the Stokes numbers of the mass dominant particles are small in the inner disk, with ${\rm{St}}\lesssim 0.001-0.01$, the local solids-to-gas volume density ratio $\epsilon=\rho_{\rm{solids}}/\rho\lesssim 0.1$ so streaming instability is precluded even for this $Z$ (see Paper II, ). In any case, the magnitude of this narrow peak in $Z$ is probably an overestimate because sublimation-condensation of organics is likely not completely reversible as we assume here (see Sec. 4.3, and Paper II for more discussion). On the other hand, the sublimation and condensation of troilite (FeS, green) is likely reversible (Arm et al., 1960), but possibly only down to temperatures of $\sim 500$ K (see, e.g., Pasek et al., 2005). When sublimated, FeS transitions to H₂S and refractory iron in the fraction 56/88 (see e.g., Paper I, ). This explains the sharp increase in refractory iron (red) around $\sim 1$ AU, inside the FeS EF. For this work we implicitly assume condensation back to FeS outside the troilite EF if refractory Fe is available there (however see Sec. 4.3).

Evaporation Fronts: In the outer disk, enhancements in supervolatiles CO (grey), CO₂ (pink) and CH₄ (yellow) can be seen outside their respective EFs, which evolve inwards as the disk cools (mostly due to decreasing stellar luminosity). A large enhancement spike in solid CO is seen at 700AU after 0.1 Myr, while all other solid species’ fractions have systematically decreased in magnitude. This enhancement is not due to inward radial drift of particles, because early on in this region of the disk where the gas density is so low, even monomers have very large Stokes numbers and are effectively immobile. Rather, this enhancement is due to the nebula gas advecting outward as the gas disk spreads. In our models, the initial gas surface density profile drops off sharply outside of $\sim 100$ AU (see Fig. 1, panel a) so that even small changes can make a significant difference in the local solid fractions in the outermost parts of the protoplanetary disk. Even though the nebula gas density remains very low, it has increased by orders of magnitude over the initial value, resulting in decreased instantaneous $Z$ of all species except CO. The CO is enhanced because vapor phase CO is advecting with the nebula gas and condenses outside the CO EF at a rate fast enough to produce the relatively large spike in solid CO fraction. The magnitude of the spike eventually decreases over time once the advected nebula gas begins to increase significantly enough beyond the CO EF as the disk continues to spread outwards. The particle Stokes numbers in this outermost region then decrease, to unity or less, at which point the region can be vacated due to rapid radial drift (e.g., see top left panel of Fig. 6).

The methane and CO₂ EFs lie inside $\sim 200$ AU and thus exhibit a different behavior because the particle Stokes numbers are $\ll 1$. Both EFs have solids enhancements that extend outwards over many AU, with the CH₄ enhancement being especially extended. This is because (like the case for CO) each region is enhanced already by the outward diffusion of sublimated vapor and subsequent recondensation outside their respective EF, as well as the outward advection of the nebula gas, but here the nebula gas can effectively advect small particles with it as well. The peak for CH₄ is broader because the Stokes numbers of the mass-dominant particles are about an order of magnitude smaller than they are at the CO₂ EF, and the gas advection velocity $v_{\rm{g}}$ is increasing sharply in this region and is higher at the methane EF. Thus, small-St particles are advected further at the methane EF.

The bottom panels of Fig. 2 indicate that there are seven EFs still on the computational grid after 0.1 Myr (all 9 are present at $t=0$), so for the times shown, both iron and Ca-Al (“CAIs”, blue) are solid everywhere in the disk. The silicate EF (brown) persists for quite some time, but has evolved inwards off our inner grid boundary at 0.5AU by 0.5 Myr. Because the location where the nebula gas velocity switches from outwards to inwards occurs at $\sim 20$ AU, the enhancements in solids outside the various EFs in the inner disk tend to be due purely to diffusion for these models, while similar EF enhancements outside 20AU are due partly to advection. It is noteworthy that the enhancements in the metallicity of the gas can also become substantial (5-8%). This could systematically increase the local gas molecular weight, affecting the pressure gradient (and viscosity) through both the gas density and the local sound speed (see, e.g. Charnoz et al., 2021). A feedback effect from the increasing molecular weight of the gas might possibly lead to a pressure bump. This becomes more of an issue for lower $\alpha_{\rm{t}}$. We will consider this effect more carefully in future work.

Fractal aggregate vs. compact particle growth: In Figure 3, we show a fractal aggregate growth model (fa3g) for the fiducial $\alpha_{\rm{t}}=10^{-3}$ for comparison with the compact particle growth simulations of Figure 2. Many of the trends described for the compact growth case are also seen here. For example, the sharp decrease in total solids abundance outside the water ice EF ($\sim 6-15$ AU) still occurs, with the solids surface density profiles overall being very similar (Figure 1). Like in model sa3g, the metallicity in the region from the H₂O EF out to $\sim 15$ AU (dotted black lines in top panels of Fig. 2 and 3) has dropped to $Z\gtrsim 0.001$, because the St in this region are similar for sa3g and fa3g ($\sim 0.03$; see Fig. 1, panel f)⁵⁵5It should be understood that growth to large ${\rm{St}}$ and decreasing $Z$ occur concomitantly. That is, the local $Z$ is low because the local ${\rm{St}}$ is large, and $Z$ and ${\rm{St}}$ are never both large enough at the same time and place that conditions for the Streaming Instability can be satisfied (see \al@Est16,Est21; \al@Est16,Est21, and also Umurhan et al. 2020).. The main difference is that now beyond $\sim 20$ AU much more material is retained because the fluffy aggregates have Stokes numbers $5-10$ times smaller than the compact particle case (orange curves, Fig. 1, panel f), though the fractal masses are several orders of magnitude larger than their compact counterparts (see Fig. 2, Paper II). Again, this is because aggregates have not yet compacted enough, so their St remain lower. The aggregates thus have much slower inward radial drifts, and a fraction of them can even be advected outwards. Indeed, where we saw much of the solids being depleted in the region $\gtrsim 40-100$ AU in Fig. 2 leaving a clear enhancement in narrow supervolatile bands of CO₂ and CH₄ solids, the outer edge of this region remains populated with particles in the fractal case and has actually moved outwards from $\sim 200$ to $\sim 300$ AU over 0.5 Myr. The supervolatile bands are still present for the fractal case, but are very much muted (see Sec. 4.1). The differences can also be seen in the $\Sigma_{\rm{solids}}$ profiles (panel b, Fig. 1). Eventually though, we expect even the fractal aggregates outside $\sim 20$ AU will continue to compact and drift into the inner disk. The vapor fractions are similar in both models with the exception that the enhancements inside the supervolatile EFs are somewhat reduced.

“Cold (non-sticky) H₂O ice”: In Figure 4 we show simulations for compact particle (sa3Qg) and fractal aggregate (fa3Qg) growth simulations for $\alpha_{\rm{t}}=10^{-3}$, which differ from the fiducial model in that we here assume that water ice is only “sticky” in a limited temperature range near the snowline (Musiolik & Wurm, 2019). We implement this in our simulations by only using the “sticky ice” $Q_{\rm{f}}=10^{6}$ (see Table 2.4) within 20 K of the water ice condensation temperature, and adopting the lower (silicate) value of $Q_{\rm{f}}$ at lower temperatures. In Fig. 4 we only show the solids fractions for the “cold ice” models at 0.5 Myr. One finds that for both fractal and compact particles, the depletion outside the snowline has been significantly reduced, because particles and aggregates cannot grow to nearly as large St further out in the disk, decreasing the mass influx of material. The region near the H₂O EF in which $Q_{\rm{f}}$ and the fragmentation threshold remain high still has larger particles, with ${\rm{St}}$ similar at their peak to the fiducial case (Fig. 1, panel f), but dropping off by nearly 2 orders of magnitude from there outwards. Enough material drifts in to prevent as large a drop in $Z$ as seen in the fiducial model, but because less material drifts across the water ice EF overall, the enhancements in the inner disk are also smaller. The fractal case fa3Qg shows slightly more retention of material in the outermost portions of the disk than the corresponding compact particle case sa3Qg, but the latter retains much more than the fiducial compact particle growth case sa3g due to the smaller ${\rm{St}}$ (Fig. 1).

Larger $\alpha_{\rm t}$ ($10^{-2}$): Simulations where we increase the turbulence parameter to $\alpha_{\rm{t}}=10^{-2}$ are shown in Figure 5. Because relative velocities between compact or aggregate particles increase with $\alpha_{\rm t}$, the fragmentation sizes and Stokes numbers are considerably smaller (at least initially) so that they remain more well coupled to the nebula gas and less affected by radial drift. As a result, the enhancements of solids and vapors at EFs are more muted overall compared to the fiducial models, and less material is lost from the outer disk. There is essentially no difference seen between compact and fractal growth in the evolution of the vapor species (2nd and 4th rows) over the course of the simulation.

New result on Fe and FeS: However, some interesting behavior is seen for the particles in the inner disk. In the compact growth model (1st row), the condensed FeS outside its EF appears to have the strongest enhancement - similar to or even larger than seen at the organics and water EFs. This FeS-enrichment comes at the expense of a dramatic depletion in iron metal near $0.9$ AU. This interesting effect, which actually begins around 0.2 Myr (not shown) in both fractal and compact growth models, occurs because the quick diffusion of a substantial amount of H₂S across the FeS EF is reacting with inward drifting, native Fe metal to produce a large spike in condensed FeS, but the Fe metal is not replaced as quickly by inward radial drift of material from further out. By 0.5 Myr and beyond though, enough Fe dust is drifting in to begin to erase the Fe-depletion signature. This localized Fe-depletion is not seen in the $\alpha_{\rm t}=10^{-3}$ fiducial model where the enhancement in FeS is absent, and instead the refractory iron is enhanced interior to the troilite EF. This Fe-depletion disappears sooner in the fractal aggregate growth model (mostly absent by 0.3 Myr) for this large $\alpha_{\rm t}$ (3rd row of Fig. 5) because, due to slower radial drift (smaller St in the inner disk overall), the enhancement of H₂S, and its feedback across the EF, does not build up sufficiently to cause a depletion in refractory iron for nearly as long. We discuss this further in section 4.1.

Global evolution with time for vapor and solids: As pointed out above, there is no notable difference between the compact and fractal cases in the vapor evolution for the $\alpha_{\rm t}=10^{-2}$ models, and at 0.1 Myr, the solids fractions in the outer disk regions beyond the snowline (1st and 3rd rows) also look quite similar. This changes gradually at later times. For both fractal and compact cases as the gas disk spreads outwards and the gas surface density increases significantly beyond $\sim 100$ AU (see panel a, Fig. 1), particles and aggregates are more easily advected outwards with the nebula gas, whose velocity $v_{\rm{g}}$ increases strongly beyond $\sim 20$ AU. Outside $\sim 90$ AU, the ${\rm{St}}$ of the fractal aggregates are smaller than their compact counterparts and are thus more coupled to the gas, so they are transported outward in greater abundance. This explains the slight difference in $Z$ between the models at large distances at later times.

The rapid gas evolution for $\alpha_{\rm t}=10^{-2}$ also decreases the gas density everywhere inside of $\sim 100$ AU significantly, and ${\rm{St}}$ generally increases (red curves in Fig. 1, panel f) to values comparable to or even larger than the lower $\alpha_{\rm{t}}$ models despite being smaller in mass in comparison (cf. Figs. 2 and 8, Paper II, ). The differences in particle mass between the compact and fractal cases is still a few orders of magnitude (Fig. 8, Paper II, ), but the Stokes numbers of the mass dominant fractal aggregates are roughly the same for the compact and fractal growth cases, between the snowline and 90AU. We note that the much higher level of coupling between particles and the gas in both cases more or less maintains $Z$ close to its initial value across the disk. This is largely true even in the inner disk, though there are slight enhancements in the solids outside EFs out to the snowline (including even the silicates at 0.1 Myr) in the compact particle growth simulation. Overall it appears the distinction between compact particle and fractal aggregate models diminishes for increasing $\alpha_{\rm{t}}$.

Smaller $\alpha_{\rm t}$: The change in overall behavior for smaller $\alpha_{\rm{t}}$ is much more dramatic, especially for compact particle growth, as can be seen for the case of $\alpha_{\rm{t}}=10^{-4}$ plotted in Figure 6. The top panel (1st row), for the compact particle growth case, shows that the outer disk is much more quickly drained of solids to the inner disk regions compared to the fiducial $\alpha_{\rm{t}}=10^{-3}$ model for compact particles, and also that material within the inner disk continues to be lost through the inner boundary and presumably to the star. The ${\rm{St}}$ values of the mass dominant particles for both fractal and compact growth models are still relatively small and similar to each other (Fig. 1), but now an order of magnitude larger than in the fiducial model (so is the gas surface density). Similar St between sa4g and fa4g also means the loss rate of particles is roughly the same for the compact and fractal aggregate growth cases (1st and 3rd rows) inside the snowline despite the aggregates being 2 orders of magnitude larger in mass than the compact ones (see Fig. 8 Paper II, ). The exception is the extended knee inwards of the snowline from $\sim 2-3$ AU seen in Fig. 1 (panel f). We note that this region inside the snowline corresponds to a sharp drop in the normalized pressure gradient, or headwind parameter $\beta$ (panel d), which would slow the radial drift rate in the region inside the snowline. For the fractal aggregates, the temperature (panel c) is hotter and has a steeper gradient in this region (which corresponds to the organics EF). Coupled with the enhancement in organics, this drop in the headwind parameter allows aggregates that decreased in mass (due to loss of water) upon crossing the snowline to grow significantly larger again (compared to the compact particle case, see black curves Fig. 1, panel f) before drifting to the organics EF where considerable material is evaporated, so that aggregate Stokes numbers and particle masses plummet. The spike in water ice outside the snow line is also lower, meaning that aggregates that cross the snowline are less ice rich compared to the compact particle case, so they lose less mass as they cross this EF. The same effect is not seen in model sa4g because the dip in $\beta$ is further inward and not as large in magnitude.

Both fractal and compact particle growth cases see large enhancements in vapors inside, and in solids outside, the organics EF giving local (solids) $Z$ as high as $\sim 0.03-0.04$. The overall $Z$ curves inside the snow line look very similar in magnitude for both models, but differences in their $\Sigma_{\rm{solids}}$ are notable (Fig. 1, panel b). The silicates EF also contributes to the local enhancement of $Z$, and persists longer in the fractal aggregate case. That is, the compact growth model cools much more quickly (see below). The enhancements in vapor (2nd and 4th rows) are considerably larger than in the fiducial model after 0.5 Myr, reaching values as high as $\sim 0.1$ and helping to produce the large enhancement in solid organics exterior to its EF. However, as noted previously, the solids fraction just outside the organics EF may not be realistic because the sublimated or decomposed CHON organics might break down into more volatile constituents such as CO or CH₄, rather than recondense as a CHON. This would not greatly affect the enhancement to the total vapor fraction, which gets quite large in this region at least until these constituents have time to diffuse. As mentioned in Sec. 3.2.1, such large values for the vapor fraction further indicate that consideration of changes in the gas’ molecular weight may become important.

Fe-FeS behavior: The depletion of refractory iron over a small radial range is even more dramatic in this $\alpha_{\rm{t}}=10^{-4}$ case than it was in the case for $\alpha_{\rm{t}}=10^{-2}$, and unlike the latter case, the depletion is not temporary. For both the compact and fractal growth cases (1st and 3rd rows), the depletion of refractory iron begins very early ($<0.1$ Myr) and the radial region over which this depletion grows increases with time as the excess diffusing H₂S gas gobbles up available Fe, maximizing the amount of troilite that can be sublimated again at the EF. The feedback creates so much excess H₂S that it diffuses further outward until it can react with more refractory iron further out in the disk. This process thus creates a situation in which troilite and H₂S coexist out to as far as $\sim 2$ AU in the compact growth case, and somewhat less in the fractal case (2nd and 4th rows), though with time we would expect the situation to extend further out, as it did in the case with $\alpha_{\rm t}=10^{-2}$. Recall this effect is not seen in the fiducial model. Although it is not clear why the $\alpha_{\rm{t}}=10^{-3}$ model does not show this effect, we suspect that it is because the inner disk remains hotter, longer than the lower and higher $\alpha_{\rm{t}}$ models (see Fig. 1, panel c). Despite the more rapidly drifting compact and aggregate particles (due to their higher St overall), Fe still cannot drift in fast enough to overcome this effect; presumably the effect could continue until most of the S from a range of radii is utilized (for more discussion see section 4.1).

Much like we saw in the fiducial case, the region outside the snowline is significantly depleted in solids, as growth there is rapid due to the higher threshold $Q_{\rm{f}}$ for ice. However, depletion in this region is especially rapid and substantial in the fractal aggregate growth model (3rd row). Even after 0.1 Myr the region between $\sim 4-15$ AU sees a large drop in the local metallicity that is absent in the compact particle case at this juncture. The more rapid growth of fractal aggregates comes about in large part due to their larger cross sections (Paper II). As aforementioned, particle Stokes numbers inside the snowline are similar between the compact and fractal models (with the exception of the knee region from $\sim 2-3$ AU), but outside the water ice EF fractal growth has led to mass-dominant aggregates with ${\rm{St}}$ about a factor of 5 larger than the compact particles. In fact, by this time, the aggregate masses are already as much as 6 orders of magnitude more massive than the compact particles in the same region (see Fig. 9 of Paper II, ). This leads to initially more rapid drift of solids into the inner disk regions and explains why, at the early time of 0.1 Myr, there is more enhancement in the inner disk for the fractal case. But the rapid depletion means that by later times the net inward flux of material is lower than in the compact case which steadily drains into the inner disk. This is because in the fractal case, the mass-dominant aggregate ${\rm{St}}$ drops below that of the compact case outside of $\sim 15-20$ AU by a factor $\sim 5-10$ (Fig. 1, panel f) so their drift rates are much slower, and considerable material remains there. On the other hand, rapid, continuous drift in model sa4g accounts for the much larger enhancement spike at the snowline in that model. By 0.5 Myr, the compact growth case has effectively lost most all of its material outside the snowline, while the fractal case retains a significant fraction.

Supervolatile bands? An effect of strong radial drift of disk material in the compact growth model case is the “trapping” of material that manifests as bands of supervolatiles centered at $\sim 55$AU (CO₂), $\sim 150$AU (CH₄) and beyond $\sim 400$AU (CO). The magnitudes of these spikes in $Z$ decrease over time, but for different reasons, as discussed for the fiducial model. We noted before that the region between $\sim 200-400$ AU is quickly vacated⁶⁶6This specific radial range is not special but will vary with the choice of $R_{0}$. A larger (smaller) value of $R_{0}$ would move this location outwards (inwards). Indeed, a smaller $R_{0}$ may better correspond to an early stage disk subject to infall., even for small particles, because the Stokes numbers are already of order unity or larger there, so radial drift is rapid. On the other hand, outside $\sim 400$ AU, Stokes numbers become increasingly large due to the sharply decreasing gas density so the particles and aggregates are not drifting much at all. The decrease in mass fraction there is solely due to the slow expansion of the nebula gas disk over time, which in this low $\alpha_{\rm{t}}$ case is too slow to diminish the CO peak after 0.5 Myr compared to the fiducial model. Also as before, the band of CH₄ is seen to be more extended in the fractal model due to the slower drift rates of the porous aggregates compared to the compact particles. It should be noted that these bands do not have a lot of mass, at least at this epoch, as can be gleaned from panel (b) of Fig. 1. For the compact growth model where a clear CO₂ band is seen, there is a peak surface density of $\sim 0.01$ g cm^-2, while for methane it is a hundred times smaller. The fractal model has a similar magnitude for CH₄, but at the CO₂ EF most species (Fig. 6, 3rd row and panel) are still present thanks to the low radial drift speeds of the porous, fluffy aggregates there (see also Fig. 8) so the fraction of CO₂ is 2 orders of magnitude larger.

Temperature variations with $\alpha_{\rm t}$: Amongst the above models, previously we mentioned that the fiducial model ($\alpha_{\rm t}=10^{-3}$) remains the hottest after 0.5 Myr (Fig. 1, panel c), and at earlier times there is variation in the cooling rates of the disk amongst the different models. This is not difficult to understand. The $\alpha_{\rm{t}}=10^{-2}$ case is cooler than the fiducial model despite having smaller particles and aggregates, which increases the opacity (Fig. 1, panel e), simply due to the order of magnitude decrease in gas density by 0.5 Myr compared to $\alpha_{\rm{t}}=10^{-3}$ (Fig. 1, panel a). The compact particle growth case with $\alpha_{\rm{t}}=10^{-4}$ is only marginally hotter than the corresponding high turbulence case because even though there is so much more nebula gas (and much more solids mass) for $\alpha_{\rm{t}}=10^{-4}$, the opacity has decreased significantly due to rapid particle growth. The fractal case for the low turbulence value $\alpha_{\rm{t}}=10^{-4}$ stays hotter longer than the compact growth case for this $\alpha_{\rm{t}}$, as seen in Fig. 6 (4th row) where the silicate EF persists up to 0.3 Myr, whereas it has already evolved inward of our computational inner boundary in the compact particle growth case by 0.1 Myr.

Extremely low $\alpha_{\rm t}$ case: Finally, in Figure 7 we present models for $\alpha_{\rm{t}}=10^{-5}$. The compact (sa5g) and fractal (fa5g) cases provide an even starker contrast to the $\alpha_{\rm t}=10^{-3}-10^{-2}$ models than did $\alpha_{\rm{t}}=10^{-4}$. For the compact particle growth case, such a low turbulent intensity leads to rapid loss of most solid material in a considerably shorter period of time than for $\alpha_{\rm{t}}=10^{-4}$. The end result is that after 0.3 Myr essentially only the banded EF-related structures remain, which only fade slightly after 0.5 Myr. The spike in H₂O ice at its EF actually increases slightly, due to condensation from a large reservoir of vapor phase water as the disk cools and the EF moves inwards. The fractal aggregate case fa5g proceeds in a similar fashion to fa4g, with rapid growth outside the snowline to aggregate masses about an order of magnitude larger even than the $\alpha_{\rm{t}}=10^{-4}$ model ($\gtrsim 10^{9}$ g, see Fig. 10, Paper II), and correspondingly larger ${\rm{St}}$. The larger ${\rm{St}}$ quickly lowers the local $Z$ to values $\lesssim 0.001$. For this $\alpha_{\rm t}$, advection and diffusion by the nebula gas are of course very small. It is clear that even for this low $\alpha_{\rm{t}}$, the solid material in the fractal growth case persists longer than for the compact particle growth case, but unlike model fa4g, much of the mass in fa5g is gone by 0.5 Myr. However, we do find for this lowest $\alpha_{\rm t}$ that there is a short period of time at $\sim 0.2$ Myr where the solids-to-gas ratio just outside the snowline in both growth models is greater than unity and, with both the mass-dominant compact and aggregate particles having large ${\rm{St}}$, may briefly satisfy conditions for the streaming instability (see Umurhan et al., 2020) and a burst of ice-rich planetesimal formation. This is discussed in detail in Section 4.1 of Paper II.