From Planetesimal to Planet in Turbulent Disks. II. Formation of Gas Giant Planets

Hiroshi Kobayashi; Hidekazu Tanaka

doi:10.3847/1538-4357/aacdf5

1. Introduction

Gas giant planets such as Jupiter, Saturn, and massive exoplanets are formed in protoplanetary disks containing solid and gas. In the core accretion scenario, once planetary embryos grow to ∼5 M_⊕, the planetary atmospheres that Mars-sized or larger embryos acquire become too massive to keep hydrostatic states, resulting in rapid gas accretion that forms gas giants (e.g., Mizuno 1980). Therefore, gas giant formation requires such massive solid cores to be formed within the disk lifetime (∼4 Myr).

In protoplanetary disks, submicron-sized dust grains accumulate to be planetesimals. Collisional growth of grains produces centimeter-sized or larger particles, which drift onto the host star on significantly short timescales due to moderate coupling with gas. If particles accumulate up to the Roche density, gravitational instability creates kilometer-sized or larger planetesimals (Goldreich & Ward 1973; Youdin & Goodman 2005; Michikoshi et al. 2012; Takeuchi & Ida 2012). On the other hand, collisional growth of dust grains naturally forms fluffy aggregates (Suyama et al. 2008; Suyama et al. 2012). The bulk density of dust aggregates becomes ∼10⁻⁴ g cm⁻³, and the collisional growth timescale of such bodies is much shorter than their drift timescale, so that kilometer-sized planetesimals are formed via direct collisional growth prior to radial drift (Okuzumi et al. 2012).

Planetesimals grow through mutual collisions. In a turbulent disk, the stirring by turbulence increases the random velocity of planetesimals v_r. For large v_r, the gravitational focusing of planetesimals is negligible and the orderly growth occurs until v_r ≳ 1.5v_esc (Kobayashi et al. 2016), where v_esc are the surface escape velocities of planetesimals. In the orderly growth, the mass-weighted average radius of planetesimals is comparable to the radius of the largest planetesimals. As planetesimals grow, v_esc increases to ∼v_r, and runaway growth occurs (Wetherill & Stewart 1989). The mass-weighted average radius of planetesimals at the onset of runaway growth, which is determined by turbulent strength (Kobayashi et al. 2016), almost remains the same after runaway growth.

Runaway growth forms a planetary embryo in each annulus of the disk. Embryos further grow mainly through collisions with surrounding planetesimals. As embryos become massive, their viscous stirring increases v_r, resulting in destructive collisions between planetesimals. Collisional fragments of planetesimals become further smaller via collisions between themselves. This collisional cascade decreases the sizes of bodies until radial drift. This process reduces the surrounding planetesimals, and embryo growth is then stalled (Kobayashi et al. 2010).

The collisional strength depends on the radius of planetesimals, r (e.g., Benz & Asphaug 1999). For r ≳ 1 km, the collisional outcome of a single collision is controlled by the self-gravity of colliders so that larger bodies are effectively stronger for collisions. Larger planetesimals, which have a longer collisional depletion timescales, contribute more to the growth of massive embryos. Embryo growth strongly depends on the mass-weighted average radius of planetesimals surrounding embryos (Kobayashi et al. 2010, 2011; Kobayashi & Dauphas 2013), which is mainly determined by the strength of turbulence (Kobayashi et al. 2016). Therefore, the growth and formation of solid cores of gas giants depends on the strength of disk turbulence.

In this paper, we investigate gas giant planet formation via core accretion in turbulent disks. In Section 2, we introduce the critical core mass from simple analysis. In Section 3, we explain the model of simulations. In Section 4, we perform simulations of collisional evolution of bodies from r = 1 km. We show the dependence of embryo growth on turbulent strength and collisional property and find the conditions for the formation of cores massive enough for the onset of rapid gas accretion to form gas giants. In Section 5, we discuss the growth timescale of embryos in turbulent disks, and the type I migration of growing embryos and the gas dispersal timescale of disks. In Section 6, we summarize our findings.

2. Critical Core Mass

The runaway growth of planetesimals produces a planetary embryo with mass M_E and radius R_E in each annulus of a protoplanetary disk. Once R_E is larger than its Bondi radius, R_B defined by ${R}_{{\rm{B}}}={{GM}}_{{\rm{E}}}/{c}_{{\rm{s}}}^{2}$ , the embryo has an atmosphere, where G is the gravitational constant and c_s is the sound velocity of gas. The density of a hydrostatic atmosphere, ρ_a, at the distance R from the center of the planetary embryo is given by (derivation in Appendix A; Mizuno 1980; Stevenson 1982; Inaba & Ikoma 2003)

$\begin{eqnarray}&&{\rho }_{{\rm{a}}}=\displaystyle \frac{\pi {\sigma }_{\mathrm{SB}}}{12\kappa {L}_{{\rm{e}}}{R}^{3}}{\left(\displaystyle \frac{{{GM}}_{{\rm{E}}}\mu {m}_{{\rm{H}}}}{k}\right)}^{4},\end{eqnarray} \tag{ 1 }$

where σ_SB is the Stefan–Boltzmann constant, κ is the opacity of the atmosphere, k is the Boltzmann constant, μ is the mean molecular weight, m_H is the mass of a hydrogen atom, and L_e is the planetary luminosity. The accretion of bodies onto planetary embryos mainly determines L_e so that

$\begin{eqnarray}&&{L}_{{\rm{e}}}=\displaystyle \frac{{{GM}}_{{\rm{E}}}}{{R}_{{\rm{E}}}}\displaystyle \frac{{{dM}}_{{\rm{E}}}}{{dt}}.\end{eqnarray} \tag{ 2 }$

The total atmospheric mass, M_A, is given by

$\begin{eqnarray}&&{M}_{{\rm{A}}}={\displaystyle \int }_{{R}_{{\rm{E}}}}^{{R}_{{\rm{B}}}}4\pi {R}^{2}{\rho }_{{\rm{a}}}{dR},\,\end{eqnarray} \tag{ 3 }$

$\begin{eqnarray}&&\,=\,\displaystyle \frac{{\pi }^{2}{\sigma }_{\mathrm{SB}}{G}^{3}{M}_{{\rm{E}}}^{3}{R}_{{\rm{E}}}}{3\kappa {\dot{M}}_{{\rm{E}}}}{\left(\displaystyle \frac{\mu {m}_{{\rm{H}}}}{k}\right)}^{4}\mathrm{ln}\left(\displaystyle \frac{{R}_{{\rm{B}}}}{{R}_{{\rm{E}}}}\right).\end{eqnarray} \tag{ 4 }$

Once M_A is comparable to M_E, the hydrostatic atmosphere is not maintained and then rapid gas accretion occurs to form gas giants. The embryo mass at the onset of rapid gas accretion is called the critical core mass, which is approximately estimated from ${M}_{{\rm{A}}}/{M}_{{\rm{E}}}\approx 1/3$ (Mizuno 1980; Stevenson 1982; Ikoma et al. 2000). Using Equation (4) with ${M}_{{\rm{A}}}={M}_{{\rm{E}}}/3$ , the critical core mass M_crit is given by

$\begin{eqnarray}&&{M}_{\mathrm{crit}}=5\,{M}_{\oplus }{\left(\displaystyle \frac{{M}_{{\rm{E}}}/{\dot{M}}_{{\rm{E}}}}{10\mathrm{Myr}}\right)}^{-3/4}{\left(\displaystyle \frac{\kappa }{0.01{\mathrm{cm}}^{2}{{\rm{g}}}^{-1}}\right)}^{3/4},\end{eqnarray} \tag{ 5 }$

where the dependence of $\mathrm{ln}({R}_{{\rm{B}}}/{R}_{{\rm{E}}})$ on M_E is ignored for this derivation.

Therefore, forming embryos with ${M}_{{\rm{E}}}\gtrsim 5\,{M}_{\oplus }$ are required for giant planet formation via core accretion. Note that the critical core mass may be smaller 5 M_⊕ because the depletion of planetesimals due to collisional fragmentation results in smaller ${\dot{M}}_{{\rm{E}}}$ (e.g., Kobayashi et al. 2011).

We perform simulations for the formation and growth of embryos below and find conditions for the formation of embryos with the critical core mass. Although we here simply estimate M_A using the radiative temperature gradient with constant κ, we calculate M_A from M_E and ${\dot{M}}_{{\rm{E}}}$ obtained in simulations, taking into account the convective temperature gradient as well as the radiative one with κ dependent on temperature (Inaba & Ikoma 2003).

3. Model for Planetary Accretion

