Making Terrestrial Planets: High Temperatures, FU Orionis Outbursts, Earth, and Planetary System Architectures

The dynamics of dust grains in a protoplanetary disk are controlled by their aerodynamical drag:

$$\begin{eqnarray}&&{\partial }_{t}{\boldsymbol{v}}=-\displaystyle \frac{{\boldsymbol{v}}-{\boldsymbol{u}}}{\tau },\end{eqnarray} \tag{ 1 }$$

where ${\boldsymbol{v}}$ is the dust grain's velocity, ${\boldsymbol{u}}$ the gas velocity at the dust grain's location, and τ the drag based stopping time. It is conventional to non-dimensionalize τ through the local Keplerian orbital frequency Ω_K to define a Stokes number:

$$\begin{eqnarray}&&\mathrm{St}\equiv \tau {{\rm{\Omega }}}_{K}.\end{eqnarray} \tag{ 2 }$$

In regions associated with the formation of the Earth, the dust grains we will consider are small enough to be in the Epstein drag regime, with midplane Stokes numbers

$$\begin{eqnarray}&&{\mathrm{St}}_{0}=\displaystyle \frac{\pi a{\rho }_{\bullet }}{2{{\rm{\Sigma }}}_{g}},\end{eqnarray} \tag{ 3 }$$

where a and ${\rho }_{\bullet }$ are the dust grain radius and solid density, and Σ_g the disk's gas surface density.

For the inner disk, we will also consider solids larger than the gas mean-free-path, in the Stokes drag regime, which depends on the gas mean molecular mass rather than the gas density. Such grains have Stokes numbers

$$\begin{eqnarray}&&{\mathrm{St}}_{0}=\displaystyle \frac{2}{9}\displaystyle \frac{{\rho }_{\bullet }{a}^{2}}{{m}_{g}}\displaystyle \frac{{\sigma }_{g}}{H},\end{eqnarray} \tag{ 4 }$$

where m_g ≃ 2.3 amu and σ_g ≃ 2 × 10⁻¹⁵ (Chapman & Cowling 1970) are the gas mean molecular mass and molecular collisional cross-section, respectively. Further,

$$\begin{eqnarray}&&{v}_{\mathrm{th}}=\langle v\rangle =\sqrt{\displaystyle \frac{8{k}_{B}T}{\pi {m}_{g}}}\end{eqnarray} \tag{ 5 }$$

is the local gas thermal speed and

$$\begin{eqnarray}&&H=\sqrt{\displaystyle \frac{\pi }{8}}\displaystyle \frac{{v}_{\mathrm{th}}}{{{\rm{\Omega }}}_{K}}\end{eqnarray} \tag{ 6 }$$

defines the local vertically isothermal pressure scale height.

Due to their radial pressure gradient, disks orbit slower than Keplerian, with a velocity deficit δv. If the pressure is a power law in radius with power s, δv is given by

$$\begin{eqnarray}&&\delta v=\displaystyle \frac{\pi }{16}s\displaystyle \frac{{v}_{\mathrm{th}}^{2}}{R}{{\rm{\Omega }}}_{K}^{-1}\simeq 50\,{\rm{m}}\,{{\rm{s}}}^{-1}.\end{eqnarray} \tag{ 7 }$$

where s ≃ 13/4 for a Hayashi MMSN (Hayashi 1981). Note that δv is independent of R in a disk with T ∝ R^−1/2. A dust grain with St ≪ 1 drifts inward at a speed of

$$\begin{eqnarray}&&{v}_{r}\simeq 2\,\mathrm{St}\,\delta v\simeq \mathrm{St}\times {10}^{4}\,\mathrm{cm}\,{{\rm{s}}}^{-1},\end{eqnarray} \tag{ 8 }$$

allowing us to define a drift timescale,

$$\begin{eqnarray}&&{\tau }_{R}=\displaystyle \frac{R}{{v}_{r}}\simeq \displaystyle \frac{R}{2{St}\delta v}=\displaystyle \frac{1}{2\pi s\,\mathrm{St}}\displaystyle \frac{{R}^{2}}{{H}^{2}}\mathrm{Orb},\end{eqnarray} \tag{ 9 }$$

where Orb is the local orbital period. Local structures such as pressure bumps alter δv, potentially reversing its sign, allowing the trapping of particles (Dittrich et al. 2013).

