© ESO 2021

1 Introduction

Pebble accretion is a version of core accretion where planets grow primarily by intercepting a stream of sub-centimeter sized particles (Ormel & Klahr 2010; Lambrechts & Johansen 2012). Previous works have mainly considered its dynamical aspects and have demonstrated that it can provide a rapid channel for growth (Morbidelli & Nesvorny 2012; Chambers 2014), especially when pebbles are settled in the mid-plane and accretion proceeds in a 2D manner (Liu & Ormel 2018; Ormel & Liu 2018). A key feature of pebble accretion is that smaller particles are naturally prevented from reaching the planet once it has grown sufficiently large to perturb the surrounding disk (Morbidelli et al. 2015; Bitsch et al. 2018; Eriksson et al. 2020). This has been suggested as a way of discerning gas and ice giant formation (Lambrechts et al. 2014) and of limiting the growth of super-Earths around low-mass stars (Liu et al. 2019; Chen et al. 2020). Its drag-based mode of capture also introduces an asymmetry that might explain the structural prograde rotation of smaller Solar System objects (Johansen & Lacerda 2010; Visser et al. 2020).

Besides these global effects, the accretion of smaller particles also has important consequences for a planet’s internal structure. So far, most inquiries have focused on compositional changes as the rapid vaporization of pebbles (Love & Brownlee 1991; McAuliffe & Christou 2006) naturally leads to the deposition of significant amounts of vapor, which can be stored in the deep interior (Iaroslavitz & Podolak 2007; Lozovsky et al. 2017). This naturally limits the size of central cores (Alibert 2017; Brouwers et al. 2018) and results in a dense, polluted interior (Iaroslavitz & Podolak 2007; Venturini et al. 2016; Lozovsky et al. 2017; Bodenheimer et al. 2018). The traditional rigid core-envelope structure of these planets is replaced by the natural emergence of a compositional gradient (Valletta & Helled 2020; Ormel et al. 2021; Vazan et al. 2020). This seems to be in line with Juno’s measurements of Jupiter’s gravitational moments of inertia, which imply such a dilute core structure (Wahl et al. 2017; Debras & Chabrier 2019; Debras et al. 2021).

Before pebbles get to these inner regions, however, they must first sediment through the tenuous outer layers, where their combined surface area can contribute to the opacity. This influence of pebbles on the outer envelope is still largely unexplored but has important thermodynamic consequences. The more opaque an envelope, the less heat is able to escape it and the less gas it can gravitationally bind. By regulating the pace of cooling, the opacity is a key variable that directly affects the outcome of planet formation. The quantitative importance of the opacity is often understated in current formation models, where a lack of physically motivated values is still a serious issue. It is common practice to adopt conveniently high ISM-like opacities or values that are arbitrarily scaled down, with important consequences for the outcome of these models. One of the reasons that envelope opacity values are currently poorly constrained is that previous grain growth models focused specifically on the formation of Jupiter at 5 AU. In the most detailed of these, Podolak (2003), Movshovitz & Podolak (2008) and Movshovitz et al. (2010) calculated grain growth with numerical models that solve the Smoluchowski equation under the assumption that grains stick when they collide. They found that grain growth effectively reduces the dust opacity as a function of depth, with the opacity declining by around three orders of magnitude from the outer layers to the inner radiative-convective boundary (RCB). Simpler analytical (Mordasini 2014) and numerical (Ormel 2014) models with a single characteristic grain size at a given height have been able to replicate this main result. However, because proto-planetary disks are thermodynamically very different across distances and masses, it is not clear that these results can be applied generally to planets of varying sizes throughout the disk.

In addition, whether a planet is predominantly accreting sub-centimeter pebbles or 100 km planetesimals is clearly important for the abundance of grains in the outer envelope. The assumption of grain growth models that planetesimals deposit large amounts of small grains in the upper layers is in conflict with impact simulations, which predict that most of their mass is released close to the planet’s central cores, far below the RCB (Brouwers et al. 2018; Valletta & Helled 2019, 2020). In contrast, Ali-Dib & Thompson (2020) and Johansen & Nordlund (2020) have recently shown that pebbles are very susceptible to erosion by small dust grains when they enter planetary envelopes, especially when they are accelerated further by convective cells. It was also pointed out by Ormel et al. (2021) that even without any size evolution, milimeter-sized pebbles can contribute a significant opacity to planetary envelopes if they accrete at a sufficiently high rate.

The total envelope opacity is a sum of the contributions of solids and gas (κ = κ_s+ κ_gas). In this study, we formulate a model for the contribution of solids, which often dominates over the gas during accretion. In Sect. 2, we develop a physical opacity model that accounts for the main processes that influence the evolution of solids in planetary envelopes. We consider a variation of the single-size approximation, where we model the populations of both small dust (κ_d) and larger pebbles (κ_peb). The characteristic size of the pebbles in our model is either determined by their growth from sticking collisions (coalescence) at low velocities or limited by fragmentation and erosion at high velocities, depending on the local thermodynamic conditions. When pebbles travel at speeds below the erosion velocity, they sweep up dust grains without experiencing mass loss (Krijt et al. 2015; Schräpler et al. 2018). We assume a constant dust size and calculate its steady state abundance between dust sweep-up and production, which can occur as a by-product of collisions between pebbles or via erosion. In Sect. 3, we then formulate simple, physically motivated analytical expressions for the opacity from solids, similar to what was done by Mordasini (2014). We apply our opacity model across a wide parameter space to explore and map out the resulting opacity as a function of distance, planet mass and accretion rate in Sect. 4. We then present the implications of these findings for the formation of gas and ice giants in Sect. 5. We discuss our model in relation to contemporary works, list potential model improvements and contrast the scenarios of pebble and planetesimal accretion in Sect. 6. Finally, we conclude our work in Sect. 7. In Appendix A, we derive an expression for the non-isothermal critical metal mass that we use to study the onset of runaway growth. In Appendix B, we verify our constant dust size by comparing how far dust grains can travel before they collide. In the final Appendix C, we vary the limiting velocity, the parameter in our model that represents the onset of fragmentation as well as erosion by micron-sized grains.

2 Model description

2.1 Two-population approach

Growing proto-planets are supplied with a flux of solid particles that enter their envelopes by two distinct mechanisms. The first is the drag-assisted capture of sub-centimeter pebbles that drift toward the central star (Ormel & Klahr 2010; Lambrechts & Johansen 2012), supplemented by the accretion of larger, km-sized planetesimals (Alibert et al. 2018; Guilera et al. 2020). The second source is the population of small, micron-sized dust grains that are coupled to the nebular gas and that enter the planet’s envelope when it is massive enough to bind one.

In two recent works, Ali-Dib & Thompson (2020) and Johansen & Nordlund (2020) have shown that the interaction between dust and pebbles is crucial in evaluating their relative abundances as high-velocity collisions between the two can significantly erode the incoming pebbles and add to the dust population (Krijt et al. 2015; Schräpler et al. 2018). This erosion is a runaway process and is only halted at the pebble size, where collisions switch to sticking and the dust abundance begins to drop. In our model, we incorporate both dust and pebbles in a simple manner, where each is represented with a single characteristic size. The pebble size is determined locally by collisions within the population (growth or fragmentation), as well as collisions with the dust (sweep-up and erosion). The size of individual dust particles is a parameter in our model, which we keep constant as motivated by the ongoing production and sweep-up of these smaller grains in the envelopes (see Appendix B). Figure 1 provides a qualitative overview of our model, which we work out in the next subsections.

2.2 Sedimentation of solids

Particles referred to as dust and pebbles are small enough to be effectively slowed down by gas drag and sink at speeds close to their terminal velocities (v_fall) relative to the local medium. We adopted a simple, continuous two-regime approximation to the drag force, where the transition from free molecular (Epstein) to continuum flow (Stokes) is set at the canonical value l_mfp ∕R_s = 4∕9 and the terminalvelocity is given by (Weidenschilling 1977a): $$v_{\text{fall}}=\frac{g{R}_{\text{s}}{\rho }_{•}}{{\rho }_{\text{g}}{v}_{\text{th}}}\text{max}\left(\frac{4R_{\text{s}}}{9l_{\text{mfp}}},1\right),$$(1)

where g and ρ_g are the local gravitational acceleration and gas density at a distance from the planet’s center r and the subscript s refers to either population of solids (pebbles with size R_peb and dust with size R_d), which share the same material density ρ_∙. We take the standard ideal gas expressions for the thermal velocity v_th and the mean free path of molecules l_mfp. The second component to the sedimentation velocity of solids is given by the downward flow of gas. Because we are mostly interested in the planet’s upper layers that only contain a fraction of the envelope mass, we can approximate this gas velocity (v_gas) from mass conservation of the gas accretion rate (Ṁ_xy) as: $$v_{\text{gas}}=\frac{{\dot{M}}_{\text{xy}}}{4\pi r^2{\rho }_{\text{g}}}.$$(2)

When the solids are small enough, as is the case for the micron-sized dust particles in our model, the local downward flow of gas dominates their total sedimentation velocity. This was described as the advection regime in the work by Mordasini (2014). For larger particles, the free-fall term instead dominates and, generally, the total sedimentation velocity relative to the planet’s core is a sum of the two: $$v_{\text{sed}}=v_{\text{fall}}+v_{\text{gas}}.$$(3)