The surface number density of planetesimals with masses m to m + dm orbiting the host star with mass M_* at the distance a, ${n}_{{\rm{s}}}(m,a){dm}$ , evolves through collisions and radial drift. The governing equation is given as (e.g., Kobayashi et al. 2016)

$\begin{eqnarray}\displaystyle \frac{\partial {{mn}}_{{\rm{s}}}(m,a)}{\partial t} & = & \displaystyle \frac{m}{2}{\displaystyle \int }_{0}^{\infty }{{dm}}_{1}{\displaystyle \int }_{0}^{\infty }{{dm}}_{2}{n}_{{\rm{s}}}({m}_{1},a){n}_{{\rm{s}}}({m}_{2},a)\\ & & \times \,K({m}_{1},{m}_{2})\delta (m-{m}_{1}-{m}_{2}+{m}_{{\rm{e}}})\\ & & -\,{{mn}}_{{\rm{s}}}(m){\displaystyle \int }_{0}^{\infty }{{dm}}_{2}{n}_{{\rm{s}}}({m}_{2},a)K(m,{m}_{2})\\ & & +\,\displaystyle \frac{\partial }{\partial m}{\displaystyle \int }_{0}^{\infty }{{dm}}_{1}{\displaystyle \int }_{0}^{{m}_{1}}{{dm}}_{2}{n}_{{\rm{s}}}({m}_{1},a){n}_{{\rm{s}}}({m}_{2},a)\\ & & \times \,K({m}_{1},{m}_{2}){\rm{\Psi }}(m,{m}_{1},{m}_{2})\\ & & -\,\displaystyle \frac{1}{a}\displaystyle \frac{\partial }{\partial a}[{{amn}}_{{\rm{s}}}(m,a){v}_{\mathrm{drift}}(m,a)],\end{eqnarray} \tag{ 6 }$

where the collisional kernel $K({m}_{1},{m}_{2})$ between bodies with masses m₁ and m₂ is given by

$\begin{eqnarray}&&K{({m}_{1},{m}_{2})=({h}_{{m}_{1},{m}_{2}}a)}^{2}\langle {{ \mathcal P }}_{\mathrm{col}}({m}_{1},{m}_{2})\rangle {\rm{\Omega }},\end{eqnarray} \tag{ 7 }$

with the reduced mutual Hill radius ${h}_{{m}_{1},{m}_{2}}=[({m}_{1}+{m}_{2})/3{M}_{* }{]}^{1/3}$ , the dimensionless mean collisional rate $\langle {{ \mathcal P }}_{\mathrm{col}}({m}_{1},{m}_{2})\rangle$ , and the Keplerian orbital frequency Ω, m_e, and ${\rm{\Psi }}(m,{m}_{1},{m}_{2})$ are the total and cumulative masses, respectively, of bodies produced by a single collision between bodies with masses m₁ and m₂, and v_drift is the drift velocity of a body due to gas drag. The collisional rate is determined according to the orbits of colliding bodies and the amount of their atmospheres: $\langle {{ \mathcal P }}_{\mathrm{col}}({m}_{1},{m}_{2})\rangle$ is given by a function of relative eccentricities and inclinations between bodies and the accretion rates of the bodies, which are summarized in Inaba et al. (2001). The enhancement of collisional cross section due to embryo's atmosphere (Inaba & Ikoma 2003) and due to gas drag for small bodies (Ormel & Klahr 2010) is taken into account. The atmospheric opacity is given by the sum of the opacities of gas and grains, given by (Inaba & Ikoma 2003)