The streaming instability (SI) is a mechanism that produces planetesimals from the minimal constituent size, significant because dust coagulation theory does not predict large dust grains. The trigger conditions for the SI are nonetheless relatively large dust grains of ${\mathrm{St}}_{2}\equiv \mathrm{St}/0.01\simeq 1$ combined with elevated local dust-to-gas mass ratios of $\tilde{\epsilon }\equiv \epsilon /3\times {10}^{-2}\simeq 1$ (Carrera et al. 2015). While that value is above the expected overall value of  = 5 × 10⁻³ for the Solar Nebula inside of the frost line (Lodders 2003), there are several ways to effectively enhance , including strong settling, or radial pressure or temperature traps (Dittrich et al. 2013; Hubbard 2016b).

Bouncing sets the upper size limit for locally coagulating reasonably compact (non-fractal) dust grains. We assume that the critical bouncing velocity v_b is significantly below the fragmentation velocity, allowing us to neglect the latter. Fragmentation velocities only modestly above the bouncing velocity would further reduce the upper size limit. The bouncing barrier and its consequences have been explored numerically (e.g., Zsom et al. 2010; Windmark et al. 2012; Garaud et al. 2013; Estrada et al. 2016) and are found to prohibit grain growth to SI triggering dust sizes.

From Güttler et al. (2010), the critical bouncing velocity for porous dry silicates is

$$\begin{eqnarray}{v}_{b} & \simeq & f\times {10}^{-1}{\left(\displaystyle \frac{m}{{10}^{-4}{\rm{g}}}\right)}^{-5/18}\,\mathrm{cm}\,{{\rm{s}}}^{-1},\\ & \simeq & f\times 3\times {10}^{-4}\,{\tilde{R}}^{5/4}{\mathrm{St}}_{2}^{-5/6}\,\mathrm{cm}\,{{\rm{s}}}^{-1},\end{eqnarray} \tag{ 10 }$$

where m is the dust grain mass, $\tilde{R}$ the orbital position in au, and f an arbitrary scaling factor which parameterizes the dust's surface parameter uncertainties. While strict sticking/bouncing velocity cut-offs such as Equation (10) are simplifications, we can estimate the rate of sticking events as the rate of collisions at velocities v < v_b. In Equation (10), following Güttler et al. (2010), the grains have volume filling factors of ϕ = 0.12 and monomer densities of ρ = 2 g cm⁻³, slightly at odds with the canonical molten solid density we will adopt later of ${\rho }_{\bullet }=3$ g cm⁻³ (Friedrich et al. 2015).

Assuming turbulent stirring, we can estimate the collision velocity scale as (Voelk et al. 1980)

$$\begin{eqnarray}&&{v}_{0}\simeq \sqrt{\alpha \mathrm{St}}\,{c}_{s}=38\sqrt{{\alpha }_{5}{\mathrm{St}}_{2}}\,{\tilde{R}}^{-1/4}\,\ \mathrm{cm}\ {{\rm{s}}}^{-1}.\end{eqnarray} \tag{ 11 }$$

In this section, to give dust coagulation the best chance, we assume that the disk is locally barely turbulent, approximating an α disk with ${\alpha }_{5}\equiv \alpha /{10}^{-5}\simeq 1$, appropriate for quiescent dead zone midplanes (Oishi & Mac Low 2009). In the limit of v_b ≪ v₀, we can estimate the rate at which a given grain collides with like grains at velocities v < v_b to be (Hubbard 2016a)

$$\begin{eqnarray}&&S\simeq \sqrt{2\pi }{{na}}^{2}{v}_{0}{\left(\displaystyle \frac{{v}_{b}}{{v}_{0}}\right)}^{4},\end{eqnarray} \tag{ 12 }$$

where the number density of the grains is

$$\begin{eqnarray}&&n=\displaystyle \frac{\epsilon {\rho }_{g}}{m}.\end{eqnarray} \tag{ 13 }$$

Equating the drift rate v_r/R from Equation (9) with the estimated sticking rate S from Equation (12) gives an extremely optimistic estimate for the maximum size dust grains can reach before radially drifting out of our zone of interest. That becomes