We consider quasi-static envelopes, where the grain transport timescale is assumed to be short relative to the timescale on which the envelope’s thermodynamic conditions change. This is a good assumption in the outer layers, provided that the total envelope mass is dominated by the planet’s polluted interior. If the solids additionally enter in a spherically symmetric manner, their volume density is given by: $${\rho }_{\text{s}}=\frac{{\dot{M}}_{\text{s}}}{4\pi r^2{v}_{\text{sed}}},$$(4)

which represents radial mass conservation. We do not include large-scale convective motions (Popovas et al. 2018, 2019) in our model, which could locally transport clumps of solids and alter their collision velocities (Ormel & Cuzzi 2007).

Fig. 1 Schematic diagram of our two-population opacity model. Depending on the velocities of sedimenting pebbles, they can experience growth (coalescence), fragmentation, or erosion. In the growth-limited regime, the pebble size in the interior is regulated by the collision and sedimentation timescales (Rpeb = Rcoal). In the velocity-limited regime, fragmentation and/or erosion restrict the pebble size below the growth potential (Rpeb = Rvlim < Rcoal). In our model, dust grains of constant size are produced in collisions between pebbles or in high-velocity pebble-dust encounters (erosion), while they are lost by sticking encounters with slow-moving pebbles (dust sweep-up), leading to a local steady state in their abundance.

2.3 Dust-pebble collisions

2.3.1 Erosion at high velocities

The sizes of pebbles that enter planetary envelopes are determined by their growth, fragmentation and radial drift in the proto-planetary disk (i.e., Liu & Ji 2020; Drazkowska et al. 2021). The observational evidence points toward typical sizes between 100 μm − 1 cm (i.e., Birnstiel et al. 2010; Kataoka et al. 2016; Hull et al. 2018; Carrasco-González et al. 2019; Ohashi et al. 2020; Tazzari et al. 2021b), which generally decrease with distance to the central star (Tazzari et al. 2021a). During their capture by the planet, pebbles are subject to additional gravitational acceleration and begin to collide with dust at increasing velocities. It was pointed out by Ali-Dib & Thompson (2020) and Johansen & Nordlund (2020) that if this velocity crosses the threshold for erosion, the dust particles begin to chip off pebble material. This process was measured by Schräpler et al. (2018) to begin at: $$v_{\text{erosion}}=2.4\text{m}{\text{s}}^{-1}{\left(\frac{R_{\text{d}}}{1\text{upmum}}\right)}^{\frac{1}{1.62}},$$(5)

with an erosive mass loss per collision Δm_erosion of: $$\frac{\Delta m_{\text{erosion}}}{m_{\text{d}}}=4.3\left(\frac{v}{10\text{m}{\text{s}}^{-1}}\right){\left(\frac{R_{\text{d}}}{1\text{upmum}}\right)}^{-0.62},$$(6)