$\begin{eqnarray}\kappa =\left\{\begin{array}{lll}0.01+4f\,{\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1} & \mathrm{for} & T\leqslant 170\,{\rm{K}},\\ 0.01+2f\,{\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1} & \mathrm{for} & 170\,{\rm{K}}\lt T\leqslant 1700\,{\rm{K}},\\ 0.01\,{\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1} & \mathrm{for} & T\gt 1700\,{\rm{K}},\end{array}\right.\end{eqnarray} \tag{ 8 }$

where f is the grain depletion factor compared to the interstellar one. The effective growth of dust aggregates depletes the opacity significantly (Okuzumi et al. 2012), so that we set $f={10}^{-4}$ .

The collisional kernel depends on orbital eccentricities e and inclinations i of colliders, which depend on mass m and distance a. The time variation of $e(m,a)$ and $i(m,a)$ via gravitational interactions between bodies depends on the mass distribution of bodies. Therefore, the evolution of $e(m,a)$ and $i(m,a)$ has to be treated simultaneously with the collisional evolution. The squares of dispersions for eccentricities and inclinations of bodies evolve according to

$\begin{eqnarray}&&\displaystyle \frac{{{de}}^{2}(m,a)}{{dt}}={\left.\displaystyle \frac{{{de}}^{2}}{{dt}}\right|}_{{\rm{g}}}+{\left.\displaystyle \frac{{{de}}^{2}}{{dt}}\right|}_{{\rm{d}}}+{\left.\displaystyle \frac{{{de}}^{2}}{{dt}}\right|}_{{\rm{c}}}+{\left.\displaystyle \frac{{{de}}^{2}}{{dt}}\right|}_{{\rm{t}}},\end{eqnarray} \tag{ 9 }$

$\begin{eqnarray}&&\displaystyle \frac{{{di}}^{2}(m,a)}{{dt}}={\left.\displaystyle \frac{{{di}}^{2}}{{dt}}\right|}_{{\rm{g}}}+{\left.\displaystyle \frac{{{di}}^{2}}{{dt}}\right|}_{{\rm{d}}}+{\left.\displaystyle \frac{{{di}}^{2}}{{dt}}\right|}_{{\rm{c}}}+{\left.\displaystyle \frac{{{di}}^{2}}{{dt}}\right|}_{{\rm{t}}},\end{eqnarray} \tag{ 10 }$

where the subscripts "g," "d," "c," and "t" indicate gravitational interaction between bodies (Ohtsuki et al. 2002), damping by gas drag (Adachi et al. 1976), collisional damping (Ohtsuki 1992), and turbulent stirring due to density fluctuation (Okuzumi & Ormel 2013) and aerodynamical friction (Völk et al. 1980), respectively. We use the damping rates by gas drag linear functions of e and i given by Inaba et al. (2001) because e and i ≪1, although the e and i damping rates by gas drag are greater for higher e and i (Kobayashi 2015). We use the damping rate and radial drift due to gas drag, taking into account three gas drag regimes (see details in Kobayashi et al. 2010).

For 10 m sized or larger bodies, the stirring by turbulence is mainly caused by turbulent density fluctuation, given by (Okuzumi & Ormel 2013)

$\begin{eqnarray}&&{\left.\displaystyle \frac{{{de}}^{2}}{{dt}}\right|}_{{\rm{t}}}\approx {f}_{{\rm{d}}}{\left(\displaystyle \frac{{{\rm{\Sigma }}}_{{\rm{g}}}{a}^{2}}{{M}_{* }}\right)}^{2}{\rm{\Omega }},\end{eqnarray} \tag{ 11 }$

$\begin{eqnarray}&&{\left.\displaystyle \frac{{{di}}^{2}}{{dt}}\right|}_{{\rm{t}}}\approx {\epsilon }^{2}{\left.\displaystyle \frac{{{de}}^{2}}{{dt}}\right|}_{{\rm{t}}},\end{eqnarray} \tag{ 12 }$

where Σ_g is the gas surface density, f_d is the dimensionless factor, and = 0.01 (see Kobayashi et al. 2016). According to magnetohydrodynamical simulations, f_d is given by (Okuzumi & Ormel 2013),

$\begin{eqnarray}&&{f}_{{\rm{d}}}=\displaystyle \frac{0.94\alpha }{{(1+4.5{H}_{\mathrm{res},0}/H)}^{2}},\end{eqnarray} \tag{ 13 }$

where H is the scale height of the disk, ${H}_{\mathrm{res},0}$ is the half vertical width of the dead zone, α is the dimensionless turbulent viscosity at the midplane. For simplicity, we set ${H}_{\mathrm{res},0}=H$ .

Dynamical stirring by large bodies increases collisional velocities of bodies. The collisional outcome is determined by m_e and ${\rm{\Psi }}(m,{m}_{1},{m}_{2})$ in Equation (6), which are mainly controlled by the fragmentation energy ${Q}_{{\rm{D}}}^{* }$ (see detailed setting in Appendix B). We use the following formula for ${Q}_{{\rm{D}}}^{* };$

$\begin{eqnarray}&&{Q}_{{\rm{D}}}^{* }={Q}_{0{\rm{s}}}{\left(\displaystyle \frac{r}{1\mathrm{cm}}\right)}^{{\beta }_{{\rm{s}}}}+{Q}_{0{\rm{g}}}{\rho }_{{\rm{s}}}{\left(\displaystyle \frac{r}{1\mathrm{cm}}\right)}^{{\beta }_{{\rm{g}}}}+{C}_{\mathrm{gg}}\displaystyle \frac{2{Gm}}{r},\end{eqnarray} \tag{ 14 }$

where ρ_s is the density of a body, ${Q}_{0{\rm{s}}},{Q}_{0{\rm{g}}},{\beta }_{{\rm{s}}}$ , β_g, and C_gg characterize the collisional strength. On the right-hand side of Equation (14), the first term is dominant for kilometer-sized or smaller bodies. The second term is important for larger bodies. The third term is determined by pure gravity, which is dominant for $r\gtrsim {10}^{7}\,\mathrm{cm}$ : we set C_gg = 10 (Stewart & Leinhardt 2009). The values of ${Q}_{0{\rm{s}}},{Q}_{0{\rm{g}}},{\beta }_{{\rm{s}}}$ , and β_g depend on the impact velocity (e.g., Benz & Asphaug 1999) and the structure of a body (e.g., Wada et al. 2013). For monolith bodies, these values are estimated from the interpolation of simulation data for the impact velocities v_imp of 3 and 5 km s⁻¹ obtained by Benz & Asphaug (1999), given by

$\begin{eqnarray}&&{Q}_{0{\rm{s}}}=1.6\times {10}^{7}{\left(\displaystyle \frac{{v}_{\mathrm{imp}}}{3\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{{\gamma }_{\mathrm{sv}}}\,{\mathrm{erg}{\rm{g}}}^{-1},\end{eqnarray} \tag{ 15 }$

$\begin{eqnarray}&&{\beta }_{{\rm{s}}}=-0.39{\left(\displaystyle \frac{{v}_{\mathrm{imp}}}{3\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{{\gamma }_{\mathrm{bsv}}},\end{eqnarray} \tag{ 16 }$

$\begin{eqnarray}&&{Q}_{0{\rm{g}}}=1.2{\left(\displaystyle \frac{{v}_{\mathrm{imp}}}{3\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{{\gamma }_{\mathrm{gv}}}\,\mathrm{erg}\,{\mathrm{cm}}^{3}{{\rm{g}}}^{-2},\end{eqnarray} \tag{ 17 }$

$\begin{eqnarray}&&{\beta }_{{\rm{g}}}=1.26{\left(\displaystyle \frac{{v}_{\mathrm{imp}}}{3\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{{\gamma }_{\mathrm{bgv}}},\end{eqnarray} \tag{ 18 }$

where γ_sv = −0.82, γ_gv = −0.31, γ_bsv = −0.080, and γ_bgv = 0.032. However, growing bodies via collisions are porous (e.g., Wada et al. 2013), and might not have monolith-like structures until melting. Therefore, bodies smaller than about 10 km have ${Q}_{{\rm{D}}}^{* }$ of porous bodies, which is given by (Wada et al. 2013),

$\begin{eqnarray}&&{Q}_{{\rm{s}}0}=4\times {10}^{7}\left(\displaystyle \frac{{v}_{\mathrm{imp}}}{100\,{\rm{m}}\,{{\rm{s}}}^{-1}}\right)\,{\mathrm{erg}{\rm{g}}}^{-1},\end{eqnarray} \tag{ 19 }$

$\begin{eqnarray}&&{\beta }_{{\rm{s}}}=0.\end{eqnarray} \tag{ 20 }$

We investigate the growth of melted bodies (monolith) and primordial bodies. For melted bodies, ${Q}_{{\rm{D}}}^{* }$ is given by Equations (15)–(18). For primordial bodies, ${Q}_{{\rm{D}}}^{* }$ is calculated from Equations (17)–(20). For melted and primordial bodies, ${Q}_{{\rm{D}}}^{* }$ are shown in Figure 1.

**Figure 1.** Fragmentation energy ${Q}_{{\rm{D}}}^{* }$ for primordial and monolith bodies with impact velocities of 50 and 100 m s⁻¹, given in Equations (14)–(19).
Download figure:
Standard image High-resolution image

4. Mass Evolution of Bodies

We perform simulations for planet formation in disks from 4.8 to 26 au via the time integration of Equations (6), (9), and (10) for M_* = M_☉, where M_☉ is the solar mass. The disk is divided into 10 annuli, and the mass distribution is described using mass bins with ratios between adjacent mass bins 1.05. We fix the bulk density of bodies at ${\rho }_{{\rm{s}}}=1\,{\rm{g}}\,{\mathrm{cm}}^{-3};$ the mass–radius relation is given by $m=4\pi {\rho }_{{\rm{s}}}{r}^{3}/3$ . The mass corresponding to the lowest mass bin is set to 4.2 × 10⁶ g (r = 1 m). We set the initial bodies at a radius of 1 km. The bodies initially have $e=i/\sqrt{\epsilon }=3{v}_{\mathrm{esc}}/a{\rm{\Omega }}\approx 9.9\times {10}^{-5}$ . The collisional growth of bodies is almost independent of the initial radius and orbits of bodies if the new-born planetesimals are smaller than those at the onset of runaway growth (Kobayashi et al. 2016).

The initial surface densities of solid and gas, ${{\rm{\Sigma }}}_{{\rm{g}},0}$ and ${{\rm{\Sigma }}}_{{\rm{s}},0}$ , respectively, are set to have power-law radial distributions,

$\begin{eqnarray}&&{{\rm{\Sigma }}}_{{\rm{g}},0}=1700{x}_{{\rm{g}}}{\left(\displaystyle \frac{a}{1\mathrm{au}}\right)}^{-1.5}\,{\rm{g}}\,{\mathrm{cm}}^{-2},\end{eqnarray} \tag{ 21 }$

$\begin{eqnarray}&&{{\rm{\Sigma }}}_{{\rm{s}},0}=30{x}_{{\rm{s}}}{\left(\displaystyle \frac{a}{1\mathrm{au}}\right)}^{-1.5}\,{\rm{g}}\,{\mathrm{cm}}^{-2},\end{eqnarray} \tag{ 22 }$

where x_g and x_s are the scaling factors, and the disk with x_g = x_s = 1 corresponds to the minimum-mass solar nebula (MMSN) model (Hayashi 1981). The solid surface density Σ_s evolves due to the radial drift of bodies, while the gas surface density artificially decreases from the formula ${{\rm{\Sigma }}}_{{\rm{g}}}={{\rm{\Sigma }}}_{{\rm{g}},0}\exp (-t/{\tau }_{\mathrm{gas}})$ with ${\tau }_{\mathrm{gas}}={10}^{7}\,\mathrm{years}$ . We discuss τ_gas and the time dependence of depletion in Section 5. In addition, we set the temperature at the disk midplane as

$\begin{eqnarray}&&T=280\left(\displaystyle \frac{a}{1\,\mathrm{au}}\right){\rm{K}},\end{eqnarray} \tag{ 23 }$

which affects R_B and the aerodynamical turbulent stirring. However, they are unimportant for the evolution of M_E and M_A.

We show the size and velocity evolution of bodies in the whole disk in Section 4.1 and focus on the size and velocity distributions around 5 au in Section 4.2. The evolution of planetary embryos dependent on α is shown in Section 4.3. We discuss the formation of cores with critical core masses from M_E and M_A obtained within the disk lifetimes by simulations in Section 4.4.

4.1. Size and Velocity Evolution of Bodies in the Whole Disk

We perform a simulation for planet formation with primordial strength in the disk with ${x}_{{\rm{s}}}={x}_{{\rm{g}}}=3$ (3MMSN) for α = 3 × 10⁻³, which results in a size and orbit evolution of bodies as shown in Figures 2–4. The surface density of solid bodies is given by ${{\rm{\Sigma }}}_{{\rm{s}}}(a)=\int {{mn}}_{{\rm{s}}}(m,a){dm}=\int {m}^{2}{n}_{{\rm{s}}}(m,a)d\mathrm{ln}m$ , so that we define the surface density of bodies in the logarithmic mass interval as

$\begin{eqnarray}&&{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r,a)\equiv {m}^{2}{n}_{{\rm{s}}}(m,a),\end{eqnarray} \tag{ 24 }$

with r corresponding to m. Figure 2 shows ${\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r,a)$ . The collisional growth of bodies occurs inside-out. Orderly growth initially occurs so that ΔΣ_s have a maximum around radii of the largest bodies in each annulus for t ≲ 1 Myr (Figures 2(a), (b)). Larger bodies tend to have greater e and i for r ≳ 10 m (Figures 3(a), (b) and 4(a), (b)). Density-fluctuation turbulent stirring and gas damping mainly control e and i for r ≳ 10 m, while the fluid-dynamical turbulent stirring is dominant instead of that by density fluctuation for r ≲ 10 m. The runaway growth of bodies produces planetary embryos inside 8 au by 4 Myr (Figures 2(c), (d)). The stirring by embryos increases e and i (Figures 3(c), (d) and 4(c), (d)), which induces collisional fragmentation due to high-speed collisions. Collisional fragments become smaller due to further collisions between themselves until radial drift, which reduces the surface density of bodies with r ≲ 10–100 km (compare Figures 2(c) and (d)). Planetary-embryo formation in the outer disk induces collisional fragmentation and radial drift, which supplies small bodies in the inner disk. However, the net flux of small bodies reduces the surface density of bodies around planetary embryos in the inner disk. This stalls the growth of embryos (e.g., Kobayashi et al. 2010, 2011). In order to see detailed mass and orbit evolution, we focus on the evolution at 5.2 au in the following paragraphs.

**Figure 3.** Orbital eccentricity distribution obtained from the same simulation as Figure 2.
Download figure:
Standard image High-resolution image

**Figure 4.** Orbital inclination distribution obtained from the same simulation as Figure 2.
Download figure:
Standard image High-resolution image

4.2. Size and Velocity Distributions around 5 au

As shown in Figures 2–4, collisional fragmentation naturally occurs during planet formation and plays an important role for planet formation (e.g., Wetherill & Stewart 1993; Inaba et al. 2003; Kobayashi et al. 2010, 2011, 2012; Kobayashi & Dauphas 2013). Therefore collisional fragmentation needs to be taken into account for planet formation. For comparison, however, we first show the size and orbit evolution of bodies in the case without collisional fragmentation (m_e = Ψ = 0 in Equation (6)). Figure 5 shows the evolution of ΔΣ_s, e, and i of bodies at 5.2 au in the 3MMSN disk with α = 3 × 10⁻³ as a function of r corresponding to m. For t ≲ 1 Myr, the surface densities ΔΣ_s(r) have single peaks, which move to large r. This is caused by orderly growth of bodies. The peak radius r_pk is approximated to be the mass-weighted average radius of bodies. For t ≳ 1 Myr, the size distribution becomes wider, and then collisional evolution produces another peak at the high-mass end, which indicates planetary embryos. This is caused by the runaway growth of bodies. After the onset of runaway growth, r_pk does not change significantly so that the peak radius at the onset of runaway growth, r_rg, is approximated to be r_pk, which is estimated to ∼100 km from the size distribution of bodies at ∼1 Myr. The onset of runaway growth occurs if v_r ≲ 1.5v_esc for bodies of r ∼ r_pk (Kobayashi et al. 2016). Since v_r/a Ω ≈ e and ${v}_{\mathrm{esc}}/a{\rm{\Omega }}\approx 6\times {10}^{-3}(r/100\,\mathrm{km})$ , we estimate r_rg ≈ 100 km from the data of e at 1.9 Myr in Figure 5, which is consistent with the mass distribution in Figure 5.

**Figure 5.** Size distribution of ΔΣ_s of bodies (top) and their eccentricities e (bottom; solid curves) and inclinations i (bottom; dotted curves) at 5.2 au in a disk with x_s = x_g = 3 for α = 3 × 10⁻³ without collisional fragmentation. The curves indicate those at 1.7 × 10⁵ years (blue), 8.7 × 10⁵ years (purple), 1.9 × 10⁶ years (magenta), and 3.2 × 10⁶ years (red). The time evolution results in the existence of large bodies, so that the curves move from left to right.
Download figure:
Standard image High-resolution image

The size distributions of bodies have single peaks at r = r_pk prior to runaway growth, t ≲ 1 Myr. Bodies with r ≲ r_pk have similar slopes of the size distributions of bodies. Even after runaway growth, bodies with r ≲ r_rg have similar slopes. The slope of the mass distribution of bodies is estimated by eye to be $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 1.8$ . For a collisional cascade with ${Q}_{{\rm{D}}}^{* }$ and ${v}_{{\rm{r}}}$ independent of m, $d\mathrm{ln}{n}_{{\rm{s}}}(m)/d\mathrm{ln}m=1/2$ (Dohnanyi 1969; Tanaka et al. 1996). The slopes obtained from simulations without collisional fragmentation are steeper than the slope of the simple collisional cascade. On the other hand, the slopes of bodies with r ranging from ≈r_rg to 10⁵ km are $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx -2.0$ in t ≳ 1 Myr. The collisional cross section is proportional to ${m}^{5/3}/{v}_{{\rm{r}}}^{2}$ due to gravitational focusing and ${v}_{{\rm{r}}}\propto {m}^{-1/2}$ due to dynamical friction, which analytically gives $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r=-2$ (Makino et al. 1998). The slopes formed via runaway growth are roughly explained by gravitational focusing and dynamical friction. In a more detailed analysis, the runaway-growth slopes bend slightly (see Morishima 2017).

Figure 6 shows the evolution of bodies at 5.2 au in the same condition as Figures 2–4 (i.e., the case with collisional fragmentation for primordial strength). For t ≲ 1 Myr, ΔΣ_s (r) has a single peak, although ΔΣ_s (r) has a long tail at the low-mass side, which is produced by collisional fragmentation. These small bodies make collisional damping effective, which reduces v_r. Runaway growth occurs slightly earlier than the case without collisional fragmentation, resulting in r_rg ≈ 60 km smaller than r_rg without fragmentation. After the onset of runaway growth (t ≳ 1 Myr), r_pk moves insignificantly. The stirring by planetary embryos created by runway growth increases e and i of the surrounding planetesimals, which induces collisional fragmentation of planetesimals of r ∼ r_rg. The solid surface density is decreased via collisional fragmentation and radial drift of the resulting fragments (Kobayashi et al. 2010, 2011). Planetary embryos grow through the accretion of planetesimals and fragments until their depletion.

The size distributions of bodies have single peaks at r = r_pk prior to runaway growth, which is similar to the case without collisional fragmentation. However, bodies smaller than the peaks have multiple slopes: the slopes are estimated by eye to be $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 2.5$ for r ≪ 30 m, $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.6$ from r ≈ 30 m to ∼0.2r_pk, and $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 2.0$ for r ≈ 0.2r_pk–r_pk (see Figure 6). The slope for large bodies of r = 0.2r_pk–r_pk is similar to that for no fragmentation, which is simply determined by collisional growth of bodies in orderly growth. The intermediate-sized bodies for r = 30 m–0.2r_pk, which are mainly produced by erosive collisions of bodies with r ∼ r_pk, have tiny e and i so that collisional fragmentation is negligible. Their collisional growth results in the slope that is determined by the collisional cascade to the positive mass direction, $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\,=1/2$ (Tanaka et al. 1996). For bodies with r ≪ 30 m, radial drift as well as the collisional cascade affects the mass distribution so that the slope is given by $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r=5/2$ .³ On the other hand, the runaway growth produces a different power-law size distribution of bodies with r ≳ r_rg; $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\,=-2.0$ , which is caused by runaway growth similar to the case without fragmentation. For r ≲ r_rg, the slope is controlled by collisional fragmentation due to large e and i. A collisional cascade therefore occurs, resulting in $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\,=(1+3p)/(2+p)$ with $3p=d\mathrm{ln}{Q}_{{\rm{D}}}^{* }/d\mathrm{ln}r-2d\mathrm{ln}{v}_{{\rm{r}}}/d\mathrm{ln}r$ (Kobayashi et al. 2010). Although $d\mathrm{ln}{Q}_{{\rm{D}}}^{* }/d\mathrm{ln}r$ depends on collisional velocities, p ≈ 0 is roughly estimated and then $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.5$ , which is similar to the slope for r ∼ 30 m–50 km. For r ≲ 30 m, e and i are so small for the bodies that collisional fragmentation no longer occurs. The slope is determined by collisional growth and radial drift, which is almost the same as the slope of small bodies prior to runaway growth.

Figure 7 shows the result of another simulation with monolith strength for the same disk condition. The evolution of the size distribution is similar to the case for primordial strength, but the onset of runaway growth occurs slightly earlier, resulting in r_rg ≈ 30 km, because the small bodies produced by collisional fragmentation prior to runaway growth result in slightly effective collisional damping. After the onset of runaway growth (t ≳ 1 Myr), forming planetary embryos induce collisional fragmentation of planetesimals, which mainly produces bodies with r ∼ 10 m–10 km through the collisional cascade. The surface density ΔΣ_s (r) has a peak at a radius of 10–100 m because the collisional cascade starting from collisional fragmentation of bodies with r ∼ r_rg is stalled by e and i damping by gas drag in the Stokes regime for bodies with r ≲ 10–100 m (Kobayashi et al. 2010). The magnitudes of such peaks depend on r_rg, ${Q}_{{\rm{D}}}^{* }$ , Σ_s, and so on. Planetary embryos mainly grow through collisions with bodies of ∼r_rg or 10–100 m (Kobayashi et al. 2010, 2011). On the other hand, ∼300 m sized bodies are quickly destroyed via collisions with small bodies because of low ${Q}_{{\rm{D}}}^{* }$ .

**Figure 7.** Same as Figure 5, but with collisional fragmentation for monolith strength. The curves indicate those at $2.1\times {10}^{5}\,\mathrm{years}$ (blue), 4.9 × 10⁵ years (purple), 1.1 × 10⁶ years (magenta), and 3.1 × 10⁶ years (red).
Download figure:
Standard image High-resolution image

The mass distribution prior to runaway growth is similar to the case with primordial strength. The slopes of small bodies are almost the same, while the population of small bodies is larger because of a high production caused by collisional fragmentation due to low ${Q}_{{\rm{D}}}^{* }$ . On the other hand, although runaway growth results in a slope for r ≳ 30 km ≈ r_rg similar to the case with primordial strength, bodies with r ≲ r_rg have a "wavy" structure, where the slopes $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r$ vary in small size ranges; $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.1$ for r ≈ 3–60 km, $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.2$ for r ≈ 0.3–3 km, $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx -0.2$ for r ≈ 30–300 m, and $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.2$ for r ≲ 30 m (t = 1.9 Myr in Figure 7). The slope controlled by the collisional cascade is analytically estimated to be $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.1$ for r ≳ 500 m and $d\mathrm{ln}{\rm{\Delta }}{{\rm{\Sigma }}}_{{\rm{s}}}(r)/d\mathrm{ln}r\approx 0.0$ for r ≈ 30–500 m. The collisional cascade may explain the slope only around r ≈ 3–60 km. The "wavy" pattern for r ≲ 3 km is formed due to large ${v}_{{\rm{r}}}^{2}/{Q}_{{\rm{D}}}^{* }$ (e.g., Campo Bagatin et al. 1994; Durda & Dermott 1997; Thébault et al. 2003; Krivov 2007; Löhne et al. 2008). The value of ${v}_{{\rm{r}}}^{2}/{Q}_{{\rm{D}}}^{* }$ becomes significant around r ≈ 500 m because of minimum ${Q}_{{\rm{D}}}^{* }$ (see Figure 1), which produces a bump at r ≈ 500 m in the size distribution.

Figure 8 shows the growth of planetary embryos obtained from the simulations shown in Figures 5–7. For t ≲ 0.4–1 Myr, embryo masses M_E grow as M_E ∝ t³ because of orderly growth (Kobayashi et al. 2016). The onset of runaway growth occurs around 0.4–1 Myr, resulting in a strong time dependence of M_E. The rapid growth occurs until M_E ∼ 10⁻²–10⁻¹ M_⊕. The growth timescale in runaway growth depends on the collisional model, which is caused by different r_rg. After the rapid growth, slow growth occurs again, which is called oligarchic growth. Planetary embryos grow through surrounding planetesimals whose typical radii are r_rg (see Figures 5–7). The viscous stirring by embryos increases e and i of planetesimals, which induces the collisional cascade. The resulting bodies of ∼10–100 m rapidly drift inward. Therefore, the oligarchic growth is stalled by the depletion of planetesimals due to collisional fragmentation and the radial drift of fragments. The growth of embryos at the oligarchic stage depends on ${Q}_{{\rm{D}}}^{* }$ for r ≲ r_rg, which controls the depletion of the solid surface density (see Figures 6 and 7). The oligarchic growth for primordial strength stalls later than that for melted materials. For primordial strength, the characteristic radii of planetesimals, r_rg, are large due to the late onset of runaway growth, which has great ${Q}_{{\rm{D}}}^{* }$ . The insignificant collisional depletion of planetesimals results in the long accretion of planetesimals onto embryos so that embryos grow to be massive. In addition, the collisional cascade grinds planetesimals down to 10–100 m. Strong gas drag damps v_r for r ≲ 100 m so that the collisional cascade stalls without destructive collisions. The radius of bodies at the low-mass end of the collisional cascade, r_cc, depends on ${Q}_{{\rm{D}}}^{* };$ the primordial strength has large r_cc. Bodies with r ∼ r_cc have so small v_r that they effective accrete onto embryos, while their drift timescales are short because of strong gas drag. The accretion efficiency of bodies with r ∼ r_cc depends on r_cc (Kobayashi et al. 2010). For primordial material, embryos effectively grow via accretion of bodies with r ∼ r_cc because of slow radial drift for large r_cc due to great ${Q}_{{\rm{D}}}^{* }$ .

4.3. Turbulent Strength Dependence of Embryo Growth

We carry out collisional-evolution simulations in 3MMSN disks with α = 3 × 10⁻⁵–3 × 10⁻³ for primordial ${Q}_{{\rm{D}}}^{* }$ , and compare this with the results without collisional fragmentation (see Figure 9). For weak turbulence (α ≲ 10⁻⁴), the growth of embryos is similar to the case without fragmentation until embryos exceed the mass of the Moon (∼10⁻² M_⊕). For α = 3 × 10⁻³, the early collision fragmentation due to strong turbulent stirring leads to effective collisional damping, inducing the onset of runaway growth earlier than in the case without collisional fragmentation. After runaway growth, large embryos control v_r instead of turbulence, which induces significant collisional fragmentation of planetesimals. The collisional cascade results in bodies of r ∼ r_cc ∼ 10–100 m. The effective accretion of such small bodies leads to the rapid growth of embryos (compare the results with or without collisional fragmentation). However, the growth is stalled by the depletion of small bodies due to collisional fragmentation and radial drift.

**Figure 9.** Time evolution of planetary-embryo masses in a 3MMSN disk with α = 3 × 10⁻⁵ (blue), 3 × 10⁻⁴ (magenta), and 3 × 10⁻³ (red) for a collisional model with primordial ${Q}_{{\rm{D}}}^{* }$ (solid) and without fragmentation (dotted).
Download figure:
Standard image High-resolution image

For small α, the runaway growth occurs early, resulting in small ${r}_{\mathrm{rg}}$ . If we ignore collisional fragmentation, the early formation and growth of embryos due to small α results in large embryos (see Figure 9). However, collisional fragmentation affects the growth of embryos. Small planetesimals tend to be destroyed via collision due to small ${Q}_{{\rm{D}}}^{* }$ for weak self-gravity (see Figure 1). The stirring by small embryos increases v_r of planetesimals moderately, which can induce collisional fragmentation of planetesimals for small r_rg. The resulting fragments initially accelerate the growth of embryos, while the depletion of the surrounding planetesimals stalls embryo growth. Weak turbulence (small α) results in low masses of embryos in the late stage (t ≳ 4 Myr in Figure 9), while large α enhances r_rg and then M_E at the late stage.

Figure 10 shows embryo growth for monolith ${Q}_{{\rm{D}}}^{* }$ . The α dependence is similar to the result for primordial ${Q}_{{\rm{D}}}^{* };$ higher α produces more massive embryos. However, the final embryos are smaller. The accretion of small bodies with r ∼ r_cc ≈ 10–100 m contributes to embryo growth. However, the depletion of bodies with r ∼ r_cc due to radial drift stalls embryo growth. Low ${Q}_{{\rm{D}}}^{* }$ for r ∼ 10–100 m makes r_cc small. Therefore, the growth for monolith strength stalls earlier than that for the primordial strength (see Figures 9 and 10).

4.4. Forming Cores with Critical Core Masses

As discussed in Section 2, gas giant formation via gas accretion requires embryos to be larger than ∼5 M_⊕ prior to significant gas depletion. Figure 11 shows embryo masses M_E at 5.2 au at 4 Myr with primordial strength. Massive embryos tend to be formed in the disks with strong turbulence (large α), while too strong turbulence cannot produce massive embryos because the onset of runaway growth is too late. In 3MMSN disks, embryos grow up to 1–2 M_⊕, which is too small to start gas accretion. However, the critical core mass depends on ${\dot{M}}_{{\rm{E}}}$ . To confirm the possibility of gas giant formation, we calculate the masses of static atmospheres around embryos, M_A, based on the analytical model for atmospheric radial density profile ignoring the gravity of atmospheres (Inaba & Ikoma 2003). The mass ratio M_A/M_E is lower than 0.2, so that runaway gas accretion does not occur within 4 Myr. In the later stage, embryos grow more massive (see Figure 9); ${M}_{{\rm{A}}}/{M}_{{\rm{E}}}\gtrsim 1/3$ at t ≈ 5 Myr and 6–7 Myr for α = 3 × 10⁻⁴ and 3 × 10⁻³, respectively. In 2MMSN, the α dependence of M_E and M_A is similar to that in 3MMSN, while M_E and M_A are smaller. Therefore, embryos cannot reach the critical core mass within the disk lifetime ≲4 Myr in disks that are less massive than 3MMSN.

**Figure 11.** Masses of largest bodies, embryo mass M_E (left panel) and the mass ratio of atmosphere to core M_A/M_E (right panel) for bodies with primordial strength at 4 Myr at 5.2 au in disks with x_g = x_s = 2 (2MMSN; green), x_g = x_s = 3 (3MMSN; blue), x_g = 2 and x_s = 8 (red), and x_g = 3 and x_s = 12 (magenta). Once M_A/M_E exceeds 1/3, runaway gas accretion occurs to form gas giants.
Download figure:
Standard image High-resolution image

The collisional evolution of fluffy dust aggregates overcomes the radial drift barrier if the dust aggregates that drift most effectively are controlled in the Stokes gas drag regime, resulting in planetesimals within ∼10 au (Okuzumi et al. 2012). During collisional evolution, the radial drift of dust aggregates induces a pile-up in the planetesimal-forming region so that the solid surface density increases by a factor 3–4 (Okuzumi et al. 2012). According to the result, we set ${x}_{{\rm{s}}}/{x}_{{\rm{g}}}=4$ .⁴ In the solid-enhanced disks with x_g = 3 and x_s = 12, embryo masses at 4 Myr exceed $5\,{M}_{\oplus }$ for α ≳ 3 × 10⁻⁴ so that their atmospheric masses M_A are much higher than M_E/3 (Figure 11). For x_g = 2 and x_s = 8, M_E for α ≲ 10⁻³ is similar to that in the case of 3MMSN, while large M_E for α ≳ 10⁻³ results in M_A/M_E ≳ 1/3. Figure 12 shows the mass evolution of embryos in the disk with x_g = 2 and x_s = 8. The growth is stalled at t ≳ 0.1 Myr for α = 3 × 10⁻⁴. Embryo growth occurs within ∼1 Myr even for α ≳ 10⁻³, which forms massive cores. Therefore, in the solid-enhanced disks, runaway gas accretion for gas giant formation occurs within the disk lifetime of ∼4 Myr.

**Figure 12.** Same as Figure 9, but in the disk with x_g = 2 and x_s for α = 3 × 10⁻⁴ (magenta), $3\times {10}^{-3}$ (red), and 3 × 10⁻² (blown).
Download figure:
Standard image High-resolution image

**Figure 12.** Same as Figure 9, but in the disk with x_g = 2 and x_s for α = 3 × 10⁻⁴ (magenta), $3\times {10}^{-3}$ (red), and 3 × 10⁻² (blown).
Download figure:
Standard image High-resolution image

For monolith strength, embryo masses are lower than the mass for primordial strength (Figure 13). In the 3MMSN disks, embryos are too small to start gas accretion. Even in the 5MMSN disk, M_A/M_E is smaller than 0.2. Therefore, gas giant formation via core accretion is difficult in a disk that is less massive than 5MMSN. In solid-enhanced disks with x_g = 3 and x_s = 12, embryo masses reach ∼1–2 M_⊕ for α ≳ 3 × 10⁻³. For α = 3 × 10⁻², M_A/M_E > 1/3. In the disk with x_g = 5 and x_s = 20, M_A/M_E > 1/3 for α ∼ 3 × 10⁻³. Therefore, melted bodies may form gas giant planets via core accretion in more massive disks compared to primordial material.

5. Discussion

5.1. Growth Timescale

Prior to the runaway growth, the mass distribution of bodies is approximated to be a single-mass population (Kobayashi et al. 2016), so that the collisional timescale at the onset of runaway growth is estimated to be

$\begin{eqnarray}&&{\tau }_{\mathrm{col}}\approx \displaystyle \frac{4i{\rho }_{{\rm{s}}}{r}_{\mathrm{rg}}}{3e{{\rm{\Sigma }}}_{{\rm{s}}}{\rm{\Omega }}}.\end{eqnarray} \tag{ 25 }$

The balance between turbulent stirring and gas drag gives i/e as prior to the onset of runaway growth (Kobayashi et al. 2016), while the gravitational interaction during runaway growth results in energy equipartition, i/e = 0.5. The embryo formation timescale via runaway growth is proportional to τ_col with i/e = 0.5 (Kobayashi et al. 2010; Ormel et al. 2010). The subsequent embryo-growth timescale τ_g, which depends on the accretion of planetesimals or small bodies (Kobayashi et al. 2010, 2011), is simply approximated to be τ_g ≈ 3τ_col,

$\begin{eqnarray}&&{\tau }_{{\rm{g}}}=4.5{\left(\displaystyle \frac{{x}_{{\rm{s}}}}{3}\right)}^{-1}\left(\displaystyle \frac{{r}_{\mathrm{rg}}}{30\,\mathrm{km}}\right){\left(\displaystyle \frac{a}{5.2\mathrm{au}}\right)}^{3}\,\mathrm{Myr}.\end{eqnarray} \tag{ 26 }$

The random velocity is determined by the balance between turbulent stirring and collisional damping prior to runaway growth. Once v_r ≈ 1.5v_esc, runaway growth occurs. The radius of bodies at the onset of runaway growth is given by (Kobayashi et al. 2016)

$\begin{eqnarray}&&{r}_{\mathrm{rg}}=32\left(\displaystyle \frac{\alpha }{3\times {10}^{-3}}\right){\left(\displaystyle \frac{{x}_{{\rm{g}}}}{3}\right)}^{2}{\left(\displaystyle \frac{{x}_{{\rm{s}}}}{3}\right)}^{-1}{\left(\displaystyle \frac{a}{5.2\mathrm{au}}\right)}^{1.5}\,\mathrm{km}.\end{eqnarray} \tag{ 27 }$

From Equations (26) and (27), the timescale of embryo growth depending on turbulence strength is given by

$\begin{eqnarray}&&{\tau }_{{\rm{g}}}=4.8{\left(\displaystyle \frac{{x}_{{\rm{g}}}}{{x}_{{\rm{s}}}}\right)}^{2}\left(\displaystyle \frac{\alpha }{3\times {10}^{-3}}\right){\left(\displaystyle \frac{a}{5.2\mathrm{au}}\right)}^{4.5}\,\mathrm{Myr}.\end{eqnarray} \tag{ 28 }$

As seen in Figure 11, growing embryos at t = 4 Myr for α = 3 × 10⁻³ are less massive than those for α = 3 × 10⁻⁴ for x_g = x_s, while embryos at t = 4 Myr increase with α for ${x}_{{\rm{s}}}/{x}_{{\rm{g}}}=4$ . That is explained by the dependence of τ_g on x_s/x_g; For x_s/x_g = 4, τ_g < 4 Myr even for α ≲ 3 × 10⁻², resulting in embryo formation within the disk lifetime.

Collisional fragmentation of planetesimals stalls embryo growth so that large ${r}_{\mathrm{rg}}$ tends to form embryos that are massive enough. The formation of cores with the critical core masses requires r_rg ≳ 10 km around 5 au for primordial strength in the disks with x_g ≈ 3. On the other hand, large r_rg results in a long embryo-growth timescale (see Equation (26)); τ_g ≪4 Myr is needed. These conditions are not satisfied for 3MMSN, while the solid-enhanced disks create massive cores for α ≳ 10⁻³ because of these conditions (see Figure 12).

5.2. Planetary Migration and Gas Dispersal

Planetary embryos migrate due to density waves caused by interaction with disks; the migration timescale is estimated to be (e.g., Tanaka et al. 2002)

$\begin{eqnarray}{\tau }_{\mathrm{mig}} & = & {\gamma }^{-1}{\left(\displaystyle \frac{{M}_{{\rm{E}}}}{{M}_{* }}\right)}^{-1}{\left(\displaystyle \frac{{{\rm{\Sigma }}}_{{\rm{g}}}{a}^{2}}{{M}_{* }}\right)}^{-1}{\left(\displaystyle \frac{{c}_{{\rm{s}}}}{{v}_{{\rm{K}}}}\right)}^{2}{{\rm{\Omega }}}^{-1},\\ & \approx & 1.2{\gamma }^{-1}{\left(\displaystyle \frac{{M}_{{\rm{E}}}}{{M}_{\oplus }}\right)}^{-1}{\left(\displaystyle \frac{{{\rm{\Sigma }}}_{{\rm{g}}}}{450{\rm{g}}{\mathrm{cm}}^{-2}}\right)}^{-1}{\left(\displaystyle \frac{{c}_{{\rm{s}}}/{v}_{{\rm{K}}}}{0.05}\right)}^{2}\\ & & \times \,{\left(\displaystyle \frac{a}{5.2\mathrm{au}}\right)}^{-1/2}{\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{3/2}\mathrm{Myr},\end{eqnarray} \tag{ 29 }$

where γ is the migration coefficient and the value of Σ_g is given from that at 5.2 au in the 3 MMSN disk. In the isothermal disk, γ ≈ 4 (Tanaka et al. 2002) and the corotation torque may reduce γ ∼ 1 (Paardekooper et al. 2011). If the migration timescale is shorter than the time required to reach the critical core mass, embryos are lost due to migration prior to gas giant formation. We here estimate the migration timescale for embryos obtained in our simulations using γ = 1.

**Figure 13.** Same as Figure 11, but for monolith strength in the disks with x_g = x_s = 3 (3MMSN; blue), x_g = x_s = 5 (5MMSN; gray), x_g = 3 and x_s = 12 (magenta), and x_g = 5 and x_s = 20 (purple).
Download figure:
Standard image High-resolution image

For primordial strength, embryos as large as the critical core mass are formed within several Myr in the disk with x_s/x_g ≈ 4 for α ≈ 10⁻³–10⁻¹ (Figure 11). Figure 14 shows 3M_A/M_E at 5.2 au for α = 3 × 10⁻⁴ and 3 × 10⁻³; M_A becomes M_E/3 at t ≈ 2 Myr and 3–4 Myr for α = 3 × 10⁻⁴ and 3 × 10⁻³, respectively, and embryos then reach the critical core masses. Using M_E obtained from the simulations and Equation (29), we calculate τ_mig/t. Embryos may grow rather than migration unless τ_mig/t ≲ 1. Therefore, embryos migrate inward prior to the formation of embryos with the critical core masses; gas giant formation is inhibited by migration.

**Figure 14.** Dimensionless migration timescales τ_mig/t (solid curves) calculated from the results of simulations with primordial strength in disks with x_g = 3 and x_s = 12 for α = 3 × 10⁻⁴ (blue) and $3\times {10}^{-3}$ (red). Dotted curves indicate 3M_A/M_E. The onset of runaway gas accretion is estimated from 3M_A/M_E = 1.
Download figure:
Standard image High-resolution image

The disk dispersion timescale is required to be comparable to or shorter than τ_mig. For τ_gas = 1 Myr, τ_mig/t > 1 is satisfied in the simulations. However, the gas surface density significantly decreases prior to the onset of gas accretion of the core (M_A = M_E/3). Therefore, a more realistic gas dispersal should be taken into consideration.

For accretion disks with constant α, for simplicity, the gas surface density is proportional to ${{\rm{\Sigma }}}_{{\rm{g}}}\propto {(1+t/{\tau }_{\mathrm{dep}})}^{-1.5}$ (Lynden-Bell & Pringle 1974), where ${\tau }_{\mathrm{dep}}={\alpha }^{-1}{({v}_{{\rm{K}}}/{c}_{{\rm{s}}})}^{2}{{\rm{\Omega }}}^{-1}/3$ at the disk radius r_cut of the exponential cutoff for the surface density. We estimate

$\begin{eqnarray}{\tau }_{\mathrm{dep}} & = & 0.63{\left(\displaystyle \frac{\alpha }{3\times {10}^{-3}}\right)}^{-1}{\left(\displaystyle \frac{{c}_{{\rm{s}}}}{0.1{v}_{{\rm{K}}}}\right)}^{-2}{\left(\displaystyle \frac{{r}_{\mathrm{cut}}}{50\mathrm{au}}\right)}^{3/2}\\ & & \times \,{\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{-1/2}\,\mathrm{Myr}.\end{eqnarray} \tag{ 30 }$

Therefore, Σ_g may decrease on a timescale of 1 Myr.

We perform the simulations assuming ${{\rm{\Sigma }}}_{{\rm{g}}}\propto {{\rm{\Sigma }}}_{{\rm{g}},0}/(1+t/{{\tau }_{\mathrm{dep}})}^{3/2}$ with τ_dep = 0.5 Myr (Figure 15), which mimics the accretion disk. Embryos reach the critical core mass at 3–4 Myr. The early gas depletion weakens the damping of eccentricity and inclination due to gas drag, and then the accretion of planetesimals on embryos is suppressed due to large e and i. However, small bodies produced by collisional fragmentation feel the Stokes gas drag, which is independent of gas density. The early gas depletion does not affect the accretion of such small bodies significantly. Therefore, embryos exceed the critical core mass within 4 Myr. On the other hand, the early gas depletion significantly affects migration. The gas depletion with ${\tau }_{\mathrm{dep}}\lesssim 1$ Myr prolongs τ_mig, resulting in τ_mig/t > 1. Therefore, embryos may start rapid gas accretion prior to migration.

**Figure 15.** Same as Figure 14, but for ${{\rm{\Sigma }}}_{{\rm{g}}}={{\rm{\Sigma }}}_{{\rm{g}},0}/{(1+t/{\tau }_{\mathrm{dep}})}^{3/2}$ with ${\tau }_{\mathrm{dep}}=5\times {10}^{5}\,\mathrm{years}$ .
Download figure:
Standard image High-resolution image

The disk with x_g = 3 initially has ∼17M_J, where we set the disk edge at 50 au and M_J is the mass of Jupiter. The disk mass decreases with $\propto {(1+t/{\tau }_{\mathrm{dep}})}^{-1/2}$ (Lynden-Bell & Pringle 1974) so that the disk mass becomes 6M_J at t = 4 Myr for ${\tau }_{\mathrm{dep}}=0.5\,\mathrm{Myr}$ . The gas accretion of embryos with critical core masses occurs not only from around embryos with the critical core masses, but also from the whole disk during disk evolution (Tanigawa & Ikoma 2007), so that embryos may acquire an atmosphere comparable to that of Jupiter. Therefore, gas accretion occurs in such a depleting disk, which saves forming gas giants from type II migration (e.g., Ida & Lin 2008). It should be noted that ${\tau }_{\mathrm{dep}}$ does not correspond to the disk lifetime inferred from infrared observations. In the self-similar solution for accretion disks (Lynden-Bell & Pringle 1974), the disk evolution timescale is prolonged for $t\gg {\tau }_{\mathrm{dep}}$ . Even for ${\tau }_{\mathrm{dep}}\lesssim 1\,\mathrm{Myr}$ , the large amount of gas remaining in several Myr may be compatible with the disk lifetime from observations.

6. Summary

We investigate planet formation in a turbulent disk. Turbulence suppresses the runaway growth of planetesimals. Once the random velocity of planetesimals is comparable to their escape velocity, runaway growth occurs (Kobayashi et al. 2016). The mass-weighted average radius during and after runaway growth is approximated to be that at the onset of runaway growth, ${r}_{\mathrm{rg}}$ . Embryos formed through runaway growth become massive through the accretion of planetesimals with $r\sim {r}_{\mathrm{rg}}$ . However, the stirring by massive embryos induces destructive collisions of planetesimals. The collisional cascade grinds bodies down to 10–100 m. Embryos grow through effective accretion of small bodies, while radial drift reduces small bodies. Eventually, the embryo growth stalls due to the depletion of bodies surrounding embryos via destructive collisions and radial drift. Therefore, the formation and growth of embryos strongly depends on the collisional properties for $r\sim {r}_{\mathrm{rg}}$ , which is controlled by turbulent strength.

We have carried out simulations of the collisional evolution of bodies with collisional strengths (Figure 1) for the formation and growth of embryos in various disks, especially taking into account the stirring by density fluctuation caused by turbulence (Equation (13)). We find the following:

1.
Strong turbulence delays the onset of runaway growth and increases ${r}_{\mathrm{rg}}$ . Embryos forming within disk lifetimes of $t\approx 4\,\mathrm{Myr}$ tend to be large for high α, while the onset of runaway growth is too late to form massive embryos within disk lifetimes for $\alpha \gtrsim {10}^{-3}$ . Cores massive enough are formed from ${r}_{{\rm{g}}}\gtrsim 10\,$ km, corresponding to $\alpha \gtrsim {10}^{-3}$ . However, the formation timescales of such cores are longer than the disk depletion time for α ≳ 10⁻³.
2.
Solid-enhanced disks are preferable for the formation of massive embryos. Such a local enhancement of solids can occur through radial drift of dust aggregates. For weak turbulence α ≲ 10⁻³, embryo masses at t = 4 Myr are similar even in the enhancement of solids. However, embryos may grow within the disk lifetime even for strong turbulence α ≳ 10⁻³. Therefore, embryos are as large as the critical core mass within disk lifetimes in the gas disk more massive than 2MMSN.
3.
The embryo growth depends on the collisional strength of bodies. We investigate collisional evolution for primordial and melted bodies using the model described in Equations (14)–(20), as shown in Figure 1. For weak turbulence disks, embryo growth is independent of collisional strength until M_E ≈ 0.01 M_⊕. However, once collisional fragmentation of bodies with $r\sim {r}_{\mathrm{rg}}$ is effective, the embryo growth strongly depends on collisional strength. The low-mass end of the collisional cascade is important for the efficiency of the accretion of bodies with r = 10–100 m. The collisional strength ${Q}_{{\rm{D}}}^{* }\gtrsim {10}^{7}\,\mathrm{erg}\,{{\rm{g}}}^{-1}$ for r ≲ 1 km, which is satisfied for primordial bodies rather than monolith material, is likely to form a massive core to be a gas giant (see Figures 11 and 13).
4.
We have estimated the timescale of planetary migration during embryo growth. We have taken into account a gas density evolution similar to the accretion disk model. For example, embryos grow to the critical core mass for gas giant formation prior to migration in a disk initially containing 3MMSN gases for a dispersal timescale of 0.5 Myr (see Figure 15). In spite of the short dispersal timescale, the disk mass at the onset of rapid gas accretion remains much higher than a Jupiter mass. Therefore, moderately massive disks with short dispersal timescales are likely to form gas giants without significant migration.

The work is supported by Grants-in-Aid for Scientific Research (26287101, 17K05632, 17H01105, 17H01103) from MEXT of Japan and by JSPS Core-to-Core Program ''International Network of Planetary Sciences''.

Appendix A: Planetary Atmosphere

The atmosphere structure is governed by (e.g., Inaba & Ikoma 2003)

$\begin{eqnarray}&&\displaystyle \frac{{dP}}{{dR}}=-\displaystyle \frac{{{GM}}_{{\rm{E}}}{\rho }_{{\rm{a}}}}{{R}^{2}},\,\end{eqnarray} \tag{ 31 }$

$\begin{eqnarray}\displaystyle \frac{{dT}}{{dR}}=\left\{\begin{array}{ll}-\displaystyle \frac{3}{4\pi {\sigma }_{\mathrm{SB}}}\displaystyle \frac{\kappa {L}_{{\rm{e}}}P}{{T}^{4}} & \mathrm{if}\quad \displaystyle \frac{3\kappa {L}_{{\rm{e}}}P}{4\pi {\sigma }_{\mathrm{SB}}{T}^{4}}\lt 1-\displaystyle \frac{1}{{{\rm{\Gamma }}}_{2}},\\ -\left(1-\displaystyle \frac{1}{{{\rm{\Gamma }}}_{2}}\right)\displaystyle \frac{{GM}{\rho }_{{\rm{a}}}}{{{PR}}^{2}T} & \mathrm{otherwise},\end{array}\right.\end{eqnarray} \tag{ 32 }$

where P is the pressure, T is the temperature, and Γ₂ is the second adiabatic exponent. In this study, we assume that the atmospheric mass M_A is much lower than M_E . However, if ${M}_{{\rm{A}}}/{M}_{{\rm{E}}}\gtrsim 1/3$ , Equation (31) is invalid, and then the rapid gas accretion forms gas giants (e.g., Mizuno 1980).

In Equation (32), the first term on the right-hand side is given by the radiative energy transfer, while the second term is determined by convective transfer. Here, we derive a solution using the radiative term, although we take into account both terms for simulations in Section 3. Dividing Equation (32) by (31) and integrating over P, we have

$\begin{eqnarray}&&{\rho }_{{\rm{a}}}=\displaystyle \frac{16\pi {\sigma }_{\mathrm{SB}}{{GM}}_{{\rm{E}}}{T}^{3}\mu {m}_{{\rm{H}}}}{3k\kappa {L}_{{\rm{e}}}}.\end{eqnarray} \tag{ 33 }$

Inserting Equation (33) into (32) and integrating it, we then obtain

$\begin{eqnarray}&&T=\displaystyle \frac{{{GM}}_{{\rm{E}}}\mu {m}_{{\rm{H}}}}{4{kR}}.\end{eqnarray} \tag{ 34 }$

Substitution of Equation (34) into (33) gives Equation (1).

Appendix B: Collisional Outcome Modeling

The collisional outcome from the collision between bodies with m₁ and m₂ is expressed by m_e and ${\rm{\Psi }}(m,{m}_{1},{m}_{2})$ , which are determined as follows:

$\begin{eqnarray}&&{m}_{{\rm{e}}}=\displaystyle \frac{\phi }{1+\phi }({m}_{1}+{m}_{2}),\end{eqnarray} \tag{ 35 }$

where ϕ is the dimensionless impact energy. Using the impact velocity ${v}_{\mathrm{imp}}$ and ${Q}_{{\rm{D}}}^{* }$ , which is the specific impact energy needed for ejection of half-bodies, $\phi ={m}_{1}{m}_{2}{v}_{\mathrm{imp}}/2{({m}_{1}+{m}_{2})}^{2}{Q}_{{\rm{D}}}^{* }$ .

$\begin{eqnarray}({m}_{1}+{m}_{2}){\rm{\Psi }}(m,{m}_{1},{m}_{2})=\left\{\begin{array}{ll}{m}_{e} & (m\gt {m}_{{\rm{L}}}),\\ {m}_{e}{\left(\displaystyle \frac{m}{{m}_{{\rm{L}}}}\right)}^{2-b} & (m\leqslant {m}_{{\rm{L}}}),\end{array}\right.\end{eqnarray} \tag{ 36 }$

where ϕ is the dimensionless impact energy and ${m}_{{\rm{L}}}$ is the largest mass of fragments produced by a single collision between bodies with m₁ and m₂, given by

$\begin{eqnarray}&&{m}_{{\rm{L}}}=\displaystyle \frac{{\epsilon }_{{\rm{L}}}}{1+\phi }{m}_{{\rm{e}}}=\displaystyle \frac{{\epsilon }_{{\rm{L}}}\phi }{{(1+\phi )}^{2}}({m}_{1}+{m}_{2}),\end{eqnarray} \tag{ 37 }$

b and ${\epsilon }_{{\rm{L}}}$ are constants. We set $b=5/3$ and ${\epsilon }_{{\rm{L}}}=0.2$ . The timescale of collisional cascades is insensitive to values of b and ${\epsilon }_{{\rm{L}}}$ (Kobayashi et al. 2010).

From Planetesimal to Planet in Turbulent Disks. II. Formation of Gas Giant Planets

Article metrics

Permissions

Author e-mails

Author affiliations

Dates

Abstract

1. Introduction

2. Critical Core Mass

3. Model for Planetary Accretion

4. Mass Evolution of Bodies

4.1. Size and Velocity Evolution of Bodies in the Whole Disk

4.2. Size and Velocity Distributions around 5 au

4.3. Turbulent Strength Dependence of Embryo Growth

4.4. Forming Cores with Critical Core Masses

5. Discussion

5.1. Growth Timescale

5.2. Planetary Migration and Gas Dispersal

6. Summary

Appendix A: Planetary Atmosphere

Appendix B: Collisional Outcome Modeling

Footnotes

From Planetesimal to Planet in Turbulent Disks. II. Formation of Gas Giant Planets

Article metrics

Permissions

Share this article

Author e-mails

Author affiliations

Dates

Abstract

1. Introduction

2. Critical Core Mass

3. Model for Planetary Accretion

4. Mass Evolution of Bodies

4.1. Size and Velocity Evolution of Bodies in the Whole Disk

4.2. Size and Velocity Distributions around 5 au

4.3. Turbulent Strength Dependence of Embryo Growth

4.4. Forming Cores with Critical Core Masses

5. Discussion

5.1. Growth Timescale

5.2. Planetary Migration and Gas Dispersal

6. Summary

Appendix A: Planetary Atmosphere

Appendix B: Collisional Outcome Modeling

Footnotes