$$\begin{eqnarray}&&\sqrt{2\pi }{{na}}^{2}{v}_{0}{\left(\displaystyle \frac{{v}_{b}}{{v}_{0}}\right)}^{4}=\displaystyle \frac{{v}_{r}}{R},\end{eqnarray} \tag{ 14 }$$

which reduces to

$$\begin{eqnarray}&&f\simeq 6.6\times {10}^{4}\times {\mathrm{St}}_{2}^{41/24}{\alpha }_{5}^{3/8}{\tilde{\epsilon }}^{-1/4}{\tilde{R}}^{-11/8}.\end{eqnarray} \tag{ 15 }$$

Equation (15) implies that our estimate of the critical bouncing velocity (Equation (10)) would need to be low by four to five orders of magnitude to allow the SI (i.e., ${\mathrm{St}}_{2}\simeq 1$). Thus, dust coagulation theory predicts that dry planetesimals will not form.

That dust coagulation theory prediction is at odds with evidence from the solar system. One obvious remedy would be poor estimates of the sticking parameters of dry silicates (Kimura et al. 2015) or of the appropriate monomer size (Arakawa & Nakamoto 2016), but Equation (15) makes it clear, however, how extreme those mis-estimates would need to be. Magnetic interactions could also increase the sticking rates of dust grains, but Hubbard (2016a) showed that the critical magnetic velocity is expected to be well below the collisional velocity v₀ from Equation (11) for St ≳ 10⁻², so magnetic interactions should die off well before the SI can be triggered.

Fischer–Tropsch-like processes could have led to the synthesis of complex organic molecules within the Solar Nebula (Hill & Nuth 2000; Flynn et al. 2013), depositing layers of sticky organic solids on the dust (Kuga et al. 2015). Those layers would have dramatically increased the dust's stickiness. However, the Earth is not merely dry, but also strongly depleted in carbon. While the degree of that depletion is uncertain, it seems likely that adding more than just a few mass percent of chondritic material would oversupply the Earth's entire carbon budget (Marty et al. 2013), to say nothing of even more strongly carbon rich organic gunk. Thus, any organic gunk model would need to demonstrate that it can operate without oversuppling carbon. Indeed, if recent laboratory work on catalyzing Fischer–Tropsch like processes is confirmed (Nuth et al. 2016), models for the formation of the Earth will likely need to explain the destruction of such layers!

Lucky grains, which undergo a series of low probability low velocity sticking collisions, can grow significantly larger than their more abundant staid counterparts, but numerical studies have found that lucky grains are too rare to trigger planetesimal formation (Estrada et al. 2016). That also rules out importing large dust grains from outside the frost line. Importing sufficient numbers of fully formed planetesimals from outside the frost line would allow planet formation to proceed, but those planetesimals would need to accrete sufficient local dry material to drop their bulk volatile abundances.

Thanks to the radial pressure gradient, planetesimals move with a speed of about δv = 50 m s⁻¹ with respect to the local gas and small-scale dust (Equation (7)). Assuming perfect, purely geometric sweep-up of ambient dust, those planetesimals would grow at a rate of about

$$\begin{eqnarray}&&{\partial }_{t}a=\displaystyle \frac{\epsilon {\rho }_{g}}{4{\rho }_{\bullet }}\delta v\simeq \displaystyle \frac{3}{8}{\tilde{R}}^{-11/4}\ \mathrm{cm}\ {\mathrm{yr}}^{-1}.\end{eqnarray} \tag{ 16 }$$

Thus, purely geometrically accreting planetesimals will only gain a few km of dry surface over a Myr, which seems unlikely to sufficiently overwhelm their volatile rich core. Pebble accretion could speed the process as long as the available grains are comparable in size to chondrules (Johansen et al. 2015), but pebble accretion's size sensitivity poses its own difficulties.

Highly viscous molten grains collide as solids rather than as liquids. As long as those effectively solid grains are not extremely sticky (f ∼ 6.6 × 10⁴, Equation (15)), the bouncing barrier remains in place at high viscosity. If the grains are instead barely viscous, and too large for surface tension to contain the collisional energy, they will splash upon collision. The viscosity condition for colliding as liquids is (Ciesla et al. 2004; Hubbard 2015)

$$\begin{eqnarray}&&\eta \lesssim 6\times {10}^{7}\left(\displaystyle \frac{a}{\mathrm{cm}}\right){\left(\displaystyle \frac{v}{\mathrm{cm}{{\rm{s}}}^{-1}}\right)}^{-1/5}\ {\rm{P}},\end{eqnarray} \tag{ 17 }$$