which requires a mass ratio below ~10⁻² (Krijt et al. 2015). The proportionality of the erosion efficiency to the collision velocity was also found in numerical investigations (Seizinger et al. 2013; Planes et al. 2017). The fragments that are produced are generally similar or smaller in size than the projectiles, which in our model are dust grains with mass m_d. Because theproduced fragments add to the dust population, the erosion of pebbles is a runaway process up to a characteristic size, at which point they slow down sufficiently due to their increased susceptibility to gas drag. In our model with linear downward sedimentation, this is set by the condition that v_{peb, fall} = v_erosion (assuming v_fall,d ≪ v_{fall, peb}): $$R_{\text{erosion}}=\{\begin{array}{l}\frac{{\rho }_{\text{g}}{v}_{\text{th}}{v}_{\text{erosion}}}{g{\rho }_{•}}(\text{Epstein})\hfill \\ \hfill \\ \frac{3}{2}{\left(\frac{{\rho }_{\text{g}}{v}_{\text{th}}{l}_{\text{mfp}}{v}_{\text{erosion}}}{g{\rho }_{•}}\right)}^{\frac{1}{2}}(\text{Stokes}).\hfill \end{array}$$(7)

In order to estimate how rapidly pebbles erode, we formulate the radial equations of mass transfer and size evolution. The typical distance that pebbles travel before they encounter a dust particle is given by (assuming R_d ≪ R_peb and v_{fall, d} ≪ v_{fall, peb}): $$l_{\text{sweep},\text{peb}}=\frac{m_{\text{d}}}{{\rho }_{\text{d}}{\sigma }_{\text{peb}}}\frac{v_{\text{sed},\text{peb}}}{v_{\text{fall},\text{peb}}},$$(8)

where ${\sigma }_{\text{peb}}=\pi R_{\text{peb}}^2$ is the pebble cross section. The equations for mass transfer and shrinkage follow as: $$\begin{array}{ll}\frac{\text{d}{\dot{M}}_{\text{peb}}}{\text{d}r}{|}_{\text{erosion}}\hfill & =\frac{{\dot{M}}_{\text{peb}}}{l_{\text{sweep},\text{peb}}}\frac{\Delta m_{\text{erosion}}}{m_{\text{peb}}},\hfill \\ \frac{\text{d}{\dot{M}}_{\text{d}}}{\text{d}r}{|}_{\text{erosion}}\hfill & =-\frac{\text{d}{\dot{M}}_{\text{peb}}}{\text{d}r},\hfill \\ \frac{\text{d}{m}_{\text{peb}}}{\text{d}r}{|}_{\text{erosion}}\hfill & =\frac{\Delta m_{\text{erosion}}}{l_{\text{sweep},\text{peb}}},\hfill \end{array}$$

where Ṁ_peb is the pebble mass flux. In this simple model, the increase in dust abundance from pebble erosion is essentially exponential and the pebbles shrink very rapidly, even without assuming any additional acceleration from convective eddies. As an example, weintegrate Eqs. (9a)–(9c) from the Bondi radius inward for a nominal planet of 5 M_⊕ at 5 AU with a dust disk-to-gas ratio of 10⁻³ (see Table 1 for all default parameters). The resulting size evolution curves are shown in Fig. 2 for three different pebble entry sizes. In all the runs, the pebbles quickly erode down to a characteristic size, at which point their velocity decreases below the threshold for erosion and they are safe. In this manner, erosion essentially cancels out the upper range of the initial pebble size distribution to more uniform pebble entry sizes, limited by their erosion at the Bondi radius.

Fig. 2 Erosive mass loss of pebbles that enter the envelope of a growing planet at its Bondi radius (rB). The parameters of the planet are shown in Table 1. The three lines correspond to initial pebble entry sizes of 0.1 cm (dotted), 1 cm (dashed) and 10 cm (solid). Due to the positive feedback of erosive mass-loss on the dust abundance, the pebbles rapidly converge to the size Rerosion (Eq. (7)) where collisions with dust grains switch from causing mass loss via erosion to sticking (sweep-up).

2.3.2 Sweep-up at low velocities

If the pebbles are reduced in size below R_erosion, the relative velocities between the dust and the pebbles drop below v_erosion and dust can stick to the pebbles without dislodging additional material. In this sweep-up regime, Eqs. (9a)–(9c) still apply with Δm = −m_d: $$\frac{\text{d}{\dot{M}}_{\text{d}}}{\text{d}r}{|}_{\text{sweep-up}}=\frac{{\dot{M}}_{\text{peb}}}{l_{\text{sweep},\text{peb}}}\frac{m_{\text{d}}}{m_{\text{peb}}},$$(10)

Now, the presence of pebbles acts to reduce the dust abundance, rather than increase it. In Sect. 3.2, we will discuss the scenario where the production of dust via erosion and pebble-pebble collisions is in equilibrium with the sweep-up of dust, providing a steady state dust abundance.

2.4 Pebble-pebble collisions

When considering straight downward sedimentation, collisions between two pebbles can occur either through Brownian motion or by differential settling (coalescence). The timescale of the latter process decreases with increasing particle size, such that coalescence easily dominates their mutual collision rate when the pebbles are larger than several microns (Boss 1998; Nayakshin 2010; Mordasini 2014; Ormel 2014). As identical particles settle with the same speed, coalescence within a population formally arises from the parameterized width in the size distribution. In our model, we follow the works by Rossow (1978) and Mordasini (2014), who set the typical mass ratio of collisions to 0.5. Although a simplification, Krijt et al. (2016) and Sato et al. (2016) show that this is a good approximation to the average mass ratio found in coagulation simulations that treat a complete size distribution. When combined with Eq. (3), it leads to collisional velocities in the range of 0.21 < v_col∕v_fall < 0.37 depending on the drag regime. In our model, we neglect this order unity difference in favor of a continuous collision velocity across drag regimes and approximate it as x_R ≡ v_col∕v_fall ≃ 1∕3.

2.4.1 Growth by coalescence at low velocities

When the collision velocities are sufficiently low, contact between pebbles leads to sticking and growth (Dominik & Tielens 1997). To work out the rate of this growth, the typical travel distance between pebble collisions (l_col,peb) can be related to their local size, abundance, and velocity as: $$l_{\text{col},\text{peb}}=\frac{m_{\text{peb}}}{{\rho }_{\text{peb}}{\sigma }_{\text{peb}}}\frac{v_{\text{sed},\text{peb}}}{v_{\text{col},\text{peb}}}.$$(11)

Equation (11) scales positively with pebble size, which indicates that smaller particles collide more often. By equating l_{col, peb} to the scale height H = k_bT_g∕(μ_gg), which represents the characteristic length scale in an envelope, it is possible to approximate the maximum sizes (R_coal) to which solidscan grow as they sediment: $$R_{\text{coal}}=\{\begin{array}{l}{\left(\frac{3x_{\text{R}}H{\dot{M}}_{\text{peb}}{v}_{\text{th}}{\rho }_{\text{g}}}{16\pi G{M}_{\text{p}}{\rho }_{•}^2}\right)}^{\frac{1}{2}}(\text{Epstein})\hfill \\ \hfill \\ \frac{3}{4}{\left(\frac{x_{\text{R}}H{\dot{M}}_{\text{peb}}{v}_{\text{th}}{l}_{\text{mfp}}{\rho }_{\text{g}}}{\pi G{M}_{\text{p}}{\rho }_{•}^2}\right)}^{\frac{1}{3}}(\text{Stokes}),\hfill \end{array}$$(12)

where we assumed that their total velocities are dominated by the terminal component rather than by inward gas flow. In the regime where growth by coalescence provides the limit to their size, the pebbles are typically large enough (R_peb ≳ 100 μm) that this is indeed a good assumption. The scaling of Eq. (12) with the pebble accretion rate shows that a greater mass flux of pebbles leads to more collisions and faster growth. Its dependence on the planet’s mass is more complicated. More massive planets have a stronger gravitational pull, which increases the distance at which they can bind gas and, therefore, extends the envelope (r_B ∝ M_p). This lowers the gravitational acceleration at the planet’s outer boundary as a function of mass ($g\propto M_{\text{p}}/r_{\text{B}}^2\propto M_{\text{p}}^{-1}$). Taken together, the collision distance scales as l_col,peb ∝ M_p, which balances with the increasing scale height H ∝ M_p and leads to a constant value of R_coal across planetary masses.

2.4.2 Collisional fragmentation at high velocities

At higher velocities, collisions between pebbles lead to bouncing or even fragmentation. Collision experiments typically find that silicate particles begin to fragment if their collision velocity is larger than 1 m s⁻¹ (Wurm et al. 2005; Schäfer et al. 2007; Güttler et al. 2010), and up to 10 m s⁻¹ for water ice grains (Gundlach & Blum 2015; Musiolik & Wurm 2019). This fragmentation criterion can be combined with the velocity for the onset of erosion (see Eq. (5)) to yield a limit on the terminal velocities of pebbles: $$v_{\text{lim}}=\text{min}\left(v_{\text{erosion}},v_{\text{frag}}/x_{\text{R}}\right).$$(13)

Depending on the dust size and the pebble’s material properties, either fragmentation or erosion can be the limiting factor to pebble growth. For simplicity, we take a default dust size of 1 μm and a fragmentation velocity of 0.8 m s⁻¹, for which the two velocity limits on the pebble’s terminal velocity are both equal to v_fall,peb = v_lim = 2.4 m s⁻¹. In Appendix C, we vary the limiting velocity within a wider range. We refer to the velocity-limited pebble size as R_vlim and it can be found from Eq. (1) as: $$R_{\text{vlim}}=\{\begin{array}{l}\frac{{\rho }_{\text{g}}{v}_{\text{th}}{v}_{\text{lim}}}{g{\rho }_{•}}(\text{Epstein})\hfill \\ \hfill \\ \frac{3}{2}{\left(\frac{{\rho }_{\text{g}}{v}_{\text{th}}{l}_{\text{mfp}}{v}_{\text{lim}}}{g{\rho }_{•}}\right)}^{\frac{1}{2}}(\text{Stokes}).\hfill \end{array}$$(14)

The dependence of R_vlim at the Bondi radius is visualized in Fig. 3, with comparison to the limit set by the growth rate. In contrast to R_coal, the value of R_vlim scales positively with the planet’s mass. Hence, the sizes of pebbles in small envelopes are typically velocity-limited via erosion or fragmentation, whereas pebbles in more massive envelopes are typically only limited by their rate of growth. Close to the central star, the disk is relatively dense and pebbles of the same size sediment more slowly. Consequently, pebbles erode or fragment down to the smallest sizes when they accrete onto low-mass planets that reside in the outer disk.

Besides limiting the sizes to which pebbles can grow, collisions between pebbles can also potentially result in production of small dust as collisional by-products. In our two-component model, we include this with a fractional dust production efficiency F, which can theoretically be between 0 (no dust production) and 1 (all growth beyond R_vlim is turned into dust). Physically, F is a measure for the number of pebble-pebble collisions needed to completely grind down a pebble. The rate at which this dust is produced is then given by the pebble collision rate as: $$\begin{array}{ll}\frac{\text{d}{\dot{M}}_{\text{d}}}{\text{d}r}{|}_{\text{frag}}\hfill & =-F\frac{{\dot{M}}_{\text{peb}}}{l_{\text{col},\text{peb}}},\hfill \\ \frac{\text{d}{m}_{\text{peb}}}{\text{d}r}{|}_{\text{frag}}\hfill & =F\frac{m_{\text{peb}}}{l_{\text{col},\text{peb}}}.\hfill \end{array}$$

We will discuss reasonable values for the parameter F in Sect. 3.3.

Fig. 3 Characteristic pebble size limits at the Bondi radius, plotted with the parameters of Table 1. The colors indicate the pebble size where vfall,peb = vlim = 2.4 m s−1, the common limit set via erosion with micron-sized grains and fragmentation. The white dashed line separates the Stokes and Epstein regimes while the dotted line indicates the regions where either growth rate or erosion/fragmentation limits the local pebble size. The plus sign marks the default planet mass (5 M⊕) and distance (5 AU).

3 Analytical opacity expressions

In this section, we formulate analytical expressions for the opacity contributions from pebbles and dust. For the pebbles, these follow from the previously identified pebble size limits in different regimes (R_coal, R_vlim). The basic expression for the Rosseland mean opacity from solids (κ_s) is: $${\kappa }_{\text{s}}=\frac{3Q_{\text{eff}}{\rho }_{\text{s}}}{4{\rho }_{•}{R}_{\text{s}}{\rho }_{\text{g}}},$$(16)

where the extinction efficiency Q_eff can be approximated as $Q_{\text{eff}}\simeq \text{min}(0.6\pi R_{\text{s}}/{\lambda }_{\text{peak}},2)$, with the peak wavelength of the emitted photons by the local gas λ_peak(cm) = 0.290∕T_g (Wien’s law). Laboratory experiments provide a more detailed temperature scaling of Q_eff, as well as a factor of ~2 between different species (Movshovitz & Podolak 2008; Bitsch & Savvidou 2021). We do not include these additional details here but do account for the most important opacity difference between species, which is that solid particles are only present in layers of planetary envelopes that are sufficiently cool for the solids to escape sublimation. In the case of silicates, which we consider as the default composition of pebbles and dust, this threshold is positioned deep inside planetary envelopes around ~2500 K.

3.1 Pebble opacity expression

The first opacity regime is the growth-limited regime, which applies when pebbles are traveling relatively slowly through planetary envelopes and experience sticking collisions, rather than erosion or fragmentation. In this case, their size can be approximated by R_coal (Eq. (12)) and their opacity contribution is independent of their accretion rate. It follows the same equation in both drag regimes: $${\kappa }_{\text{coal}}=\frac{Q_{\text{eff}}}{x_{\text{R}}H{\rho }_{\text{g}}},$$(17)

which was the main finding of Mordasini (2014). As indicated by Fig. 3 at the Bondi radius, this growth-limited regime is mainly applicable in the envelopes of more massive planets, whose larger envelopes allow pebbles to sediment more slowly. The required planetary size to enter this regime at the same accretion rate is an increasing function of distance from the central star, as pebbles in the tenuous outer disk sediment faster than those closer in. When the conditions are such that the terminal velocity of pebbles exceeds v_lim, the size to which solids can grow becomes limited to R_vlim (Eq. (14)) by either fragmentation or erosion. In this regime, the pebble opacity is instead given by: $$\begin{array}{ll}{\kappa }_{\text{vlim}}\hfill & =\frac{3Q_{\text{eff}}{\rho }_{\text{peb}}}{4{\rho }_{•}{R}_{\text{vlim}}{\rho }_{\text{g}}}\hfill \\ \hfill & =\{\begin{array}{l}\frac{3Q_{\text{eff}}{\dot{M}}_{\text{peb}}g}{16\pi r^2{\rho }_{\text{g}}^2{v}_{\text{th}}{v}_{\text{lim}}(v_{\text{lim}}+v_{\text{gas}})}(\text{Epstein})\hfill \\ \hfill \\ \frac{Q{\dot{M}}_{\text{peb}}}{8\pi r^2(v_{\text{lim}}+v_{\text{gas}})}{\left(\frac{g}{{\rho }_{•}{\rho }_{\text{g}}^3{v}_{\text{th}}{v}_{\text{lim}}{l}_{\text{mfp}}}\right)}^{\frac{1}{2}}(\text{Stokes}).\hfill \end{array}\hfill \end{array}$$

In combination, the pebble opacity can be approximated analytically from the expressions above by determining the appropriate regime based on the characteristic sizes, which have to be evaluated locally in the envelope: $${\kappa }_{\text{peb}}=\{\begin{array}{l}{\kappa }_{\text{vlim}}\text{if}{R}_{\text{vlim}}<R_{\text{coal}}\hfill \\ \hfill \\ {\kappa }_{\text{coal}}\text{if}{R}_{\text{vlim}}>R_{\text{coal}}.\hfill \end{array}$$(19)

For completeness, we also include the simplest scenario where solids enter a planetary envelope with a constant size that remains unchanged during their sedimentation. While we follow the more physically motivated expressions from Eq. (19) in the rest of this work, the scenario of a constant pebble size might be applicable if pebbles are both sufficiently small to escape fragmentation and bounce rather than stick upon contact, as is the case if the pebbles are modeled as small molten chondrules, rather than dust agglomerates. Under this assumption, the opacity follows from Eqs. (16), (3) as (see also the work by Ormel et al. 2021): $${\kappa }_{\text{peb},\text{cst}}=\frac{3Q_{\text{eff}}{\dot{M}}_{\text{peb}}}{16\pi r^2{\rho }_{•}{\rho }_{\text{g}}{R}_{\text{peb}}{v}_{\text{sed},\text{peb}}}.$$(20)

3.2 Steady state between dust replenishment and sweep-up

In previous grain growth models (Ormel 2014; Mordasini 2014), a large influx of small dust grains in the envelope’s outer layers was found to quickly diminish due to mutual sticking collisions. Hence, even when a significant fraction of the pebble mass is transferred to dust grains upon entry by efficient erosion (Sect. 2.3.1), the smallest particles soon disappear from the envelope, growing to the same limiting sizes as pebbles that we discussed in the previous section. The diminishing of the small grains abundance is further hastened by the sweep-up of pebbles below their erosion size. In order to maintain a dust population, therefore, it must be supplied by either continued erosion of pebbles that grow beyond the erosion limit, or by fragmentation. If that happens, the dust abundance has both a source and a sink term, generating a steady state when they are equal and opposite: $$\frac{\text{d}{\dot{M}}_{\text{d}}}{\text{d}r}{|}_{\text{sweep-up}}=-\frac{\text{d}{\dot{M}}_{\text{d}}}{\text{d}r}{|}_{\text{frag}},$$(21)

from which their respective volume densities follow as $$\frac{{\rho }_{\text{d}}}{{\rho }_{\text{peb}}}=F{x}_{\text{R}}.$$(22)

Equation (22) implies that in a steady state between fragmentation/erosion and sweep-up, most of the radial mass flux is generated by the larger sedimenting pebbles. But if F is near unity, the volume density of grains can nevertheless be comparable to that of the pebbles due to their slower sedimentation.

3.3 Dust opacity in steady state

The dust opacity in steady state between dust production and sweep-up is given by: $$\begin{array}{ll}{\kappa }_{\text{d}}\hfill & =\frac{3Q_{\text{eff},\text{d}}{\rho }_{\text{d}}}{4{\rho }_{•}{R}_{\text{d}}{\rho }_{\text{g}}}\hfill \\ \hfill & ={\kappa }_{\text{peb}}\frac{Q_{\text{eff},\text{d}}}{Q_{\text{eff},\text{peb}}}\frac{R_{\text{peb}}}{R_{\text{d}}}F{x}_{\text{R}}.\hfill \end{array}$$

Equation (23b) is proportional to the pebble opacity, only differing from its trends due to the additional dependence on the pebble size. Because the opacity in planetary envelopes is generally far more variable than the pebble size, the opacity from dust in steady state generally follows a very similar trend to that of the pebbles. This makes it possible to evaluate the dust opacity, even though the parameter F is largely unconstrained. In principle, it can vary between 0−1, with a value around unity more appropriate in the erosion-limited regime where any pebble growth beyond the erosion limit is converted into dust. Lower values, where collisions convert a smaller fraction of their mass into dust, are likely more appropriate in cases where either fragmentation or growth (both pebble-pebble interactions) limits the pebble size. We show in the next section that with our default value of F = 0.1 the dust opacity is typically comparable to the pebble contribution. The biggest positive change with respect to the default model, therefore, occurs for F = 1. For lower F,the opacity is not affected much because the pebble opacity then dominates and remains unchanged.

4 Envelope opacity trends

In this section, we apply our model for the opacity of solids to a broad parameter space in order to investigate the its trends in envelopes of planets throughout the disk. We first detail our envelope model and then evaluate the opacity of dust and pebbles as a function of depth (Figs. 4–6) for a range of planetary masses, distances and pebble accretion rates. After that, we consider the opacity at the boundary between radiative and convective zones to visualize the same trends in a single graph (Fig. 7).

Fig. 4 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). The different lines indicate a range of planet masses, with the arrow indicating a logarithmic progression from 0.5− 20 M⊕. The triangles indicate the location of the RCB, while the stars indicate the depth where the ambient temperature exceeds the sublimation temperature (2500 K) and the opacity from solids vanishes. The colors in the top panel show the ratio of the pebble size in the model relative to Rvlim (Eq. (14)). The colors in the two lower panels show the relative value of the indicated opacity to the total opacity.

4.1 Envelope structure

We focus our modeling efforts on the outer envelope down to the polluted region, ulterior to which no significant sublimation occurs and the opacity of solids is relevant. We refer to an accompanying paper by Ormel et al. (2021) and independent studies by Bodenheimer et al. (2018) and Valletta & Helled (2020) for detailed numerical models of polluted envelope interiors. The structure of the outer envelope is simple by contrast as it shares its gaseous composition with the surrounding disk and is unaffected by self-gravity prior to the onset of runaway accretion. In quasi-hydrostatic equilibrium, its structure equations read: $$\begin{array}{ll}\frac{\partial m}{\partial r}\hfill & =4\pi r^2{\rho }_{\text{g}},\hfill \\ \frac{\partial P_{\text{g}}}{\partial r}\hfill & =-\frac{G{M}_{\text{p}}{\rho }_{\text{g}}}{r^2},\hfill \\ \frac{\partial T_{\text{g}}}{\partial r}\hfill & =\frac{\partial P_{\text{g}}}{\partial r}\frac{T_{\text{g}}}{P_{\text{g}}}\text{min}({\nabla }_{\text{rad}},{\nabla }_{\text{ad}}),\hfill \end{array}$$

where $${\nabla }_{\text{rad}}=\frac{3\kappa L{P}_{\text{g}}}{64\pi \overline{\sigma }G{M}_{\text{p}}{T}_{\text{g}}^4}$$(25)

is the radiative temperature gradient, which contains the gravitational and Stefan-Boltzmann constants ($G,\overline{\sigma }$) and is a function of the luminosity L and the total Rosseland mean opacity κ. In our model, we assume a constant (global) accretion luminosity equal to L = GM_cṀ_peb∕r_c as is commonly done. However, we note the caveat issued by Ormel et al. (2021), that this term can in reality vary substantially, as refractory material can be absorbed before sinking to the core and a significant portion of the gravitational energy can be processed in envelope mixing or used to heat the surrounding gas.

The total envelope opacity is defined as the sum of the gas, dust, and pebble contributions (κ = κ_gas + κ_d + κ_peb). The gas opacity is often taken from lookup tables (i.e., Freedman et al. 2008, 2014) that contain contributions from different gaseous species. We not not model this contribution here but choose to take a simple reference value of κ_gas as an analytical scaling of the molecular opacity from Bell & Lin (1994): $${\kappa }_{\text{mol}}={10}^{-8}{\rho }_{\text{g}}^{\frac{2}{3}}{T}_{\text{g}}^3{\text{cm}}^2{\text{g}}^{-1},$$(26)

which was also used in the first two papers in this series. A simple molecular opacity scaling has the advantage that trends in the opacity from solids, which we seek to characterize here, can be more easily isolated in the results.

We integrate the envelope structure Eqs. (24a)–(24c) from the planet’s outer edge, which is taken as the minimum between the Hill and Bondi radii. The outer boundary conditions are equal to the local disk environment, for which we assume the simple temperature and density relations from the Minimum Mass Solar Nebula (see Table 1) (Weidenschilling 1977b; Hayashi 1981). The total opacity is explicitly included in the integration of the structure equations, and no iteration of the model is needed. We use the ideal gas equation to relate the local density to the pressure and temperature, which is justified in the upper layers of the atmosphere that are most important for the opacity calculation.

Fig. 5 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). This figure is the same as Fig. 4 but now the mass is fixed at the default 5 M⊕ and the planet’s distances from the star are varied logarithmically from 0.1−30 AU.

Fig. 6 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). This figure is the same as Fig. 4 but now the mass is fixed at the default 5 M⊕ and the pebble accretion rates are varied logarithmically from 10−7−10−4 M⊕ yr−1.

4.2 Trend with planet mass

In Fig. 4, we plot the pebble sizes (panel a) along with their opacity (b), and the dust opacity (c) of our steady state model for a range of masses that we vary logarithmically between 0.5 and 20 M_⊕. The other parameters in the runs share are the default values from Table 1. We first discuss the top panel, which shows two separate trends depending on the process that limits the pebble size. In the velocity-limited regime (R_peb = R_vlim, yellow color range), the pebble size is limited by either the erosion with small dust grains or by collisions with other pebbles. Because the Bondi radius is positioned closer to the core in low-mass planets, pebbles of the same size are able to sediment faster ($g\propto M_{\text{p}}^{-1}$) and the same terminal velocity limit leads larger pebbles in more massive planets. Once the planet becomes massive enough that R_vlim > R_coal, the pebble size ceases to be velocity-limited and instead becomes only limited by the rate of growth. In this growth regime, the pebble size is nearly invariant with a further increase in the planet mass (see Sect. 2.4.1).

Regardless of the regime, we find that the opacity from solids declines steeply when planets grow more massive. This can be explained with two competing processes. As planets accrete material, the surface area of their Bondi sphere increases quadratically with mass and the same pebbles obscure a smaller fraction of the envelope. Part of this trend is compensated by the slower sedimentation of pebbles in more massive planets but the net result is that the pebble volume density scales as ${\rho }_{\text{peb}}\propto M_{\text{p}}^{-1}$ at the Bondi radius. If the pebble size is instead limited via erosion and fragmentation, in which regime their sizes also scale positively with planet mass, the opacity from solids is an even steeper declining function. The pebble and dust opacity trends are almost identical because they are nearly proportional to one another (Eq. (23b)).

The difference in opacity between the low- and high-mass planets begins at around three orders of magnitude at the Bondi radius and these differences increase further with depth. Most of this increase is due to the transition from Epstein to Stokes drag that occurs in the more massive envelopes, where larger pebbles enter a denser medium. With this, the difference at the Bondi radius is extended further, to over six orders of magnitude between the plotted values (0.5−20 M_⊕) at the RCB.

4.3 Trend with orbital separation

Next, we evaluate the importance of the planet’s distance to the central star in Fig. 5, where we plot five lines that indicate semimajor axes between 0.1 and 30 AU. The distance to the central star determines the local thermodynamic conditions of the surrounding disk. In our model, we use the Minimum Mass Solar Nebula, where the temperature and density scale as T_disk ∝ d^−1∕2, ρ_disk ∝ d^−11∕4 (Weidenschilling 1977b; Hayashi 1981). While the local conditions will also generally differ between different disk types and models, the general trend of lower densities and temperatures in the outer regions that receive less light is universal. This trend is important, because a more tenuous medium allows pebbles of the same size to sediment faster. The velocity-limited sizes, therefore, decrease with distance from the star and pebbles experience fewer collisions, reducing their rate of growth as well. This trend is reflected in Fig. 5, where the opacity of solids at the Bondi radius is seen to increase with orbital separation, regardless of the regime.

The opacity trend with distance becomes more complex at greater depth, where the opacity of solids at the RCB is seen to converge to similar values for the parameters in Fig. 5. In these runs, the pebble sizes of planets in the outer disk increase with depth as the surrounding medium becomes denser and the same velocity limit allows for larger pebbles. As a result, the opacities also decline with depth until they are surpassed by the molecular contribution near the RCB, where the local pebble sizes become limited by growth with the parameters plotted here. For smaller planets in the outer disk, or those that are accreting pebbles at a higher rate, the pebble sizes at the RCB instead remain limited by fragmentation and erosion. We will look into the distance trend of this velocity-limited opacity regime in Sect. 4.5.

4.4 Trend with pebble accretion rate

Finally, we examine the opacity trend with the pebble accretion rate in Fig. 6, where it is varied between 10⁻⁷ and 10⁻⁴ M_⊕ yr⁻¹. This probes both the variation in solid mass flux and accretion luminosity, which are proportional (L = GM_cṀ_peb∕r_c) in our model. The effect on the pebble sizes (top panel) can be divided into two regimes, depending on whether the pebble size is limited by growth (orange color range) or velocity (yellow color range). In the former case, the pebble size scales positively with the accretion rate as an increased pebble mass flux increases the number of collisions during their sedimentation and allows for faster growth. In this regime, the pebble opacity (per Eq. (17)) is not explicitly dependent on the mass flux, aside from the extinction efficiency which is higher for larger pebbles in these conditions. The second effect is that a higher luminosity increases the local radiative temperature gradient, which alters the conditions of the local medium, ultimately affecting the opacity as well.

When the accretion rate is increased beyond ~10⁻⁵ M_⊕ yr⁻¹ for the plotted default parameters, the pebble size becomes velocity-limited at the fragmentation-erosion barrier. At this point, any additional material no longer increases the pebble size and directly increases their volume density, allowing the opacity to rise sharply. The increased opacity in these runs with high pebble accretion rates turns the envelopes almost entirely convective (see the upper two curves), with RCB locations at or close to the Bondi radius.

Fig. 7 Compilation grid of 106 runs, where the colors indicate the opacity at the RCB as a function of the planet’s semimajor axis and mass at four different pebble accretion rates. The white dashed lines divide the zones where different opacity contributions dominate. The red dotted lines mark the parameter space where the entire envelope is convective due to high opacities in the velocity-limited regime.

4.5 Description of three opacity regimes

In order to visualize the broad opacity trends more clearly, we also evaluate a large 2D grid of envelope opacities at the RCB (κ_rcb) as a function of their distance to the star and planetary mass. While the opacity throughout the entire radiative zone is relevant for the envelope’s structure, its value at the boundary between the radiative and convective zones is generally considered the most important. This is both because most of the radiative portion of the envelope’s mass is contained near the RCB and because it controls the rate at which the convective interior cools (i.e., Lee & Chiang 2015; Ginzburg et al. 2016).

The results of these runs are shown in Fig. 7, ordered into four panels that correspond to different pebble accretion rates, which typically increase as planets grow and capture pebbles more efficiently (Ormel & Klahr 2010; Lambrechts & Johansen 2012; Liu et al. 2019). The opacity trends reflect the discussions of the previous subsections but also contrast the opacity of solids and molecules, which generally follow a dichotomy based on the planet’s distance and mass. Hot molecules contribute most of the envelope’s opacity in the warm and dense inner disk, where solids sediment slowly and coalesce to form larger agglomerates. Because the molecular opacity increases with both temperature and density,this provides a clear contrast based on orbital separation. The molecular opacity at the RCB is not very sensitive to planet mass, as can be seen from the nearly vertical opacity contours. In contrast, the opacity from solids declines steeply as planets grow more massive and particles sediment more slowly due to the higher densities.

Figure 7 shows that the total envelope opacity can be divided into three regimes. The first, molecular opacity regime applies to higher-mass planets in the inner disk that are accreting pebbles at a low rate. In this regime, the opacity at the RCB is almost independent of planet mass and steeply declines with orbital separation. The second, growth-limited solids opacity regime applies to higher-mass planets in the outer disk that are accreting pebbles at a low rate. In this regime, the RCB opacity is almost independent from distance and declines steeply with planet mass. The opacity is invariant to an increase in the mass flux as any additional mass just adds to pebble growth. However, the RCB opacity scales positively with the accretion luminosity. The third, velocity-limited solids opacity regime applies to lower-mass planets that are accreting pebbles at a sufficiently high rate. It is typically characterized by fully convective envelopes, with the plotted RCB located at the Bondi radius. A higher pebble accretion rate increases the volume density of solids in this regime and greatly increases the parameter range where it applies (to include larger M_p and smaller d).

One caveat to this third regime is that if the entire outer envelope is convective, the envelope density remains relatively low and as pointed out by Johansen & Nordlund (2020), the envelope can then locally become radiative at the depth where solids sublimate (2500 K in our model). At this point, which roughly coincides with the dissociation temperature of molecular hydrogen, only the molecular contribution to the opacity remains. If there are sufficient free electrons to ionize the hydrogen and produce H⁻, this becomes thedominant opacity component (Lee & Chiang 2015). The situation in this regime is further complicated by the balancing effect of a compositional gradient that also begins to form at these temperatures (Bodenheimer et al. 2018; Müller et al. 2020), in which the Ledoux criterion (Ledoux 1947) has to be used to evaluate stability against convection, rather than the Schwarzschild criterion. The thermodynamic consequences of this transition region are an active topic of investigation (Valletta & Helled 2020; Ormel et al. 2021) and are outside the scope of this opacity study, where we focus specifically on the opacity of solids in the outer envelope. We provide a broader discussion of this potential zone in Sect. 6.

5 Implications for giant planet formation

Planets transition from slow to runaway gas accretion when they reach a critical point, where their self-gravity begins to exceed the envelope’s pressure support. A transparent envelope is characterized by more efficient radiative energy transport in the outer regions, which results in faster cooling and more rapid gas accretion. As a result, both analytical and numerical works (e.g., Mizuno 1980; Stevenson 1982; Pollack et al. 1996; Movshovitz et al. 2010; Lee & Chiang 2015) have shown that the mass at which planets transition to runaway growth is positively linked to the envelope’s opacity. In this section, we apply our model to examine this trend across the proto-planetary disk.

5.1 Trends in the critical metal mass

In a previous work (Brouwers & Ormel 2020), we showed that the traditional criterion of a critical core mass becomes meaningless if the growth of the core is limited. For this reason, we introduced and derived an analogous criterion called the “critical metal mass”, which measures the total metal content of a planet at the onset of runaway gas accretion. It shares the opacity dependence with previous expressions for the critical core mass but is additionally an explicit function of the core mass itself. One limitation of our previous work was that we modeled envelopes with an isothermal radiative layer, which is a good assumption for grain-free envelopes but invalid for envelopes with higher opacities due to the presence of solids. In order to correctly account for thermodynamic changes in the radiative part of the envelope, we slightly modify the analytical structure model of Brouwers & Ormel (2020) in Appendix A to be applicable to envelopes with non-isothermal radiative regions. The critical metal mass, defined as the total mass of solids accreted at the onset of runaway accretion, then becomes an explicit function of both T_rcb and κ_rcb: $$\begin{array}{ll}{M}_{\text{z},\text{crit}}\hfill & \approx 5.5M_{\oplus }{\left(\frac{{\kappa }_{\text{rcb}}}{0.01\text{g}{\text{cm}}^{-2}}\right)}^{\frac{1}{6}}{\left(\frac{d}{\text{AU}}\right)}^{\frac{7}{108}}{\left(\frac{T_{\text{vap}}}{2500\text{K}}\right)}^{\frac{8}{27}}\hfill \\ \hfill & {\left(\frac{{\dot{M}}_{\text{peb}}}{{10}^{-6}{M}_{\oplus }{\text{yr}}^{-1}}\right)}^{\frac{1}{6}}{\left(\frac{M_{\text{c}}}{M_{\oplus }}\right)}^{\frac{1}{2}}{\left(\frac{T_{\text{rcb}}}{T_{\text{disk}}}\right)}^{-\frac{126}{972}}.\hfill \end{array}$$(27)

In line with the findings by Ormel et al. (2021), we assume core growth to be uniformly limited to 2 M_⊕, with the rest of the solids being absorbed in the deep envelope. We plot the resultant critical metal masses as a function of orbital separation for three accretion rates in Fig. 8. The black lines include the opacity from dust and pebbles, contrasted with the gray lines that assume entirely grain-free envelopes.

The general shape of the curves in Fig. 8 reflects the opacity trends described in the previous section. Molecules dominate the opacity of envelopes in the hot inner disk, where the critical mass is seen to decrease with orbital separation. In the grain-free curves, this downward trend continues toward the outer disk, where the critical mass reaches very low values of only a few Earth masses. Accounting for the opacity of pebbles and dust breaks this trend in the intermediate disk, where solids begin to dominate the opacity of the envelopes. This halts the decline of the critical mass, which then instead flattens toward the outer disk. Besides its variation with distance, the critical metal mass is also positively dependent on the solids accretion rate. This is both due to the increased luminosity and the increased opacity from solids, although the two are linked thermodynamically (see Sect. 4.4). We note that in reality, a pebble-accreting planet would not experience a single pebble accretion rate during its evolution but one that generally increases as it gains in mass and is able to capture pebbles more efficiently.

Physically, the trends in the critical metal mass translate to similar trends in the occurrence rates of gas giant planets. Data from the Kepler mission have provided good exoplanet abundance statistics up to about 1 AU in semimajor axis (Borucki et al. 2010, 2011; Batalha et al. 2013), which show a general increase in the abundance of giant planets with distance in this range (Howard et al. 2012; Dong & Zhu 2013; Santerne et al. 2016). When combined with data from radial-velocity measurements (Mayor et al. 2011), it shows a break in the giant planet occurrence rate between 2 and 3 AU, with a declining power law beyond this value (Fernandes et al. 2019). While other factors are certainly at play in planet formation across the proto-planetary disk, these trends are generally well matched by the critical metal mass predicted by our opacity model.

Fig. 8 Trends in the critical metal mass as a function of distance. The different lines correspond to a variation in solids accretion rates: 10−7 M⊕ yr−1 (dotted), 10−6 M⊕ yr−1 (dashed), 10−5 M⊕ yr−1 (solid), with proportional gas accretion Ṁxy∕Ṁpeb = 5 and the remaining parameters set by Table 1. Planets with grain-free envelopes (gray lines) exhibit a downward trend with orbital separation because the molecular opacity scales positively with temperature. The opacity generated by pebbles and dust (black lines) increases the critical mass significantly in the outer disk, especially at higher accretion rates.

5.2 Implications for the formation of Uranus and Neptune

Far out in the disk, one of the key theoretical challenges is to explain how Uranus and Neptune accreted their substantial masses of ~ 15 M_⊕ without entering runaway growth (Helled & Bodenheimer 2014; Helled et al. 2020). In the model developed by Lambrechts et al. (2014), it is suggested that the accretion of solids continued throughout the disk lifetime and that the heat from impactors kept the envelope stable during its evolution. However, the planet’s cold birth environments and low molecular opacity means that their accretion luminosity is easily radiated away, causing grain-free envelopes to implode into gas giants before they have grown beyond a few Earth masses. In order to prevent runaway accretion, the planets must either form when most of the nebula is gone (Lee & Chiang 2016; Ogihara et al. 2020), or their envelopes must still be sufficiently opaque that their heat does not escape so quickly. The common approach to this problem is to assume that the envelopes contain a large ISM-like opacity from small dust grains. In line with previous works on coagulation (Movshovitz & Podolak 2008; Movshovitz et al. 2010; Mordasini 2014; Ormel 2014), our results indicate that such a high opacity is not realistic when grain growth determines the pebble sizes. It is possible, however, for planetary envelopes in the outer disk to attain these high opacities if the pebble sizes are velocity-limited by fragmentation or erosion. As we show in Fig. 7, this requires the planet to continuously accrete pebbles at a sufficiently high rate.

Alternatively, planets begin to carve partial gaps around their orbit when they grow beyond several Earth masses (Paardekooper & Mellema 2004, 2006). When they continue to grow, they locally invert the pressure gradient in the disk and stop accreting pebbles, a criterion that is known as the pebble isolation mass (Morbidelli et al. 2015; Morbidelli & Nesvorny 2012; Ataiee et al. 2018). As shown by Bitsch et al. (2015, 2018), it is important to account for the cooling of the Solar nebula and the reduced disk scale height over time. In these simulations, the pebble disk scale height is found to reduce sufficiently in the first few megayears that both Uranus and Neptune can reach their isolation mass and stop accreting pebbles. This introduces a new issue however, as ceasing the supply of accretion heat can again easily trigger cooling and rapid gas accretion. Unlike the close-in super-Earths and mini-Neptunes, which can be prevented from accreting gas by entropy advection (Ormel et al. 2015; Cimerman et al. 2017; Kurokawa & Tanigawa 2018; Ali-Dib & Thompson 2020; Moldenhauer et al. 2020), the ice giants will rapidly accrete gas due to their Kelvin-Helmholtz contraction. The solution offered by Alibert et al. (2018) and Guilera et al. (2020) is that the reduction in pebble accretion can be compensated by the accretion of planetesimals, which could naturally form at the surrounding pressure bump. In a similar vein, Lambrechts & Lega (2017), Chen et al. (2020), and Bitsch & Savvidou (2021) reason that while centimeter-sized pebbles are easily blocked by a local pressure bump, small dust (≤ 200 μm) is sufficiently coupled to the gas to pass through (Pinilla et al. 2016; Bitsch et al. 2018; Weber et al. 2018; Haugbølle et al. 2019) and provide the required opacity to prevent rapid cooling. Our results indicate that this solution essentially faces the same problem as continuous pebble accretion: regardless of whether the mass is supplied by pebbles or planetesimals, the opacity from solids is insufficient to prevent runaway growth without a high accretion rate. In fact, the required planetesimal accretion rate is even higher than the equivalence in pebbles because the abundance of small solids is insufficient to reach the velocity-limited regime and the envelope opacity thus remains low.

The implication is that regardless of whether the ice giants became isolated from the pebble flow, their formation requires the constant accretion of material at high rates if they formed in-situ. This introduces a timescale problem as they must then have formed in the window where they had enough time to accrete their observed masses, but not too early such that they accreted too much material and became gas giants. Accurately estimating the duration of this formation window requires more detailed models of the ice giant’s interiors, which are currently not yet well constrained (Helled et al. 2020; Vazan et al. 2020). In addition, it is important to account for the cooling of the solar nebula over time. The cooling of the surrounding disk has the same general effect as increasing the planet’s distance from the central star (reduced T_disk, ρ_disk), which we showed increases the envelope opacity contribution of solids. In a future work, we plan to quantitatively estimate the required accretion rate and formation window of these ice giants.

6 Discussion

6.1 Comparison with contemporary works

In this work, we focused on the growth and destruction of solids as they travel through planetary envelopes, ignoring their evolution prior to accretion. A contrasting approach was taken in a contemporary work by Bitsch & Savvidou (2021), who modeled thesize distribution of solids in the disk similar to Savvidou et al. (2020), with the assumption that their distribution and opacity in the disk mid-plane also apply to the envelope’s interior. We present evidence that the opacities at the RCB are actually very different from the disk due to a typically much higher gas temperature and density, as well as the significant size evolution of solids in the envelope. Larger pebbles fragment or face erosion when they enter planetary envelopes, while the lower end of the accreted size distribution experiences significant growth during sedimentation. Effectively, we predict much smaller envelope opacities in the growth-limited regime than Bitsch & Savvidou (2021), although high opacities remain possible in the velocity-limited regime. Nevertheless, while we argue that the sizes of solids in the envelope are disconnected from their size distribution in the disk, modeling their size evolution in the disk remains important for physical estimates of their accretion rates and efficiency.

An alternative approach to modeling the envelope opacity is presented by Chachan et al. (2021), which is based on the work by Lee & Chiang (2015). The authors argue that the sublimation of dust leads to a second radiative zone at the sublimation front, with an inner boundary around 2500 K, where dissociation of H₂ provides a new opacity source in the form of H⁻. The appearance of this inner radiative zone is based on the assumption that the envelope’s outer layers are convective due to a high dust opacity. The requirement of free elections to form H⁻ means that the opacity studied by Chachan et al. (2021) is proportional to the local metalicity of the gas, which the authors take equal to the local dust-to-gas ratio of the disk, yielding an opacity minimum at intermediate distances (1–10 AU). While this is an interesting result, we note that the interior region around 2500 K is highly complex due to the compositional gradient that forms at this temperature, which affects the criterion for convective stability (Ledoux 1947; Müller et al. 2020). The precise value of the metallicity at 2500 K is, therefore, hard to constrain without a detailed interior model that accounts for the saturated vapor pressure of the sublimated solids, which determines the envelope absorption of metals and is highly variable with temperature and composition (Ormel et al. 2021). Besides the compositional complexity, the appearance of a second radiative region furthermore requires a high dust opacity in the outer envelope, which we show is significantly reduced when coalescence is accounted for. Nevertheless, the work by Chachan et al. (2021) hints at a more complex scenario for the fully convective envelopes that we encounter in the velocity-limited opacity regime. The details of a potential deeper RCB are likely crucial in any model where it appears. It is, therefore, warranted to engage in a detailed study that also accounts for the important effects in this temperature range, including the varying composition of the envelope.

6.2 Pebble versus planetesimal accretion

Our work also highlights the difference between the formation of planets with pebbles or planetesimals. Whereas pebbles naturally provide an abundance of small grains in the outer envelope, planetesimals are much larger (>1 km, typically several hundred km, (i.e., Li et al. 2019; Klahr & Schreiber 2020; Rucska & Wadsley 2021) and only produce smaller grains when they experience significant mass loss or a dynamical fracture, where the pressure on the front of the planetesimal exceeds its strength and the planetesimal bursts apart. For this to occur, the planetesimal must be moving at very high velocity through a region with sufficient density, a condition that is only fulfilled neara planet’s central core and far below the RCB (Podolak et al. 1988; Mordasini et al. 2015; Pinhas et al. 2016; Brouwers et al. 2018; Valletta & Helled 2019). Similarly, mass ablation through frictional heating is inversely effective with planetesimal mass and only becomes important for planetesimal-sized objects close to the planet’s central core. Besides this, the localized nature of planetesimal impacts alsomakes it more difficult to effectively obscure an envelope, allowing unaffected parts to continue radiating away energy.

6.3 Model caveats and improvements

The main simplification of our opacity model is the assumption that the distribution of the solids in a planetary envelope can be modeled with a characteristic local pebble size, supplemented with a constant-size dust population. This assumption has been shown to provide good agreement with coagulation simulations (Mordasini 2014) and we argue that it is also justified when the pebbles are limited by a velocity threshold and pile up at a single size. Modeling a more complex size distribution self-consistently during sedimentation would require a computationally intensive approach with a large number of size bins and experimental data on the redistribution function of grains produced by both erosion and collisions. Furthermore, grains of different compositions are characterized by different strengths and sublimation temperatures, meaning that the maximum grain size will in reality fluctuate depending on the temperature regime where different species can survive. Finally, we considered quasi-static envelopes and did not model any convective overshoot from the envelope’s deeper regions, which can potentially return a portion of the particles to the surface (Popovas et al. 2018, 2019) after they experience collisional grind-down from the additional acceleration (Ali-Dib & Thompson 2020; Johansen & Nordlund 2020).

7 Summary and conclusions

In the scenario of pebble accretion, planets grow by the subsequent accretion of solids (micron-sized dust + larger pebbles) and gas. The amount of gas that a planet is able to bind at a certain mass depends on its cooling rate. This, in turn, is set by the envelope’s opacity. In this work, we designed an opacity model that incorporates the main physical processes that influence the evolution of solids in planetary envelopes. Our approach can be considered an extension to previous works by Mordasini (2014) and Ormel (2014) as we model two populations of solids instead of one and add the effects of pebble erosion, fragmentation, and dust sweep-up by pebbles to a grain growth prescription. We formulate convenient analytical expressions for the pebble opacity in different regimes (Eqs. (17)–(20)) and we estimate the dust opacity in a steady-state between dust sweep-up and production (Eq. (23b)). Finally, we apply our model across a wide parameter space to map out the resulting opacity trends as a function of distance, planet mass and accretion rate. Our main findings are that:

1. The initial size distribution of solids that a planet accretes is not an accurate representation of their sizes within the envelope. At the upper end of the distribution, we corroborate the recent findings by Ali-Dib & Thompson (2020) and Johansen & Nordlund (2020) that larger pebbles are prone to erosion by dust grains when they accrete onto planetary envelopes. This erosion occurs as a result of the additional gravitational acceleration by the planet and does not require violent convective gas motion. The size limit at entry varies from a few microns in envelopes of small planets in the outer disk to several centimeters for more massive planets that are located closer to the central star.

2. Within the envelope, the pebble sizes are limited by their rate of growth or by the velocity at which they begin to erode or fragment. The growth-limited regime applies to larger planets that are accreting pebbles at a low rate. As found by Mordasini (2014) and Ormel (2014), the opacity in this regime drops far below ISM values, especially at greater depth or in the envelopes of more massive planets. Any additional pebbles that a planet accretes in the growth-limited opacity regime only lead to larger pebbles and do not increase the opacity, other than through any added luminosity. Nevertheless, the contribution of solids typically still dominates over the molecules, with the exception of higher-mass planets in the inner disk (see Fig. 7). Smaller planets at greater distance from the star or those that are accreting pebbles at a sufficiently high rate enter the velocity-limited regime, where pebble growth is restricted by erosion or fragmentation. In this regime, a higher pebble accretion rate further increases the pebble abundance and the opacity is able to remain high throughout the envelope, typically turning it entirely convective.

3. One proposed solution for the fact that Uranus and Neptune accreted their substantial masses of ~15 M_⊕ without entering runaway gas accretion is that they had a high envelope opacity during formation. We find that the commonly made assumption of an ISM-like opacity from small dust grains is not realistic when grain growth is accounted for. High envelope opacities remain possible if the pebble accretion rate is sufficiently high and the pebble sizes become limited by fragmentation or erosion, rather than by their rate of growth. However, this introduces a timescale problem as the planets must have formed at the right time to grow to their current sizes at a sufficiently rapid pace without becoming gas giants.

Acknowledgements

Part of the work presented here is based on discussions conducted at the second ISSI ice giants Meetings in Bern, March 2020. M.G.B. acknowledges the support of a Royal Society Studentship, RG 160509. A.B. is grateful to the Royal Society for a Dorothy Hodgkin Fellowship. We thank an anonymous referee for comments that helped to substantially improve the clarity of this work.

Appendix A Expression of the critical metal mass with a non-isothermal radiative layer

In this appendix, we slightly modify the analytical structure model of Brouwers & Ormel (2020) to derive an expression for the critical metal mass that is more appropriate to study the effects of envelope opacity. The original model consists of an outer radiative region, followed by an intermediate convective region with solar composition and a polluted interior convective region that consists primarily of silicate vapor. Because that model was built with a grain-free envelope in mind, the radiative region was assumed to be isothermal. This ceases to be a good approximation when the opacity from solids is included and the radiative temperaturegradient becomes important.

Fortunately, it is straightforward to modify the model to be applicable to a scenario with higher envelope opacities. We shift the temperature boundary of the intermediate convective region from T_disk to T_rcb and recalculate the thermodynamic structure of the intermediate layer. This amounts to a single variable change in Eqs. (9a)–(9c) of Brouwers & Ormel (2020): $${r^{\prime }}_{\text{B}}=\frac{{\gamma }_{\text{xy}}-1}{{\gamma }_{\text{xy}}}{r}_{\text{B}}\text{to}{r^{\prime }}_{\text{B}}=\frac{{\gamma }_{\text{xy}}-1}{{\gamma }_{\text{xy}}}\frac{T_{\text{disk}}}{T_{\text{rcb}}}{r}_{\text{B}}$$(A.1)

which in effect lowers the density of the intermediate region and increases its temperature. The relation between the opacity and density at the RCB then becomes: $${\kappa }_{\text{rcb}}=\frac{64\pi \overline{\sigma }{T}_{\text{rcb}}^4{r^{\prime }}_{\text{B}}}{3{\rho }_{\text{rcb}}L}.$$(A.2)

The corresponding structure equations of the polluted region remain unchanged. The only modification we have to make is at the outer boundary of the polluted region, where the density changes from $${\rho }_{\text{vap}}={\left(\frac{T_{\text{vap}}}{T_{\text{disk}}}\right)}^{\frac{1}{{\gamma }_{\text{xy}}-1}}\frac{{\mu }_{\text{g}}}{{\mu }_{\text{xy}}}{\rho }_{\text{rcb}}\text{to}{\rho }_{\text{vap}}={\left(\frac{T_{\text{vap}}}{T_{\text{rcb}}}\right)}^{\frac{1}{{\gamma }_{\text{xy}}-1}}\frac{{\mu }_{\text{g}}}{{\mu }_{\text{xy}}}{\rho }_{\text{rcb}}.$$(A.3)

In our model for a polluted envelope, most of the mass is contained in the inner regions. In order to derive the critical metal mass, we substitute the modified variables into Eq. (22b) of Brouwers & Ormel (2020) and follow the same steps as in Sect. 4.2 of that work. In practice, this comes down to substituting T_rcb for T_d in their Eq. (27). The critical metal mass is approximated to occur at the crossover mass where M_xy = M_z, which yields: $$\begin{array}{ll}{M}_{\text{z},\text{crit}}\hfill & \approx 5.5M_{\oplus }{\left(\frac{{\kappa }_{\text{rcb}}}{0.01\text{g}{\text{cm}}^{-2}}\right)}^{\frac{1}{6}}{\left(\frac{d}{\text{AU}}\right)}^{\frac{7}{108}}{\left(\frac{T_{\text{vap}}}{2500\text{K}}\right)}^{\frac{8}{27}}\hfill \\ \hfill & {\left(\frac{{\dot{M}}_{\text{peb}}}{{10}^{-6}{M}_{\oplus }{\text{yr}}^{-1}}\right)}^{\frac{1}{6}}{\left(\frac{M_{\text{c}}}{M_{\oplus }}\right)}^{\frac{1}{2}}{\left(\frac{T_{\text{rcb}}}{T_{\text{disk}}}\right)}^{-\frac{126}{972}}.\hfill \end{array}$$(A.4)

The only difference in the final expression is the appearance of an inverse scaling with T_rcb∕T_disk. This is due to a combination of two opposing effects. The first effect is that a higher RCB temperature leads to a lower RCB density at the same pressure and, therefore, results in an increased critical mass. The second trend is that a higher RCB temperature moves the location of the RCB inward to a higher pressure and an increased density at the same opacity, leading to a reduced critical mass. In combination, these two trends yield a slight inverse scaling in the final expression. We show a comparison with the expression from Brouwers & Ormel (2020) in Fig. A.1, where the non-isothermal expression is shown in black as opposed to gray for the original expression. As expected, the difference is greatest when there is a large opacity from solids. Grain-free envelopes remain well approximated by the isothermal model. This is because their opacity scales positively with temperature and density, such that their radiative temperature gradients are steep functions of depth. We also plot the critical metal mass for a run where we keep the pebble size constant at 100 μm. In agreement with the numerical results from Ormel et al. (2021), we find that this yields a flat line. The quantitative values from Eq. (A.4) are generally a factor ~2 reduced compared to detailed numerical simulations with non-ideal equations of state performed by Ormel et al. (2021) but the simple derived scaling captures the main physical trends well.

Fig. A.1 Comparison of the critical metal mass expressions of Brouwers & Ormel (2020) (isothermal, gray) and this work (non-isothermal, black), both examined with different opacities. The dotted lines indicate grain-free envelopes, whereas the solid lines assume the accretion of small pebbles that can grow, fragment and erode according to our model. The dash-dotted line corresponds to the accretion of 100 μm pebbles thatare kept at a constant size during their sedimentation.

Appendix B Justification for neglecting dust growth in production-sweep-up steady state

In this Appendix, we investigate whether dust growth can be neglected in our steady-state model of dust production and sweep-up. To begin, we estimate the typical distance that dust grains travel before mutually colliding (l_col,d), which is analogous to Eq. (11) for pebbles: $$l_{\text{col},\text{d}}=\frac{m_{\text{d}}}{{\rho }_{\text{d}}{\sigma }_{\text{d}}}\frac{v_{\text{sed},\text{d}}}{v_{\text{col},\text{d}}}.$$(B.1)

We compare this quantity with the typical distance that these dust grains can travel before they are swept up by a sedimenting pebble (l_sweep,d), which in turn is analogous to Eq. (8) for pebbles: $$l_{\text{sweep},\text{d}}=\frac{m_{\text{peb}}}{{\rho }_{\text{peb}}{\sigma }_{\text{peb}}}\frac{v_{\text{sed},\text{d}}}{v_{\text{fall},\text{peb}}},$$(B.2)

where we assumed that v_fall,peb ≫ v_fall,d. As discussed in Sect. 3.2, the dust and pebble volume densities in a steady state between dust sweep-up sweep-up and production are related through Eq. (22). This gives the ratio of the relevant length scales as: $$\frac{l_{\text{col},\text{d}}}{l_{\text{sweep},\text{d}}}=F^{-1}{x}_{\text{R}}^{-2}\frac{R_{\text{d}}{v}_{\text{fall},\text{peb}}}{R_{\text{peb}}{v}_{\text{fall},\text{d}}}$$(B.3) $$\le \frac{9}{F},$$(B.4)

where the second line substitutes the scaling of Epstein drag, the drag regime with the most dust collisions. It is evident that even in this regime, l_col,d > l_sweep,d and the dust grains are typically swept up by pebbles before they collide with other dust grains. This remains valid with the highest possible value of the dust production constant F = 1, so we conclude that significant dust growth is not expected in a steady statebetween dust production and sweep-up.

Appendix C Variation of the limiting velocity

The most important free parameter in our model is the limiting velocity (Eq. (13)), which represents the upper limit to the terminal velocity of pebbles during their sedimentation. It combines the velocity limits from erosion (Eq. (9a)) and fragmentation. In the main text, we assumed a limiting velocity of 2.4 m s⁻¹, corresponding to the onset of erosion by micron-sized dust grains (Schräpler et al. 2018). Because smaller dust grains are more effective at removing mass, the erosion velocity scales inversely with grain size. In Fig. C.1, we vary the limiting velocity between 0.4–40 m s⁻¹, which corresponds to erosion from grain sizes between 0.05 and 100 μm.

Figure C.1 illustrates that the value of the critical velocity has important consequences for the opacity. If the limiting velocity is reduced below our default value of 2.4 m s⁻¹, the pebbles are velocity-limited to smaller sizes and the opacity can significantly increase. Due to this, a broader range of envelopes become fully convective when they experience a given pebble accretion rate. The opposite trend is also true when the limiting velocity is increased. A higher velocity limit means that pebbles are more often limited by growth, which in turn means that higher accretion rates are required to keep envelopes of the same mass convective. Within an envelope whose sedimenting pebbles are already limited by growth, a further increase of the limiting velocity has no effect.

Fig. C.1 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). This figure is the same as Fig. 4, but now the mass is fixed at the default 5 M⊕ and the limiting velocity (Eq. (13)) is varied a logarithmically between 0.4− 40 m s−1.

References

All Tables

All Figures

Fig. 1 Schematic diagram of our two-population opacity model. Depending on the velocities of sedimenting pebbles, they can experience growth (coalescence), fragmentation, or erosion. In the growth-limited regime, the pebble size in the interior is regulated by the collision and sedimentation timescales (Rpeb = Rcoal). In the velocity-limited regime, fragmentation and/or erosion restrict the pebble size below the growth potential (Rpeb = Rvlim < Rcoal). In our model, dust grains of constant size are produced in collisions between pebbles or in high-velocity pebble-dust encounters (erosion), while they are lost by sticking encounters with slow-moving pebbles (dust sweep-up), leading to a local steady state in their abundance.

In the text

Fig. 2 Erosive mass loss of pebbles that enter the envelope of a growing planet at its Bondi radius (rB). The parameters of the planet are shown in Table 1. The three lines correspond to initial pebble entry sizes of 0.1 cm (dotted), 1 cm (dashed) and 10 cm (solid). Due to the positive feedback of erosive mass-loss on the dust abundance, the pebbles rapidly converge to the size Rerosion (Eq. (7)) where collisions with dust grains switch from causing mass loss via erosion to sticking (sweep-up).

In the text

Fig. 3 Characteristic pebble size limits at the Bondi radius, plotted with the parameters of Table 1. The colors indicate the pebble size where vfall,peb = vlim = 2.4 m s−1, the common limit set via erosion with micron-sized grains and fragmentation. The white dashed line separates the Stokes and Epstein regimes while the dotted line indicates the regions where either growth rate or erosion/fragmentation limits the local pebble size. The plus sign marks the default planet mass (5 M⊕) and distance (5 AU).

In the text

Fig. 4 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). The different lines indicate a range of planet masses, with the arrow indicating a logarithmic progression from 0.5− 20 M⊕. The triangles indicate the location of the RCB, while the stars indicate the depth where the ambient temperature exceeds the sublimation temperature (2500 K) and the opacity from solids vanishes. The colors in the top panel show the ratio of the pebble size in the model relative to Rvlim (Eq. (14)). The colors in the two lower panels show the relative value of the indicated opacity to the total opacity.

In the text

	Fig. 5 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). This figure is the same as Fig. 4 but now the mass is fixed at the default 5 M⊕ and the planet’s distances from the star are varied logarithmically from 0.1−30 AU.