where P stands for Poise, the cgs unit of dynamical viscosity.

While the condition for avoiding splashing is not certain, from Jacquet & Thompson (2014) we have

$$\begin{eqnarray}&&\eta \gtrsim 4\times {10}^{-12}{\left(\displaystyle \frac{a}{\mathrm{cm}}\right)}^{3}{\left(\displaystyle \frac{v}{\mathrm{cm}{{\rm{s}}}^{-1}}\right)}^{5}\ {\rm{P}},\end{eqnarray} \tag{ 18 }$$

where we have assumed a surface tension of 400 dyn cm⁻¹ (0.4 N m⁻¹) for molten chondritic material (Susa & Nakamoto 2002). Combining Equations (17) and (18) we find the condition

$$\begin{eqnarray}&&\left(\displaystyle \frac{a}{\mathrm{cm}}\right){\left(\displaystyle \frac{v}{\mathrm{cm}{{\rm{s}}}^{-1}}\right)}^{2.6}\lt 4\times {10}^{9}\end{eqnarray} \tag{ 19 }$$

for there existing a viscosity for which liquid, non-splashing collisions are possible.

At temperatures near T ≃ 10³ K, the thermal speed is approximately c_s ≃ 2 × 10⁵ cm s⁻¹ and the turbulent collision speed is about (Equation (11))

$$\begin{eqnarray}&&v\simeq 2\sqrt{\mathrm{St}}\times {10}^{4}\ \mathrm{cm}\ {{\rm{s}}}^{-1},\end{eqnarray} \tag{ 20 }$$

where we have assumed α ∼ 10⁻². In a Hayashi MMSN (Hayashi 1981), with ${{\rm{\Sigma }}}_{g}=1700\,{\tilde{R}}^{-1.5}$, assuming T ≃ 10³ K, we can combine Equations (3), (4), (19), and (20) to find the largest possible grain sizes:

$$\begin{eqnarray}&&{a}_{E}\lesssim 5.7\,{\tilde{R}}^{-0.85}\ \mathrm{cm},\end{eqnarray} \tag{ 21 }$$

$$\begin{eqnarray}&&{a}_{S}\lesssim 6.3\,{\tilde{R}}^{0.54}\ \mathrm{cm},\end{eqnarray} \tag{ 22 }$$

for the Epstein (dominant outwards of 0.93 au) and Stokes (dominant inwards of 0.93 au) drag regimes, respectively.

Equations (21) and (22) imply largest possible Stokes numbers of

$$\begin{eqnarray}&&{\mathrm{St}}_{E}\lesssim 0.016\,{\tilde{R}}^{0.65},\end{eqnarray} \tag{ 23 }$$

$$\begin{eqnarray}&&{\mathrm{St}}_{S}\lesssim 0.015\,{\tilde{R}}^{-0.42},\end{eqnarray} \tag{ 24 }$$

respectively. Those sizes require viscosities of

$$\begin{eqnarray}&&{\eta }_{E}\simeq 7.7\times {10}^{7}{\tilde{R}}^{-0.9}\ {\rm{P}},\end{eqnarray} \tag{ 25 }$$

$$\begin{eqnarray}&&{\eta }_{S}\simeq 8.5\times {10}^{7}{\tilde{R}}^{0.6}\ {\rm{P}},\end{eqnarray} \tag{ 26 }$$

to reach. Thus, the liquid bouncing/fragmentation barriers prevent growth to St ≫ 10⁻², with critical viscosities to reach St ≃ 10⁻² of η ∼ 5–10 × 10⁷ P. The precise viscosities of Solar Nebula solids are unknown, but those values would have required temperatures near T = 10³ K (Hubbard 2015).

Being molten is not a panacea for the formation of terrestrial planets, and the collision of heated solids does not allow direct coagulation to St ≫ 10⁻². However, it does seem likely that temperatures in a relatively narrow range near T ∼ 10³ K permit nearly perfect sticking up to St ∼ 10⁻², conspicuously the critical value for triggering the SI. Unfortunately, we lack detailed laboratory experiments measuring the effective viscosity of chondritic melts of centimeter grains at the relevant temperatures, with thermal histories, and we lack detailed zero-g, low pressure splashing experiments. That limits how thoroughly molten grain growth can be modeled beyond the simple estimates above.

Neglecting bouncing and fragmentation, the growth rate of a given dust grain is the same as sweep-up up to the details of the velocity (Equation (16)):

$$\begin{eqnarray}&&\displaystyle \frac{\partial a}{\partial t}=\displaystyle \frac{\epsilon {\rho }_{g}v}{4{\rho }_{\bullet }}.\end{eqnarray} \tag{ 27 }$$

From Equations (3) and (4), we can see that the Stokes number's quadratic dependence on a means that the Stokes number increases more rapidly in the Stokes drag regime, and we need only consider the growth timescales in the Epstein regime where

$$\begin{eqnarray}&&{St}=\sqrt{\displaystyle \frac{\pi }{8}}\displaystyle \frac{a{\rho }_{\bullet }{{\rm{\Omega }}}_{K}}{{\rho }_{g}{v}_{\mathrm{th}}}.\end{eqnarray} \tag{ 28 }$$

gives the Stokes number at all altitudes.

Combining Equations (27) and (28), we find

$$\begin{eqnarray}&&\displaystyle \frac{\partial \mathrm{St}}{\partial t}=\sqrt{\displaystyle \frac{{\pi }^{3}}{32}}\epsilon \left(\displaystyle \frac{v}{{v}_{\mathrm{th}}}\right){\mathrm{Orb}}^{-1}\simeq \epsilon \left(\displaystyle \frac{v}{{v}_{\mathrm{th}}}\right){\mathrm{Orb}}^{-1}.\end{eqnarray} \tag{ 29 }$$

We can insert Equation (11) into (29) to find

$$\begin{eqnarray}&&2\sqrt{\mathrm{St}(t)}=\sqrt{\alpha }\epsilon \displaystyle \frac{t}{\mathrm{Orb}}+2\sqrt{\mathrm{St}(0)}.\end{eqnarray} \tag{ 30 }$$

Neglecting St(0), we arrive at

$$\begin{eqnarray}&&\mathrm{St}(t)\simeq \displaystyle \frac{\alpha {\epsilon }^{2}}{4}{\left(\displaystyle \frac{t}{\mathrm{Orb}}\right)}^{2}=2.25\times {10}^{-6}\,{\alpha }_{2}{\tilde{\epsilon }}^{2}{\left(\displaystyle \frac{t}{\mathrm{Orb}}\right)}^{2},\end{eqnarray} \tag{ 31 }$$

or alternatively, the time to reach a given St is

$$\begin{eqnarray}&&t(\mathrm{St})\simeq \displaystyle \frac{2}{\epsilon }\sqrt{\displaystyle \frac{\mathrm{St}}{\alpha }}\mathrm{Orb}\simeq 67\displaystyle \frac{1}{\tilde{\epsilon }}\sqrt{\displaystyle \frac{{\mathrm{St}}_{2}}{{\alpha }_{2}}}\mathrm{Orb},\end{eqnarray} \tag{ 32 }$$

where we recall that α and Stare both normalized to 10⁻²: ${\mathrm{St}}_{2}=\mathrm{St}/{10}^{-2}$ and ${\alpha }_{2}=\alpha /{10}^{-2}.$

Comparing Equation (9) to (32), we find that dust growth outpaces radial drift as long as

$$\begin{eqnarray}&&\mathrm{St}\lt {\left[\displaystyle \frac{\alpha {\epsilon }^{2}}{16{\pi }^{2}{s}^{2}}\displaystyle \frac{{R}^{4}}{{H}^{4}}\right]}^{1/3}\simeq 0.16\,{\alpha }_{2}^{1/3}{\tilde{\epsilon }}^{2/3}\,{\tilde{R}}^{-1/3}.\end{eqnarray} \tag{ 33 }$$

Dust growth can outpace radial drift up to SI triggering dust grain sizes for reasonable dust-to-gas mass ratios even for modest turbulence and relatively thick disks.