In the text

	Fig. 6 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). This figure is the same as Fig. 4 but now the mass is fixed at the default 5 M⊕ and the pebble accretion rates are varied logarithmically from 10−7−10−4 M⊕ yr−1.
In the text

Fig. 7 Compilation grid of 106 runs, where the colors indicate the opacity at the RCB as a function of the planet’s semimajor axis and mass at four different pebble accretion rates. The white dashed lines divide the zones where different opacity contributions dominate. The red dotted lines mark the parameter space where the entire envelope is convective due to high opacities in the velocity-limited regime.

In the text

Fig. 8 Trends in the critical metal mass as a function of distance. The different lines correspond to a variation in solids accretion rates: 10−7 M⊕ yr−1 (dotted), 10−6 M⊕ yr−1 (dashed), 10−5 M⊕ yr−1 (solid), with proportional gas accretion Ṁxy∕Ṁpeb = 5 and the remaining parameters set by Table 1. Planets with grain-free envelopes (gray lines) exhibit a downward trend with orbital separation because the molecular opacity scales positively with temperature. The opacity generated by pebbles and dust (black lines) increases the critical mass significantly in the outer disk, especially at higher accretion rates.

In the text

Fig. A.1 Comparison of the critical metal mass expressions of Brouwers & Ormel (2020) (isothermal, gray) and this work (non-isothermal, black), both examined with different opacities. The dotted lines indicate grain-free envelopes, whereas the solid lines assume the accretion of small pebbles that can grow, fragment and erode according to our model. The dash-dotted line corresponds to the accretion of 100 μm pebbles thatare kept at a constant size during their sedimentation.

In the text

	Fig. C.1 Pebble growth tracks (a), their resulting pebble opacity (b), and produced dust opacity (c) for a standard set of model runs at 5 AU (see Table 1). This figure is the same as Fig. 4, but now the mass is fixed at the default 5 M⊕ and the limiting velocity (Eq. (13)) is varied a logarithmically between 0.4− 40 m s−1.
In the text