Boley et al. (2014) proposed the inner disk, with T ∼ 1500 K, as a region where partially molten grains could coagulate well beyond the conventional bouncing and fragmentation limits, envisioning coagulation up to planetesimal sizes. As we showed in Section 3.1 though, new versions of the bouncing and fragmentation barriers show up as dust grains grow beyond about St ≃ 0.01. However, that approximate Stokes number of St = 0.01 is sufficiently large to allow the SI (Johansen et al. 2007) to trigger, albeit in a narrow radial annulus where the temperature is modestly above T = 10³ K and the corresponding viscosity of the molted grains several times 10⁷ P. That would naturally occur at around R = 0.1 au, and the temperature range is particularly interesting as it is close to the temperature required for sufficient thermal ionization to allow MRI activity (Gammie 1996). Thus, the region should occur near the inner edge of the dead zone, a location theorized to concentrate dust particles (Lyra et al. 2009).

We have then the potential for rapid planetesimal formation through the SI if that concentrated dust can be brought into regions of the correct temperature. That would occur either if those temperatures are the normal background temperature at the edge of the dead zone, or through secular evolution of the disk's radial profile. Thus, the picture of Boley et al. (2014), with terrestrial planetesimals forming in the hot inner disk, survives in a modified fashion. Significantly, it can potentially explain the recent observation of systems with tightly spaced inner planets (or STIPs) representing a significant fraction (more than 10%) of stellar systems.

However, that inner disk picture is a poor match for the formation of the Earth, and due to the large amount of migration that would be required, solar system-like planetary architectures in general. The Earth shows a smooth volatility dependent depletion pattern across a broad condensation temperature range of T ∼ 700–1400 K (Palme 2000; McDonough 2003), which is ill-fit by planetesimal formation in a narrow temperature range. Hubbard & Ebel (2014) put forth a model to explain the Earth's abundance pattern through time-dependent heating and cooling from an FU Orionis type event in the early solar system (Hartmann & Kenyon 1996).

FU Orionis type events are dramatic, long-lived accretion events, associated with luminosity increases of 4–6 mag and timescales of 50–100 yr (Hartmann & Kenyon 1996; Green et al. 2016). FU Orionis events (FUors) can raise the temperature at Earth's orbital position to about 1350 K (Hubbard & Ebel 2014), and drive the frost line out to about 40 au (Cieza et al. 2016). An FUor that strong would raise the temperature above 1000 K out to Mars's orbit, while a weaker, 4 magnitude FUor would still raise the temperature above 1000 K out to about Venus's orbit. A long-lived, 100 yr FUor outburst can satisfy the timescale constraint from Equation (32) out to a bit past Earth's orbit, while a 50-year FUor would still satisfy Equation (32) out to Venus's. While the frequency of FU Orionis events is not yet certain, statistics are consistent with a large fraction of protoplanetary disks undergoing one or several FU Orionis outbursts (Hartmann & Kenyon 1996).

We therefore propose that molten collisions during an FU Orionis event were the solution to the bouncing and fragmentation barriers for dry silicates at the locations of the terrestrial planets in our solar system. Thermal processing and molten collisions could explain both the existence of dry planetesimals at Earth's position, and the Earth's volatile abundance trend (Hubbard & Ebel 2014). The scenario suggests an in situ explanation for Mars's small size (Chambers 2014): the temperature and growth timescales required mean that the mechanism ceases to operate somewhere between Earth's and Mars's orbital positions (see Equation (32)), naturally reducing the number density of dry planetesimals beyond Earth's orbit. By cycling the temperature high enough to burn off organic surface layers, the model would also help explain why the Earth is carbon-poor even if carbon deposition is expected (Marty et al. 2013; Nuth et al. 2016). In our picture, molten collisions and evaporation/recondensation mechanisms (Hubbard 2017), rather than conventional dust coagulation, put the solar system's architecture in place very early in the Solar Nebula's existence. In addition to terrestrial planets, the duration and heating of an FUor leads to evaporation/recondensation cycles that might have driven icy planet formation in the outer disk (Hubbard 2017).

The briefer (duration of up to a few years) accretion events associated with Ex Lupi and Exors in general (Herbig 2007) are too brief for any given eruption to allow molten collisions to drive significant growth far beyond the inner edge of the dead zone, and are generally weaker than FUors. If Exors survive as episodic accreters long enough to accumulate decades of integrated outburst time, then they might drive intermittent molten dust growth to SI triggering sizes. Note, however, that dust drifts radially during quiescent periods as well as during active ones, so the radial drift condition (Equation (33)) becomes more stringent.