# Pebble-isolation mass: Scaling law and implications for the formation of super-Earths and gas giants

Bertram Bitsch<fnref n="1" />, Alessandro Morbidelli<fnref n="2" />, Anders Johansen<fnref n="1" />, Elena Lega<fnref n="2" />, Michiel Lambrechts<fnref n="1" />, and Aurélien Crida<fnref n="2" />,<fnref n="3" />

Received 11 September 2017 / Accepted 4 January 2018

## ABSTRACT

The growth of a planetary core by pebble accretion stops at the so-called pebble isolation mass, when the core generates a pressure bump that traps drifting pebbles outside its orbit. The value of the pebble isolation mass is crucial in determining the final planet mass. If the isolation mass is very low, gas accretion is protracted and the planet remains at a few Earth masses with a mainly solid composition. For higher values of the pebble isolation mass, the planet might be able to accrete gas from the protoplanetary disc and grow into a gas giant. Previous works have determined a scaling of the pebble isolation mass with cube of the disc aspect ratio. Here, we expand on previous measurements and explore the dependency of the pebble isolation mass on all relevant parameters of the protoplanetary disc. We use 3D hydrodynamical simulations to measure the pebble isolation mass and derive a simple scaling law that captures the dependence on the local disc structure and the turbulent viscosity parameter $\alpha$. We find that small pebbles, coupled to the gas, with Stokes number $\tau_f < 0.005$ can drift through the partial gap at pebble isolation mass. However, as the planetary mass increases, particles must be decreasingly smaller to penetrate the pressure bump. Turbulent diffusion of particles, however, can lead to an increase of the pebble isolation mass by a factor of two, depending on the strength of the background viscosity and on the pebble size. We finally explore the implications of the new scaling law of the pebble isolation mass on the formation of planetary systems by numerically integrating the growth and migration pathways of planets in evolving protoplanetary discs. Compared to models neglecting the dependence of the pebble isolation mass on the $\alpha$-viscosity, our models including this effect result in higher core masses for giant planets. These higher core masses are more similar to the core masses of the giant planets in the solar system.

**Key words.** accretion, accretion discs – planets and satellites: formation – protoplanetary discs – planet-disc interactions

## 1. Introduction

Protoplanetary discs consist of gas and approximately 1% dust grains. These grains collide and grow to millimeter and even centimeter sizes (Brauer et al. 2008; Güttler et al. 2010). These particles are often referred to as pebbles. Pebbles interact with the gas disc through gas drag and drift inwards (Weidenschilling 1977; Brauer et al. 2008). Pebbles can become concentrated in pressure bumps and through the streaming instability (see Johansen et al. 2014 for a review), leading to planetesimal formation by gravitational collapse of the filaments. Planetesimals formed by the streaming instability have characteristic sizes of 100 km (Johansen et al. 2015; Simon et al. 2015).

In classical planet formation models, the cores of the giant planets form through mutual collisions between these planetesimals (Pollack et al. 1996). However, to achieve a core mass high enough (≈10 $M_E$) to attract a gaseous envelope, a surface density of planetesimals of a few times the Minimum Mass Solar Nebula (MMSN) is needed. Additionally, the growth timescale increases steeply with orbital distance, making the formation of the ice giants in the solar system basically impossible to achieve with planetesimal accretion alone. Gravitational stirring of the planetesimals by a set of growing protoplanets decreases the growth rates even more (Levison et al. 2010).

In recent years, a new paradigm of solid accretion has emerged: pebble accretion (Johansen & Lacerda 2010; Ormel & Klahr 2010; Lambrechts & Johansen 2012; Morbidelli & Nesvorny 2012). When a pebble enters the planetary Hill sphere, it is subject to gas drag, which robs the pebble of angular momentum, resulting in an inward drift of the pebble onto the planet. When the largest planetesimals have grown to a few hundred kilometers in size by accreting other planetesimals, rapid pebble accretion allows further growth to cores of ten Earth masses well within the lifetime of the protoplanetary disc (Ida et al. 2016; Visser & Ormel 2016; Johansen & Lambrechts 2017).

A growing planet opens a partial gap in the protoplanetary gas disc, which influences the motion of solids in the disc (Paardekooper & Mellema 2006; Rice et al. 2006). Pebble accretion stops when the gap carved by the planet generates a pressure maximum outside of its orbit, which stops the inward flux of pebbles (Morbidelli & Nesvorny 2012; Lambrechts et al. 2014). This is referred to as the pebble-isolation mass. As the influx of pebbles is stopped, the planet's gas envelope loses its hydrostatic support, and the envelope can then contract to form a planet with an extensive gaseous atmosphere (Lambrechts et al. 2014). In the solar system, pebble isolation is a potential mechanism to explain the dichotomy between the ice and gas giants, with ice giants never reaching the pebble-isolation mass and hence not able to

<footnote n="1">¹ Lund Observatory, Department of Astronomy and Theoretical Physics, Lund University, 22100 Lund, Sweden
e-mail: bitsch@mpia-hd.mpg.de</footnote>

<footnote n="2">² University Nice-Sophia Antipolis, CNRS, Observatoire de la Côte d'Azur, Laboratoire LAGRANGE, CS 34229, 06304 Nice Cedex 4, France</footnote>

<footnote n="3">³ Institut Universitaire de France, 103 Boulevard Saint-Michel, 75005 Paris, France</footnote>
undergo gas accretion within the lifetime of the protoplanetary disc (Lambrechts et al. 2014; Venturini & Helled 2017; Frehlik & Murray-Clay 2017).

Lambrechts et al. (2014) used hydrodynamical simulations to infer a pebble-isolation mass given by

$$M_{\text{iso}}^{\text{L14}} \approx 20 \left( \frac{H/r}{0.05} \right)^3 M_{\text{E}}, \qquad (1)$$

where $H/r$ is the disc aspect ratio. Lambrechts et al. (2014) found a weak dependence on the viscosity parameter. This is not surprising, since a dependence of gap opening on viscosity was also reported by Crida et al. (2006). However, because the dependence on the viscosity was not explored in detail in Lambrechts et al. (2014), an explicit mapping of the pebble-isolation mass as a function of viscosity is of crucial importance. Additionally, the headwind felt by the particles depends on both the disc aspect ratio, as evident in Eq. (1) above, and on the radial, initially unperturbed, pressure gradient $\partial \ln P / \partial \ln r$.

Here, we investigate the dependence of the pebble-isolation mass on all local disc parameters, namely the disc aspect ratio $H/r$, the viscosity $\nu$, and the pressure gradient of the disc $\partial \ln P / \partial \ln r$. In order to probe this large parameter space, we adopt 3D isothermal simulations executed with the FARGOCA code (Lega et al. 2014; Bitsch et al. 2014). As particles with different sizes are coupled in different ways to the gas disc, we additionally integrate the trajectories of single pebbles with various sizes (and therefore Stokes numbers) in the gas disc to probe which particle sizes can be trapped in the pressure bump as a function of planet mass.

The paper is organised as follows. In Sect. 2 we present our hydrodynamical set-up and discuss the different parameters that influence the pebble-isolation mass, to which we provide a fit in absence of turbulent diffusion. In Sect. 3 we integrate the trajectories of single pebbles in discs with embedded planets and infer the pebble sizes that are trapped in the pressure bumps. We then discuss turbulent diffusion of dust particles through the pressure bump. We also present the fitting formula for the pebble-isolation mass including turbulent diffusion, which is useful for planet formation simulations involving pebble accretion, in Sect. 3.3. In Sect. 4 we show the influence of the new-found pebble-isolation mass on simulations of planet formation, where we compare our results to Bitsch et al. (2015b). We additionally discuss implications of our results in Sect. 5 and finally summarise in Sect. 6.

# 2. Hydrodynamic simulations

## 2.1. Simulation set-up

In order to simulate the 3D disc-planet interaction, we used the 3D hydrodynamical code FARGOCA (Lega et al. 2014; Bitsch et al. 2014) in a locally isothermal configuration, where the radial temperature profile remains fixed throughout the simulation. We used the locally isothermal configuration because it allows a fast probing of parameter space in $\alpha$, $H/r$, $\partial \ln P / \partial \ln r$, and planetary masses, which is needed to constrain the pebble-isolation mass. Here $\alpha$ is related to the viscosity through $\nu = \alpha H^2 \Omega_K$ (Shakura & Sunyaev 1973), where $\Omega_K$ denotes the Keplerian rotation. Additionally, locally isothermal simulations allow an easier probing of the aspect ratio because it can be set as an input parameter in contrast to simulations with heating and cooling, where the aspect ratio is set ultimately by the opacity profile that determines the cooling. Nevertheless, in Sect. 2.7 we test

the predictions made with the isothermal simulations against simulations with heating and cooling.

We used for our simulations a 3D grid in spherical coordinates $(r, \phi, \theta)$ with 315, 720, and 32 grid cells. Our grid ranged from 0.4 to 2.5 in radius, where the planet is located at 1, and spanned the full azimuthal range. We used evanescent boundary conditions for the radial boundaries to damp out the spiral waves exerted by the planet in order to avoid disturbances due to reflections of the spiral waves caused by the planet. We simulated different values of the viscosity parameter $\alpha$, the aspect ratio $H/r$, the pressure gradient $\partial \ln P / \partial \ln r$, and a range of planetary masses (5–120 $M_E$).

The planetary potential was modelled with a cubic potential (Kley et al. 2009). We used a smoothing length of $r_{sm} = 0.6r_H$, where $r_H$ denotes the planetary Hill radius. This is the same smoothing length as in Lambrechts et al. (2014). The smoothing length has no influence on our results because the pressure bump generated by the planet outside of its orbit lies well beyond the smoothed zones. We tested different smoothing lengths, for instance, $r_{sm} = 0.8r_H$ and $r_{sm} = 0.4r_H$, and the pressure bump outside of the planetary orbit did not change compared to our standard smoothing length of $r_{sm} = 0.6r_H$.

## 2.2. Measurement of the pebble-isolation mass

Pebbles in protoplanetary discs are subject to radial drift that is due to the headwind they feel from the gas (Weidenschilling 1977; Brauer et al. 2008). The gas orbits at a slightly sub-Keplerian speed because of the force exerted by the radial pressure gradient in the protoplanetary disc. This velocity difference is expressed as

$$v_{\text{gas},\phi} = v_K(1-\eta) = v_K - \eta v_K = v_K - \Delta v, \qquad (2)$$

where

$$\eta = -\frac{1}{2} \left(\frac{H}{r}\right)^2 \frac{\partial \ln P}{\partial \ln r}, \qquad (3)$$

and $P$ is the pressure in the protoplanetary disc. In isothermal discs the pressure is given by $P = c_s^2 \rho_g$, with $c_s$ being the isothermal sound speed $c_s = H \Omega_K$ and $\rho_g$ the gas volume density.

If $\eta$ is lower than 0, the azimuthal gas velocity becomes higher than the Keplerian velocity and thus the particles feel a net outwards acceleration that will stop the inward motion of pebbles. As the planet grows, it carves a (partial) gap in the gas distribution around it by pushing material away from its orbit. This will eventually accelerate the gas outside of the planetary orbit to super-Keplerian velocities (Lambrechts et al. 2014). We calculate in the following an azimuthally averaged value of $\eta$. The azimuthally averaged $\eta$ quantity gives a good handle on the generation of the pressure bump (see Appendix A). The planetary mass at which the created pressure bump stops the radial inward flow of pebbles through radial drift is called the pebble-isolation mass $M_{\text{iso}}^{\text{iso}}$ without diffusion.

In Fig. 1 we display the $\eta$ parameter outside of the planetary orbit (the planet is fixed at $r = 1$) for several planet masses in a disc with $H/r = 0.05$, $\alpha = 0.001$ and $\Sigma_g \propto r^{-0.5}$, where $\Sigma_g$ denotes the gas surface density of the disc. A negative value of $\eta$ means that the gas velocity is super-Keplerian, and inwards-drifting pebbles are stopped and cannot reach the planet any more. The pressure gradient is calculated for an azimuthally averaged pressure in the protoplanetary disc. In this case, a pebble-isolation mass of $M_{\text{iso}}^{\text{iso}} \approx 25 M_E$ is found, in rough agreement with Lambrechts et al. (2014).
<figure><img src="image_1.png" /><figcaption>**Fig. 1.** Pressure gradient parameter $\eta$ as a function of orbital distance from the planet for $\alpha = 0.001$, $H/r = 0.05$, and different planetary masses. The location of the planet is fixed at $r = 1$. A negative $\eta$ parameter indicates the formation of a pressure bump and with it super-Keplerian gas velocities that stop the inwards drift of pebbles. Here a mass of about 25 Earth masses is needed to generate the pressure bump outside of the planetary orbit. Inside of the planetary orbit, $\eta$ is always positive.</figcaption></figure>

<figure><img src="image_2.png" /><figcaption>**Fig. 2.** Pebble-isolation mass as a function of $\alpha$ and for different aspect ratios $H/r$. The pebble-isolation mass is fitted through different $\alpha$ values with a simple fit that also scales with $(H/r)^3$.</figcaption></figure>

## 2.3. Dependence on viscosity and aspect ratio

In Fig. 2 we present the pebble-isolation mass determined by 3D hydrodynamical simulations with different viscosities and aspect ratios. We define the pebble-isolation mass as the planetary mass at which the pressure bump outside of the planetary orbit becomes large enough to turn $\eta$ negative. The pebble-isolation mass increases with $\alpha$ and with $H/r$, as predicted by Lambrechts et al. (2014). We show a fit for our obtained data that scales with the aspect ratio, similar to Lambrechts et al. (2014),

$$ M_{\text{iso}}^{\dagger} (H/r) \propto \left( \frac{H/r}{0.05} \right)^3 . \qquad (4) $$

Additionally, the fit includes a dependency on $\alpha$, which varies the pebble-isolation mass by a factor of 2–3 between low and high $\alpha$ values. Our fit is therefore also a function of $\alpha$ in the following way:

$$ M_{\text{iso}}^{\dagger} (\alpha) \propto \left( 0.34 \left( \frac{\log(\alpha_3)}{\log(\alpha)} \right)^4 + 0.66 \right), \qquad (5) $$

where $\alpha_3 = 0.001$.

In contrast, the original formula from Lambrechts et al. (2014), which is at the base of the new, refined formula proposed here, is similar to the expression of the critical mass for the wake to shock in a disc (Goodman & Rafikov 2001; Rafikov 2002). This highlights the fact that opening a deep gap around the orbit of the planet and creating a small gap that just reverses the pressure gradient are not similar processes and obey different physics.

To emphasise this effect, we calculated the depth of the gap following the formula for giant planet gap depths by Crida & Morbidelli (2007) given as

$$ G(\mathcal{P}) = \begin{cases} \frac{\mathcal{P}-0.541}{4}, & \text{if } \mathcal{P} < 2.4646 \\ 1.0 - \exp\left(-\frac{\mathcal{P}^{3/4}}{3}\right), & \text{otherwise.} \end{cases} \qquad (6) $$

The parameter $\mathcal{P}$ is given by Crida et al. (2006) as

$$ \mathcal{P} = \frac{3 H}{4 r_{\text{H}}} + \frac{50}{qR} \le 1. \qquad (7) $$

Here $q$ is the star-to-planet mass ratio, $r_H$ the planetary Hill radius, and $\mathcal{R}$ the Reynolds number given by $\mathcal{R} = r_H^2 \Omega_P / \nu$. In our simulations we have checked that the left term of Eq. (7) varies from 2 to 10 for planets that have reached the pebble-isolation mass. Thus Eq. (6) predicts a gap depth of between 15% and 60%, while our 3D simulations yield gap depths of only 10–20%. The partial gap opened at the pebble-isolation mass thus corresponds to a different regime compared to the gap opening mass for giant planets of Crida et al. (2006) and the gap depth of Crida & Morbidelli (2007). Hence, the expression of Crida & Morbidelli (2007) for the depth of the gap should not be extrapolated to the regime of large $\mathcal{P}$ (see Eq. (7)), low planet masses, and shallow gaps, and cannot be used to estimate the pebble-isolation mass.

## 2.4. Global pressure gradient

By changing the background gradient in surface density, the global pressure gradient changes. We varied the background surface density gradient $\Sigma_g \propto r^{-s}$ in the disc from $s = 0.5$ to $s = -1.5$ and determined the pebble-isolation mass for discs with $H/r = 0.05$ and $\alpha = 0.001$ with the same method as above. A global inversion of the gas surface density gradient does not already imply a pressure bump in 3D simulations, in contrast to 2D simulations, because a pressure gradient inversion in 3D isothermal discs can only be reached by an inversion of the volume density gradient, where $\rho \propto r^{s-1}$ for a radially constant $H/r$. The difference between 2D and 3D simulations regarding the pebble-isolation mass is discussed in more detail in Appendix B. In Fig. 3 we show the pebble-isolation mass as a function of the background value of $\partial \ln P / \partial \ln r$ of the unperturbed disc.

The dependence on the background gradient of surface density (and thus on the background value of $\partial \ln P / \partial \ln r$) is not very strong. We approximated the dependency of the pebble-isolation mass on $\partial \ln P / \partial \ln r$ with the expression

$$ M_{\text{iso}}^{\dagger} \left( \frac{\partial \ln P}{\partial \ln r} \right) \propto \left( 1 - \frac{\frac{\partial \ln P}{\partial \ln r} + 2.5}{6} \right), \qquad (8) $$

where the reference value $\partial \ln P / \partial \ln r = -2.5$ corresponds to an unperturbed disc with $H/r = \text{const.}$ and $\Sigma_g \propto r^{-0.5}$.
<figure><img src="image_3.png" /><figcaption>Fig. 3. Pebble-isolation mass as a function of $\partial \ln P / \partial \ln r$ of unperturbed discs with different surface density gradients. Steeper surface density slopes $s$ result in more negative $\partial \ln P / \partial \ln r$ values. All simulations have been performed for planets in discs with $H/r = 0.05$, $f = 0$ and $\alpha = 0.001$.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 4. Pressure gradient parameter $\eta$ as a function of orbital distance in a disc set-up with heating and cooling for two different planetary masses. The planet is placed at $r = 1$. The pebble-isolation mass is reached at $\approx 28.3 M_E$ in the simulations, in good agreement with Eq. (10).</figcaption></figure>

## 2.5. Flared discs

Using $\Sigma_g \propto r^s$ and $H/r \propto r^f$, one can derive a dependency of $\eta$ on the orbital distance $r$ (using $P = c_s^2 \rho = H^2 \Omega^2 \rho$),

$$ \eta = -\frac{1}{2} \left(\frac{H}{r}\right)^2 \frac{\partial \ln P}{\partial \ln r} = -\frac{1}{2} \left(\frac{H_0}{r_0}\right)^2 r^{2f} (f + s - 2). \quad (9) $$

Here $H_0/r_0$ indicates the aspect ratio at $r = 1$. This equation indicates that $\eta$ only varies radially in discs with non-constant $H/r$. As the pressure bump generated by the planet is located about $2H$ outside of the planet position, one could imagine that a change in $\eta$ with orbital distance might influence how the pressure bump is generated. We therefore tested the influence of the flaring index on the pebble-isolation mass in isothermal discs with $\alpha = 0.001$, $\Sigma_g \propto r^{-0.5}$ and $H/r = 0.05 r^f$, where $f$ spans from $-0.42$ to $+0.42$. We did not find any dependence on the pebble-isolation mass in discs with different flaring index. The pebble isolation mass in this case is determined only by the local unperturbed $\partial \ln P / \partial \ln r$ value, $\alpha$ and $H/r$ at the location of the planet, but not by the flaring of the disc itself.

## 2.6. Pebble isolation mass without diffusion

To summarise the results of Sects. 2.3–2.5, we find that the pebble-isolation mass without diffusion $M_{iso}^\dagger$ is given by

$$ M_{iso}^{\dagger} = 25 f_{\text{fit}} M_E, \quad (10) $$

where

$$ f_{\text{fit}} = \left[ \frac{H/r}{0.05} \right]^3 \left[ 0.34 \left( \frac{\log(\alpha_3)}{\log(\alpha)} \right)^4 + 0.66 \right] \left[ 1 - \frac{\frac{\partial \ln P}{\partial \ln r} + 2.5}{6} \right], \quad (11) $$

with $\alpha_3 = 0.001$.

## 2.7. Radiative simulations

In reality, discs have complex radial temperature profiles and a non-isothermal vertical structure. Different heating sources (viscous heating, stellar heating) are balanced by radiative cooling, which can alter the disc structure quite severely compared to

simple power laws (Bitsch et al. 2015a). In order to test the predictions of the pebble-isolation mass (Eq. (10)), we studied the pebble-isolation mass in discs with heating and cooling.

In the adiabatic (and radiative) case, the sound speed changes by a factor of $\sqrt{\gamma}$ compared to the isothermal sound speed. This leads to a difference in the scale height of the protoplanetary disc for the isothermal and adiabatic configuration, which are related in the following way:

$$ H_{\text{adi}} = \sqrt{\gamma} H_{\text{iso}}. \quad (12) $$

In the radiative configuration, we therefore used the adiabatic scale height to estimate the pebble-isolation mass in a radiative disc. The disc set-up was similar to before, where we now additionally included radiative cooling and viscous heating (with $\alpha = 6 \times 10^{-3}$) as described in Kley et al. (2009). At the planet location, $H_{adi,pla} = 0.0414$ and $\partial \ln P / \partial \ln r = -3.26$ (the flaring index of the disc is $f = -0.38$ at this location), which leads to a pebble-isolation mass of $M_{iso}^\dagger \approx 28.6 M_E$ according to Eq. (10). The results of our 3D simulations in discs with heating and cooling are shown in Fig. 4. Clearly, Eq. (10) matches the 3D simulations of discs with heating and cooling well.

## 2.8. Application of the new fitting formula

The pebble-isolation mass depends not only on the disc aspect ratio $H/r$, but also on the viscosity and the radial pressure gradient of the protoplanetary disc (Eq. (11)). For a fixed $H/r$, a change in $\alpha$ from $10^{-4}$ to $10^{-2}$ increases the pebble-isolation mass by a factor of $\approx 3$ (Fig. 2), while an increase in $\partial \ln P / \partial \ln r$ from $-3.5$ to $-1.5$ decreases the pebble-isolation mass by about $\approx 30\%$. Making use of Eq. (10), we calculated the pebble-isolation mass in a disc with $H/r = 0.05$ for different values of $\alpha$ and $\partial \ln P / \partial \ln r$ and show the resulting pebble-isolation mass in Fig. 5.

Clearly, high values of viscosity ($\approx 10^{-2}$) increase the pebble-isolation mass significantly, where the pebble-isolation mass can reach over $\approx 50 M_E$ for our nominal $\partial \ln P / \partial \ln r = -2.5$. Nevertheless, the strongest dependence of the pebble-isolation mass is on the disc aspect ratio $H/r$. The disc aspect ratio is determined by the heating of the disc, either through viscosity or stellar irradiation. As the disc evolves in time, the aspect ratio decreases in
<figure><img src="image_5.png" /><figcaption>Fig. 5. Pebble-isolation mass $M_{iso}^\dagger$ as a function of the pressure gradient $\partial \ln P / \partial \ln r$ and $\alpha$ in a disc with a constant $H/r = 0.05$. The two black lines mark 20 and 50 $M_E$. Clearly, higher values of viscosity and $\partial \ln P / \partial \ln r$ result in significantly higher pebble-isolation masses.</figcaption></figure>

the inner part of the disc as a result of reduced viscous heating and in the outer parts as a result of a decreasing stellar luminosity (Bitsch et al. 2015a). These effects reduce the pebble-isolation mass in time and only discs with high viscosities can maintain a high pebble-isolation mass in the outer parts of the disc, as the disc evolves in time.

### 3. Drift of small pebbles through the bump

Small particles ($\tau_f \ll 1$) are strongly coupled and move with the radial gas accretion flow, while larger particles ($\tau_f \gg 1$) are only weakly affected by gas drag. Here $\tau_f$ denotes the Stokes number of the pebbles. The acceleration of a pebble in a gas disc is given by

$$ \frac{dv_{\text{peb}}}{dt} = -\frac{GM_*}{r^3}r - 2\Delta v \Omega_K - \frac{1}{t_f}(v_{\text{peb}} - v_{\text{gas}}), \quad (13) $$

where $r$ denotes the vector between the central star and the pebble, and $t_f$ is the friction time, which is related to the Stokes number with $\tau_f = t_f \Omega_K$ and $\Delta v = \eta v_K$. The variables $v_{\text{peb}}$ and $v_{\text{gas}}$ are the pebble and gas velocities, respectively.

In Fig. 6 we show the radial and azimuthally averaged gas velocities in locally isothermal discs with embedded planets. The pressure bump generated by the planet outside of its orbit is clearly visible in the azimuthal velocity pattern, where the gas can reach speeds higher than the Keplerian value. Particles entering this pressure bump can be trapped, depending on their size.

A negative radial velocity indicates an inward flow of the gas, while a positive radial velocity indicates an outward movement of the gas. The planet generates a radial outward flow of gas close to its vicinity, but limited to the region in front of the pressure bump. This outward flow is related to the gap-opening process, where the planet pushes the material away from its orbit. This material can then move upwards to maintain hydrostatic equilibrium and again falls in from the top regions of the disc onto the planet. The same meridional flow was also observed for gap-opening giant planets (Morbidelli et al. 2014). In this region, inflowing particles might in principle be trapped as well.

After the hydrodynamical simulations shown in Sect. 2 reached an equilibrium state, we integrated the movement of test particles in the steady-state gas distribution to determine the dependence of the pebble-isolation mass on the Stokes number of the particles. The steady-state surface density profile and $\eta$ profile of the disc are shown in Appendix A, where the pressure bump is clearly visible for a 25 $M_E$ planet.

<figure><img src="image_6.png" /><figcaption>Fig. 6. Midplane gas velocities as a function of orbital distance in discs with $H/r = 0.05$, $\alpha = 0.001$, and $\Sigma \propto r^{-0.5}$. Top plot: radial velocity in units of the sound speed at $r = 1.0$. A negative velocity indicates an inward flow of the gas. The red line indicates the radial velocity of the unperturbed disc. Bottom plot: azimuthal velocity in units of the Keplerian velocity. If $v_\theta$ is larger than 1, the gas orbits super-Keplerian, indicating the pressure bump in the disc.</figcaption></figure>

We integrated the pebble trajectories in 2D planes of the protoplanetary disc, where we mainly focused on integration in the disc midplane. Pebbles injected at higher altitudes in the disc are additionally subject to vertical settling, which moves them quickly towards the midplane. Integrating the pebble trajectories in a 2D plane above the disc midplane revealed, however, that pebbles are stopped at all altitudes when the planet has reached isolation mass. This implies that at higher altitudes pebbles cannot drift through the generated pressure bump either.

In Fig. 7 we show the trajectories of pebbles in the midplane with constant $\tau_f = 1.0$ in the gas velocity field generated by planets with 10 $M_E$ and 25 $M_E$. As the 10 $M_E$ planet does not generate a pressure bump outside of its orbit (Fig. 1), the pebbles drift through towards the inner disc. The 25 $M_E$ planet generates a pressure bump in the disc (Fig. 1) where the pebbles can be trapped.

In the top panel of Fig. 8, we show the time evolution of the orbital distance of integrated pebble trajectories with different Stokes number $\tau_f$ in a disc with an embedded 25 $M_E$ planet. Pebbles with $\tau_f > 0.005$ are trapped in the pressure bump and do not drift inwards any more. Particles with $\tau_f < 0.005$ are well enough coupled with the gas to move through the pressure bump towards the system interior to the 25 $M_E$ planet. The pressure
<figure><img src="image_7.png" /><figcaption>Fig. 7. Trajectories of pebbles with $\tau_f = 1.0$ that started at (1.5; 0) moving in the gas velocity field generated by planets with $M_P = 10 M_E$ and $M_P = 25 M_E$. We plot the pebble trajectories in the rotating frame, and the planet position is marked by the black dot. The 10 $M_E$ planet does not generate a pressure bump, so the pebble drifts through all the way towards the inner boundary of the computational domain. The 25 $M_E$ planet, on the other hand, generates a pressure bump outside of its orbit and thus blocks the flow of pebbles.</figcaption></figure>

bump generated by the 25 $M_E$ planet is even quite weak (Fig. 1), which explains why small particles can still drift through.

The pressure bump generated by the planet increases in strength with planetary mass (Fig. 1), which allows the trapping of pebbles with smaller $\tau_f$ in it for higher planetary masses (bottom panel in Fig. 8 for 30 $M_E$ planet). Increasing the planetary mass by 20% reduces the Stokes number of particles that can still drift through the generated pressure bump by more than an order of magnitude. Now particles with $\tau_f > 5 \times 10^{-4}$ are trapped inside the pressure bump and can no longer reach the inner system.

Very small particles, however, with $\tau_f \approx 1 \times 10^{-4}$, are so strongly coupled with the gas that they completely follow the gas flow (Brauer et al. 2008) and are thus no longer blocked by the pressure bump located at $r = 1.13$ (see Sect. 3.1). Instead, the pebbles drift through the pressure bump, but are then caught just outside of the planetary orbit because of the radial outward flow of the gas (Fig. 6), which prevents further inward drift. For the 25 $M_E$ planet, the pebbles also drift through the radial outward flow of the gas because it was not strong enough to keep the pebbles from the inner disc (see Sect. 3.1).

### 3.1. Dependence on the radial gas velocity

The radial velocity of the gas in an $\alpha$ accretion disc is determined directly by the viscosity of the disc,

$$v_r = -\frac{3}{2} \frac{\nu}{r} = -\frac{3}{2} \frac{\alpha H^2 \Omega_k}{r} \qquad (14)$$

However, there is a strong debate in the literature about the causes of the turbulence and about the size of its magnitude

<figure><img src="image_8.png" /><figcaption>Fig. 8. Evolution of the orbital distance of pebbles as they drift through a disc with an embedded 25 $M_E$ planet (top) and an embedded 30 $M_E$ planet (bottom). Even though the planet has generated a pressure bump outside of its orbit (maximum at $r = 1.123$, see Fig. 1), small pebbles with $\tau_f < 0.005$ can drift through the pressure bump for the 25 $M_E$ planet, while the Stokes number of the particles has to be an order of magnitude smaller to drift through the pressure bump generated by the 30 $M_E$ planet.</figcaption></figure>

(Turner et al. 2014). As the drift velocity of the particles depends on the gas velocity, we investigate in this section how a change in the radial gas velocity influences particle drift through the disc. We are particularly interested in how a change in the radial gas velocity allows or hinders particles from drifting through the pressure bump generated by the azimuthal gas velocity changes induced by the planet outside of its orbit (Fig. 6). We artificially modified the radial velocity pattern to a fixed value, but kept the azimuthal gas velocity profile of a disc perturbed by a planet (Fig. 6). In this way, we mimicked the effects of different levels of turbulence without simulating discs with magneto-rotational instability (MRI) or vertical-shear instability turbulence.

In Fig. 9 we show the trajectories of pebbles embedded in discs with fixed radial gas velocities, but with azimuthal gas velocity profiles that correspond to Fig. 6 for the 25 $M_E$ planet. In the top panel we show the trajectories of pebbles in discs with radial gas velocities lower than in Fig. 8, while in the bottom panel we show the trajectories of pebbles in discs with radial gas velocities higher than in Fig. 8. Clearly, a lower radial gas velocity allows for more efficient trapping of smaller pebbles compared to a higher radial gas velocity.
<figure><img src="image_9.png" /><figcaption>**Fig. 9.** Evolution of the orbital distance of pebbles as they drift through a disc with an embedded 25 $M_E$ planet, where the radial velocity is fixed to $-1 \times 10^{-5} c_{s,0}$ (top) and $-3 \times 10^{-5} c_{s,0}$ (bottom). As the drift speed of the particles depends on the gas velocities, particles of different sizes cannot be blocked or drift through compared to the nominal case shown in Fig. 8. In particular, a slower gas flow allows for a more efficient particle trapping in the pressure bump generated by the planet.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Fig. 10.** Velocity perturbation $\Delta v$ as a function of planet mass at the location of the pressure bump. A negative $\Delta v$ value indicates that the gas speed is super-Keplerian. The blue circles mark the lowest Stokes number of particles that can be stopped at the pressure bump assuming a radial gas velocity of 1 cm s$^{-1}$ and following Eq. (16).</figcaption></figure>

The reason is that the radial pressure gradient is negative everywhere in the disc except at the centre of the outer pressure bump generated by the planet, where it is zero for a planet that has just reached pebble-isolation mass. Hence any tiny radial gas velocity can transport dust particles of any size across the pressure bump. Pebbles are only safe when the planet is more massive and the inner edge of the pressure bump has a positive radial gradient in pressure. As the radial gas velocity is determined by the viscosity in an $\alpha$-accretion disc, the movement of particles is determined by viscosity for discs with high viscosity and by drift in discs with low viscosity (de Juan Ovelar et al. 2016).

Pebbles with a (positive) terminal velocity

$$v_{r,t} = \tau_f \frac{1}{\rho} \frac{\partial \ln P}{\partial \ln r} = 2\tau_f \Delta v, \qquad (15)$$

high enough to compensate for the (negative) radial gas velocity $v_{r,g}$ can be trapped in the pressure bump. This can be reformulated as

$$\tau_f > \frac{v_{r,g}}{2\Delta v}. \qquad (16)$$

Using this equation, we can estimate the Stokes number of particles that are blocked by the pressure bump generated by the

planet by just looking at the velocity fields. The radial gas velocity is $\approx 1$ cm s$^{-1}$ in our disc model, and the measured $\Delta v$ is shown in Fig. 10. This indicates that pebbles with $\tau_f > 0.005$ should be blocked by a planet of 25 $M_E$ and particles with $\tau_f > 0.0003$ should be blocked by a planet of 30 $M_E$, in agreement with our simulations (Fig. 8). We have marked the minimal Stokes number of particles that can be blocked by the pressure bump generated by the planet with blue circles in Fig. 10 for a radial velocity of 1 cm s$^{-1}$. The lowest Stokes number depends linearly on the radial gas velocity, which is slightly different for the higher planetary masses (Fig. 6). For a 50 $M_E$ planet, the radial gas velocity is roughly $\approx 2$ cm s$^{-1}$ (Fig. 6), higher than for the 25 $M_E$ planet, because the planet influences the velocity pattern of the disc. This implies that pebbles with $\tau_f > 1.2 \times 10^{-4}$ can be stopped by a 50 $M_E$ planet, but we note that the blue dots in Fig. 10 correspond to $v_g = 1$ cm s$^{-1}$. This is also in agreement with our simulations. These results do not take diffusion into account, which we discuss in the next section.

## 3.2. Diffusion of dust particles

Our hydrodynamical simulations do not include any turbulent motion of the gas, such as those seen in magnetohydrodynamics (MHD) simulations (Bai 2017) or in simulations with the vertical shear instability (Nelson et al. 2013; Stoll & Kley 2016). These turbulent motions in the gas velocities can act on the movements of the pebbles, giving them random kicks. Several authors have considered the effects of diffusion on particles embedded in discs with planets. Paardekooper & Mellema (2006) included diffusion into the motion of dust particles in gas discs in the presence of 30 $M_E$ planets and estimated this effect to be of the order of 1%, indicating that diffusion of dust particles does not play a role in opening a gap in dust distribution of protoplanetary discs. Pinilla et al. (2016), on the other hand, studied dust filtration by giant planets in the context of transition discs. Giant planets open deep gaps in protoplanetary discs that prevent dust from drifting through. However, in their 2D simulations, the authors found that a planet of 1 $M_{Jup}$ still does not stop all dust particles and small dust grains ($\tau_f \approx 10^{-3}$) can drift through the gap generated by the planet at 20 AU. Their disc set-up would lead to a pebble-isolation mass of 27.5 $M_E$ according to Eq. (10), where
the pebble-isolation mass is reached in 2D disc simulations at lower masses (see Appendix B) than in 3D simulations.

Pinilla et al. (2016) considered turbulent mixing of dust particles, where the dust diffusivity follows the prescriptions by Youdin & Lithwick (2007), which depend on the Stokes number and the gas diffusivity (assumed to be equal to the disc viscosity). If then the pressure gradient is not steep enough, the particles can be released from the pressure bump, where particles with $\tau_f > \alpha$ are trapped. Smaller particles are diffused out of the pressure bump and dragged by the gas. Without diffusion, the trapping of particles is more efficient (Pinilla et al. 2012). This mechanism allowed the small particles to move across the pressure bump generated by the planet in Pinilla et al. (2016), while our simulations show an effective trapping of small particles as a result of the lack of diffusion. The difference of Pinilla et al. (2016) to Paardekooper & Mellema (2006) is probably related to different prescriptions of diffusion. We therefore estimate the effects of diffusion in the following for the pebble-isolation mass.

The equilibrium between radial advection of dust particles and turbulent diffusion is achieved when (as also stated in Pinilla et al. 2012)

$$v_{r,p} \rho_p - D \rho_g \frac{d\epsilon}{dr} = 0, \qquad (17)$$

where $v_{r,p}$ is the radial velocity of the pebbles, $\rho_p$ the pebble density, $\rho_g$ the gas density, $D$ the diffusion coefficient parameterised by $D = \alpha c_s H_g$, and $\epsilon = \rho_p/\rho_g$ is the dust-to-gas ratio. The radial velocity of the pebbles is given by

$$v_{r,p} = -2\tau_f \Delta v. \qquad (18)$$

To efficiently trap dust in a pressure bump despite turbulent diffusion, the equilibrium flux stated above must be obtained for a radial dust scale-length ($H_p$) that is roughly equal to the extent of the pressure bump ($H_b$). As the latter has an approximate width of one gas scale-height, this leads to

$$H_p \sim H_g. \qquad (19)$$

We estimate the scale-height of the dust in the pressure bump from Eq. (17) to obtain

$$\frac{D}{2\tau_f \Delta v} \sim H_g. \qquad (20)$$

Expanding now $D$, we arrive at

$$\frac{\alpha c_s H_g}{2\tau_f \Delta v} \sim H_g, \qquad (21)$$

which leads to

$$\tau_f \sim \frac{\alpha}{2\Pi}, \qquad (22)$$

where $\Pi = \Delta v/c_s$. Particles with Stokes numbers larger than this (Eq. (22)) are still trapped in a pressure bump, particles smaller than this can diffuse through the pressure bump. We show this critical Stokes number in Fig. 11 for a disc with $H/r = 0.05$.

We show $\Delta v/c_s$ as a function of planetary mass in a disc with $\alpha = 0.001$ and $H/r = 0.05$ in Fig. 12. We can now compare directly to Fig. 11 to determine when and how pebbles inside the pressure bump are affected by diffusion. For a 30 $M_E$ planet in a disc with $\alpha = 0.001$, the minimum Stokes number of particles that are not affected by diffusion is $\approx 0.025$, while for the 50

<figure><img src="image_11.png" /><figcaption>Fig. 11. Critical Stokes number of particles that are affected by turbulent diffusion as a function of absolute $\Delta v/c_s$ in the pressure bump and $\alpha$. Particles with Stokes numbers lower than the critical value easily diffuse and can drift through a pressure bump. The black lines denote particles with Stokes numbers of 10, 1, 0.1, 0.01, and 0.001 from left to right. The white dots correspond to different planetary masses, where we have taken the maximum value of the pressure bump generated by the planet for simulations with different $\alpha$. The planetary masses are marked with the white numbers next to the dots. The horizontal dashed black line marks our nominal $\alpha = 0.001$.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 12. $\Delta v/c_s$ as function of orbital distance in discs with $\alpha = 0.001$ and $H/r = 0.05$ and embedded planets of different masses. Note that $\Delta v/c_s = \eta/(H/r)$. The $\Delta v/c_s$ values shown here correspond to the white dots at $\alpha = 0.001$ (dashed line) in Fig. 11.</figcaption></figure>

$M_E$ planet, this particle size is $\approx 0.004$, which corresponds to the particle size stopped by the 25 $M_E$ planet in the case of no diffusion (Fig. 8). This indicates that the pebble-isolation mass for these Stokes numbers in case of diffusion is a factor of $\approx 2$ higher than in the case without diffusion for $H/r = 0.05$ and $\alpha = 0.001$, depending on $\tau_f$.

Inside the pressure bump generated by a Jupiter-mass planet in 2D simulations of a disc with $\alpha = 0.001$ and $H/r = 0.05$, the maximum $\Delta v/c_s = 0.66$. This leads to a critical Stokes number of $\tau_f = 7.5 \times 10^{-4}$, in agreement with previous studies including diffusion, see for example Pinilla et al. (2016), where very small particles can drift through the pressure bump generated by large planets. Adding the effects of particle diffusion in 2D disc simulations, Ataiee et al. (2018) found that diffusion can increase the pebble-isolation mass, in agreement with Pinilla et al. (2012) and our estimates.

The difference between the pebble-isolation mass derived from pure hydrodynamical simulations compared to simulations taking diffusion into account also depends on the level of turbulence in protoplanetary discs. For example, blocking particles
with $\tau_f = 0.001$ requires for $\alpha = 2 \times 10^{-4}$ an increase in pebble-isolation mass by about a factor of 1.5, while for $\alpha = 0.001$ an increase of much more than a factor of 2 is needed (a Jupiter-mass planets blocks pebbles with $\tau_f > 7.5 \times 10^{-4}$ in the case of diffusion).

However, the level of turbulence in discs is not very well constrained, where $\alpha$ values from $10^{-4}$ to $10^{-1}$ can be reached in simulations of the MRI (see Turner et al. 2014 for a review). Simulations with hydrodynamical instabilities in regions of the disc that are not subject to MRI-driven turbulence show $\alpha$ values of a few times $10^{-4}$, which probably sets a lower limit on turbulence (Klahr & Bodenheimer 2003; Nelson et al. 2013; Stoll & Kley 2016). Recent simulations indicate that disc winds could be the main driver of disc accretion. In these simulations, the mid-plane regions only have a very low level of turbulence. For such a low level of turbulence, the effect of turbulent diffusion on the pebble-isolation mass is quite low compared to discs with high viscosity.

The dominant particle size in the planet formation simulations of Bitsch et al. (2015b) with $\alpha = 0.0054$ presented below is of the order of 0.1 (Fig. 13), and about 0.01–0.2 when taking fragmentation into account (see below). From Fig. 11 it can be inferred that taking turbulent diffusion into account, the pebble-isolation mass for these Stokes numbers is higher than predicted in the case without diffusion by about a factor of 2.

### 3.3. Pebble-isolation mass including diffusion

The pebble-isolation mass does not only depend on $H/r$, $\alpha$ and $\partial \ln P/\partial \ln r$, but also, as shown in the previous subsection, on the turbulent diffusion of particles. The critical pressure gradient parameter $\Pi_{crit}$ to block particles of Stokes number $\tau_f$ is given by

$$ \Pi_{crit} = \frac{\alpha}{2\tau_f}. \qquad (23) $$

From our hydrodynamical simulations we can measure how $\Pi = \Delta v/c_s$ in the pressure bump generated by a planet that has already reached $M_{iso}^{\dagger}$ changes with planetary mass as

$$ \Pi = \lambda(M_p/M_E - M_{iso}^{\dagger}/M_E), \qquad (24) $$

where

$$ \lambda \approx 0.00476/f_{fit}, \qquad (25) $$

where $f_{fit}$ is defined in Eq. (11). This fit only applies to planets that have already reached the pebble isolation mass without diffusion $M_{iso}^{\dagger}$, because $\lambda$ gives the slope of the change of $\Pi$ inside the pressure bump generated by the planet, where $M_{iso}^{\dagger}$ is the minimum mass needed to invert the radial pressure gradient $\partial \ln P/\partial \ln r$ in the disc. We note that $\lambda$ is only valid until $M_p \approx 2.5 M_{iso}^{\dagger}$, when $\lambda$ changes, because the growing planet slowly transitions into the gap depth regime predicted by Crida & Morbidelli (2007). When setting $\Pi_{crit} = \Pi$, we can define the pebble isolation mass with diffusion $M_{iso}^{\dagger}$ as

$$ M_{iso} = M_{iso}^{\dagger} + \frac{\Pi_{crit}}{\lambda} M_E. \qquad (26) $$

Using Eq. (26), we study in the next section the effect of this new-found pebble isolation mass on the formation of planetary systems and the core masses of the formed planets.

## 4. Influence on planet formation

The pebble isolation mass determines the final mass of the planetary core in the pebble accretion scenario because the reduced accretion luminosity facilitates the accretion of gas (Lambrechts et al. 2014; Lambrechts & Lega 2017) and the planet can eventually grow to become a gas giant. The formation pathway of the growing planet is determined by the growth rate and size of the planetary core because it influences its gas accretion rates (Piso et al. 2015) and migration behaviour (Baruteau et al. 2014). By reaching a different pebble isolation mass, the planet can undergo a different formation pathway.

We therefore investigate in this section the influence of the pebble isolation mass on planet growth by comparing planet growth simulations to the new pebble isolation mass (Eq. (26)) with simulations with the pebble isolation mass measured by Lambrechts et al. (2014), who only inferred the dependence on $H/r$. For this we make use of the planet growth simulations including planet migration and disc evolution presented in Bitsch et al. (2015b).

### 4.1. Planet growth and migration model

The planet growth and migration model is described in great detail in Bitsch et al. (2015b), therefore we only repeat the essential points here. The planet growth and migration rates strongly depend on the disc structure. We used here the disc structure model of Bitsch et al. (2015a). This semi-analytical disc model features bumps and dips in the inner disc structure caused by transitions in the opacity at the water ice line, which can act as planet traps for low-mass planets (Bitsch et al. 2015a; Bitsch & Johansen 2016) and evolves in time. We used a disc lifetime of 3 Myr.

The growth rate of the planet depends on the pebble surface density $\Sigma_{peb}$ at the planet location of the protoplanetary disc

$$ \dot{M}_c = 2 \left( \frac{\tau_f}{0.1} \right)^{2/3} r_H v_H \Sigma_{peb}, \qquad (27) $$

where $r_H$ is the planetary Hill radius and $v_H$ the Hill speed at which the particles enter, given by $v_H = r_H \Omega_K$. In the drift-limited growth of dust particles to pebbles (Birnstiel et al. 2012), the pebble surface density depends on the pebble flux $\dot{M}_{peb}$ (Lambrechts & Johansen 2014) in the following way:

$$ \Sigma_{peb} = \sqrt{\frac{2 \dot{M}_{peb} \Sigma_g}{\sqrt{3} \pi \epsilon_P r_P v_K}}. \qquad (28) $$

Here $\epsilon_P = 0.5$ (Lambrechts & Johansen 2014). We note that the nominal pebble flux used in Bitsch et al. (2015b) was overestimated by a factor of $\approx 10$ and that we used here a modified pebble growth model presented in Bitsch et al. (2018).

After the planet has reached its pebble-isolation mass, it can contract a gaseous envelope (Piso & Youdin 2014), and as soon as the mass of the gaseous envelope is higher than the core mass of the planet, runaway gas accretion can start (Machida et al. 2010).

Growing planets interact with their natal protoplanetary disc and migrate in it. Low-mass planets do not perturb the disc significantly and migrate in type I migration, which depends mainly on the disc viscosity and on the radial gradients of surface density, temperature, and entropy (Paardekooper et al. 2011). The disc structure is therefore of crucial importance in determining the migration rates. Planets growing further (e.g. as a result of
rapid gas accretion) push the gas away from their orbit (or accrete it, Crida & Bitsch 2017) and open a gap in the protoplanetary disc and migrate in type II migration. This migration rate depends on the viscosity of the protoplanetary disc. For a review on planet migration, see for example Baruteau et al. (2014).

To calculate the torque Γ exerted by the disc on the planet, we followed Paardekooper et al. (2011) for type I migration and the viscous evolution for type II migration. The orbital migration time $t_m$ is given as

$$t_m = -\frac{J}{(dJ/dt)}, \qquad (29)$$

where J is the angular momentum of the protoplanet. Migration here is thus defined in terms of the total torque exerted on the orbit. The migration time so defined is positive when the total torque is negative. For constant eccentricity, the time taken to migrate to the centre is $t_m/2$ (Papaloizou & Larwood 2000). This factor of 2 was absent in the original simulations of Bitsch et al. (2015b) and was now added here.

## 4.2. Stokes numbers in the planet formation model

The Stokes numbers of the pebbles in our model are shown in Fig. 13 as a function of distance and time. The blue lines indicate the growth tracks of the planets shown in Fig. 14 as the planets grow and migrate. The solid lines correspond to solid accretion, which stops when the planet reaches pebble-isolation mass (marked as a dot), while the dashed lines indicates gas accretion.

Typically, the pebbles accreted by planets in our model have Stokes numbers in the range of $0.05 < \tau_f < 0.4$, as those are the Stokes numbers of the pebbles dominated by radial drift (Birnstiel et al. 2012; Lambrechts & Johansen 2014). As shown in Sect. 3.2, only small pebbles can drift through the pressure bump (with $\tau_f < 0.01$ when the planet just reached pebble-isolation mass in the absence of diffusion, $M_{iso}^\dagger$). In our model, $\alpha = 0.0054$, indicating a strong turbulent diffusion efficiency (Fig. 11), but even turbulent diffusion cannot carry the particles of the dominant Stokes number across the pressure bump because of the large pebble sizes (Fig. 13). Therefore the pebble-isolation mass is only slightly increased by the effects of turbulent diffusion in our model. This means that the bulk of the pebbles and therefore also the bulk of the solid mass is blocked by a planet that has reached the pebble-isolation mass if the particle size is dominated by radial drift. Additionally, the blocking of pebbles becomes more efficient as the planet grows and starts to accrete a gaseous envelope.

The final pebble sizes in our model were determined by radial drift alone, where we did not take the effects of fragmentation (Brauer et al. 2008; Birnstiel et al. 2011), bouncing (Zsom et al. 2010), or condensation (Ros & Johansen 2013; Schoonenberg & Ormel 2017) into account. Fragmentation and bouncing can lead to smaller pebble sizes than in the drift-limited case, while condensation of volatiles onto already existing pebbles can increase their size.

The sizes of pebbles determined by the fragmentation limit are given in Birnstiel et al. (2015) as

$$a_f = \frac{2 \sum_g v_f}{3\pi \rho_{\text{peb}} \alpha c_s}, \qquad (30)$$

where $\rho_{\text{peb}}$ is the density of the pebble itself (set to $1.5 \text{ g cm}^{-3}$ for water ice) and $v_f$ is the fragmentation velocity limit, where water-ice particles have a higher fragmentation velocity of $10 \text{ m s}^{-1}$

<figure><img src="image_13.png" /><figcaption>Fig. 13. Stokes number of the pebbles in the simulations presented in Figs. 15 and 16 at a given time t and orbital distance r in the disc. The black lines mark Stokes numbers of 0.05 to 0.4. Note that each planet growth trajectory stars at a given point in r-t and the planet moves to higher time values, indicating that the Stokes number of the pebbles accreted by the planet increases. Additionally, the planet migrates in the disc to smaller orbital distances, which also increase the Stokes number of the accreted particles. The over-plotted blue lines correspond to the growth tracks shown in Fig. 14, where the solid line corresponds to solid (pebble) accretion, the dot marks the pebble-isolation mass, and the dashed line corresponds to the gas accretion phase. The planetary growth track moves upwards to increasing time.</figcaption></figure>

than silicate grains (Gundlach & Blum 2015). When we use only this fragmentation limit for icy particles, the Stokes numbers of the particles in our disc model are 0.01–0.2, which is slightly smaller than in the drift-limited case (shown in Fig. 13). Lower Stokes numbers will result in higher pebble-isolation masses (Eq. (26)) and thus higher core masses. However, condensation at ice lines could lead to even larger particles in these regions. Future models of planet formation have to take these effects into account in order to calculate more realistic grain sizes.

In the following, we use Eq. (26) to include the effects of diffusion in the calculations of the pebble-isolation mass and to simulate the growth of planets through pebble accretion, where we use the Stokes numbers of the drift-limited solution shown in Fig. 13.

## 4.3. Growth tracks

Growth and migration depend on the structure of the protoplanetary disc, where we follow the disc model of Bitsch et al. (2015a), and use a dust metallicity of $Z_{\text{dust}} = 0.5\%$ to set the disc opacity. The disc viscosity is $\alpha = 0.0054$. The planetary growth rate depends crucially on the amount of available pebbles that can be accreted by the planet. In the remainder of the paper, we use $Z_{\text{peb}} = 1.0\%$ as in Bitsch et al. (2015b). We set the disc lifetime to be 3 Myr.

In Fig. 14 we show the growth tracks of planetary seeds starting at different locations in a disc that is already 1.5 Myr old, meaning that the planets evolve for 1.5 Myr to reach a disc lifetime of 3 Myr. For each orbital distance we ran two simulations, where the only difference was the final pebble-isolation mass, determined either by Lambrechts et al. (2014), Eq. (1), or by the new-found pebble-isolation mass, Eq. (26). This means that the initial growth is the same for the two simulations, until the planet in one simulation reaches the pebble-isolation mass and gas accretion starts. This is generally the case for simulations
<figure><img src="image_14.png" /></figure>

**Fig. 14.** Growth tracks of planets starting at several different initial positions in a disc that is already 1.5 Myr old for pebble-isolation masses $M_{\text{iso}}^{\text{L}^{14}}$ given by Lambrechts et al. (2014) in Eq. (1) and for the new-found pebble isolation mass $M_{\text{iso}}$ in Eq. (26). The initial growth phase is the same for both isolation masses, so the growth tracks diverge only when the pebble-isolation mass (Eq. (1)) is reached, marked by the dots. The solid lines indicate pebble accretion, while the dashed lines mark gas accretion.

following the pebble-isolation mass of Eq. (1), which is generally lower than Eq. (26), especially for the given disc model with $\alpha = 0.0054$.

The final core mass also determines the contraction phase of the envelope, where $\dot{M}_{\text{env.gas}} \propto M_{\text{core}}^{11/3}$ (Piso & Youdin 2014), which in turn determines how fast the planet transitions into runaway gas accretion and can then open a gap and transition into the slower type II migration phase.

In the inner parts of the protoplanetary disc, the pebble-isolation mass is low because of the low $H/r$, while the pebble-isolation mass is high in the outer parts of the disc because $H/r$ is high. For planets forming in the inner regions of the disc ($r < 15$ AU), the difference between Eq. (1) and Eq. (26) is not that large, allowing planets to arrive at similar total masses and orbital distances. However, the core masses of planets formed using Eq. (26) are higher. In the outer disc, the pebble flux is quite low, so that the differences in the pebble-isolation mass result in a slower growth of the planets, where $M_{\text{iso}}$ is determined by Eq. (26), because pebble accretion is slower than gas contraction for these pebble densities. This also results in further inward migration before the planet opens a deep gap and transitions into type II migration. In total, the differences regarding the final orbital position and the final planetary mass seem not very great. However, the differences in the core masses can be up to 30% (see Fig. 16), which is crucial for the formation of the ice giants in our solar system, which had low core masses in Bitsch et al. (2015b).

## 4.4. Global picture

We now extend the approach of the growth tracks to probe the planetary growth for starting positions of the planetary seeds from $r_0 = 0.2$ AU to 50 AU and from $t_0 = 100$ kyr to 3 Myr. In Fig. 15 we show the final total planetary mass of planets as a function of $r_0$ and $t_0$ for pebble-isolation masses following Eq. (1) (top) and Eq. (26) (bottom). The white crosses mark the growth tracks shown in Fig. 14.

<figure><img src="image_15.png" /></figure>

**Fig. 15.** Final masses of planets (total mass $M_P = M_c + M_{\text{env}}$) as a function of formation distance $r_0$ and formation time $t_0$ in the disc. Planets that are below the dark blue line have reached pebble-isolation mass and can accrete gas. All planets that are below the white line have $M_c < M_{\text{env}}$, indicating that they have undergone runaway gas accretion. The black lines indicate the final orbital distance $r_f$ of the planet, namely 0.1, 0.5, 1.0, 5.0, 10.0, and 20.0 AU. The top plot corresponds to the pebble-isolation mass found in Lambrechts et al. (2014), while the bottom plot corresponds to the pebble-isolation mass in this work (Eq. (10)). The higher pebble-isolation mass suppresses the formation of gas giants in the very outer parts of the disc at late times.

At first glance, the difference in final orbital position $r_f$ and final planetary mass $M_P$ is not that large compared for the different pebble-isolation masses, in agreement with Fig. 14. For $0.1 \text{ Myr} < t_0 < 2.0 \text{ Myr}$ and $r < 20$ AU, the formation of close-in planets that have reached the inner edge of the disc at $0.1 \text{ AU}$ is triggered. These planets form too close to the central star, so that planetary migration drives them towards the inner edge of the disc during the lifetime of the protoplanetary disc for our migration rates in discs with high viscosity. As the pebble-isolation mass is higher, the planetary cores with low-mass gaseous envelopes (that have not reached runaway gas accretion) become too large to be contained in the region of outward migration (which can only hold planets of a few Earth masses after about 1 Myr; Bitsch et al. 2015a). These planets then drift inwards as rock-dominated planets (bottom panel of Fig. 15 in contrast to the top panel of Fig. 15).

Planets forming in the outer part of the protoplanetary disc reach higher pebble-isolation masses owing to the higher aspect ratio. In the top panel of Fig. 15, the pebble-isolation mass is reached earlier, but the core masses are high enough ($\approx 10 M_E$)
to allow a transition into runaway gas accretion and thus gas giant formation. However, in the bottom panel of Fig. 15, the pebble-isolation mass is higher, which prolongs the core formation timescale and thus results in gas accretion at later stages. The overall differences in $r_f$ and $M_p$ are not very large, however.

Clearer differences can be seen with respect to the core masses of these planets (Fig. 16), where Eq. (26) delivers core masses that are about ≈30% higher, compared to Eq. (1) in our model with $\alpha = 0.0054$. Figure 2 shows the clear dependence on $H/r$ for the pebble-isolation mass, explaining why the core masses increase with orbital distance of the formed planets, because $H/r$ increases outwards in a stellar irradiated disc (Bitsch et al. 2013, 2015a). However, the 10 $M_E$ core mass line seems to be constant as a function of $t_0$ at early times in the inner regions of the disc ($r < 15$ AU) of Fig. 16: the growth time in the outer parts of the disc is much longer because of the lower pebble flux and the larger pebble scale height. The planetary growth additionally competes with planetary migration (driving the planet inwards to parts of the disc with lower $H/r$ and thus lower pebble-isolation mass) and disc evolution, where the disc aspect ratio decreases with time (Bitsch et al. 2015a).

These higher core masses are more consistent with the planetary structure of the solar system<fnref n="1" />, where Uranus and Neptune have not reached the pebble-isolation mass (Lambrechts et al. 2014) and thus stayed at 15–20 $M_E$ without accreting a large gaseous envelope. The final core mass could additionally be increased if more pebbles are available, allowing a faster growth of the planets.

# 5. Discussion

## 5.1. Planet migration versus pebble-drift speeds

Planets embedded in protoplanetary discs interact gravitationally with the disc and move through the disc (Paardekooper & Mellema 2006; Baruteau & Masset 2008; Kley et al. 2009; Paardekooper et al. 2011). The migration timescale $\tau_{mig}$ of low-mass planets is estimated in Tanaka et al. (2002) and given as

$$ \tau_{\text{mig}} = C \frac{M_{\odot}}{M_{\text{pl}}} \frac{M_{\odot}}{\Sigma_g(r_{\text{pl}})r_{\text{pl}}^2} \left(\frac{H}{r}\right)_{\text{pl}}^2 \Omega_K^{-1} \\ \approx 7.8 \times 10^5 \left(\frac{M_{\text{pl}}}{M_E}\right)^{-1} \text{yr.} \qquad (31) $$

Here, $r_{\text{pl}}$ and $(H/r)_{\text{pl}} = 0.05 = \text{const.}$ are the orbital distance and the aspect ratio at the planet location. The constant $C$ reflects the migration speed through the disc surface density profile and disc temperature profile, given by $C = 1/(2.5 + 1.7\beta_T - 0.1\alpha_\Sigma)$ (Paardekooper et al. 2011), where $\alpha_\Sigma$ is given by $\Sigma_g = \Sigma_0 r^{-\alpha_\Sigma}$ with $\Sigma_0 = 350 \text{ g/cm}^2$ and $\beta_T$ by $T \propto r^{-\beta_T}$. For our standard disc with $\alpha_\Sigma = 0.5$ and $\beta_T = 1$, the pre-factor is $C = 0.24$.

We can now compare this with the radial drift speed of particles (Brauer et al. 2008), which is given by

$$ v_{\text{d,rad,tot}} = v_{\text{r,d}} + \frac{v_{\text{r,gas}}}{1 + \tau_f^2}. \qquad (32) $$

<figure><img src="image_16.png" /><figcaption>Fig. 16. Final core masses of planets $M_c$ as a function of formation distance $r_0$ and formation time $t_0$ in the disc. The white lines correspond to core masses of 10, 15, and 20 $M_E$ (top to bottom) and are marked by white numbers. The black and blue lines have the same meaning as in Fig. 15. The top plot corresponds to the pebble-isolation mass stated in Lambrechts et al. (2014), while the bottom plot corresponds to the pebble-isolation mass in this work (Eq. (10)). Clearly, Eq. (10) allows for higher core masses, which is crucial for the formation of ice giants in our own solar system.</figcaption></figure>

The radial speed of the gas $v_{\text{r,gas}}$ in an $\alpha$ disc is estimated by Takeuchi & Lin (2002) as

$$ v_{\text{r,gas}} = -3\alpha \frac{c_s^2}{v_K} \left(\frac{3}{2} - \alpha_\Sigma\right). \qquad (33) $$

The quantity $v_{\text{r,gas}}$ that describes the radial drift of individual dust particles is given by Weidenschilling (1977) as

$$ v_{\text{r,d}} = - \frac{2\Delta v}{\tau_f + 1/\tau_f}, \qquad (34) $$

where $\Delta v$ is the maximum drift velocity, which is calculated as

$$ \Delta v = \frac{c_s^2}{2v_K} \left( \alpha_\Sigma + \frac{7}{4} \right). \qquad (35) $$

In Fig. 17 we show the radial drift speed of particles as a function of Stokes number, the radial gas velocity, and the migration speed of planets with 25 $M_E$ (corresponding directly to the pebble-isolation mass without diffusion) and 50 $M_E$ in a disc with $\alpha = 0.001$. Clearly, the particles drift faster than the planet migrates when the pebble-isolation mass is reached. Even for

<footnote n="1">¹ The findings in Bitsch et al. (2015b) and Bitsch & Johansen (2016) produced core masses around 10 $M_E$, approximately a factor 2 lower than Uranus and Neptune.</footnote>
<figure><img src="image_17.png" /><figcaption>**Fig. 17.** Total inward velocity of a dust particle (red) as a function of Stokes number $\tau_f$ and the corresponding gas velocity for a disc with $\alpha = 0.001$ (magenta). Over-plotted are also the radial speeds of a 25 and a 50 $M_E$ planet in the same discs, assuming pure type I planetary migration. Even for the 50 $M_E$ planets, particles with $\tau_f > 0.001$ drift faster than the planet and would thus accumulate at the pressure bump outside of the planetary orbit.</figcaption></figure>

the 50 $M_E$ planet, particles with $\tau_f > 0.001$ drift faster than the planet migrates, indicating that particles drifting inwards from the outer disc will be trapped in the pressure bump outside of the planetary orbit.

However, planetary migration is more complicated than the simple estimate provided in Eq. (31). The corotation torque can change the migration speed and also the direction of migration (Paardekooper & Mellema 2006; Baruteau & Masset 2008; Kley et al. 2009; Paardekooper et al. 2011). This can result in regions of the disc where the planet does not migrate at all, so-called zero-torque regions (Bitsch et al. 2013, 2014, 2015a; Baillié et al. 2015). However, the corotation torque is prone to saturation, which depends on the local viscosity of the protoplanetary discs (Paardekooper et al. 2011), where a lower viscosity allows an easier torque saturation, preventing outward migration. Even in these cases, however, the pebbles with $\tau_f > 0.001$ drift faster than the planet migrates and will thus accumulate outside of the planetary orbit in the generated pressure bump.

## 5.2. Mass loading in the pressure bump

As soon as the planet reaches the pebble-isolation mass, pebbles drifting inwards from the outer disc are stopped in the pressure bump generated by the planet. As the flux of pebbles from the outer disc continues, pebbles accumulate in the pressure bump and the pebble-to-gas ratio increases. However, an increased pebble-to-gas ratio will trigger the streaming instability (Bai & Stone 2010; Carrera et al. 2015), transforming the pebbles into planetesimals. For the streaming instability to occur, a vertically integrated pebble-to-gas ratio of a few percent is needed.

Pebbles can also accumulate in vortices generated outside of the gap carved by the planet, where they would form planetesimals. Raettig et al. (2015) showed that the accumulated pebbles could destroy a vortex in a disc, but in our case, the vortex is fed by the presence of the planet itself, which was not taken into account in their work. Auffinger & Laibe (2018) studied the linear growth regime of the streaming instability in pressure bumps in discs, and they found that streaming instability can occur within the pressure bump. The accumulated pebbles

inside the pressure bump therefore turn into planetesimals, which do not affect the gas velocities and thus do not disrupt the pressure bump. However, how the presence of a planet influences the streaming instability in a pressure bump is still subject to investigation.

Nevertheless, the current evidence occurrence of planetesimal formation inside the pressure bump before mass loading with pebbles influences the gas dynamics of the pressure bump itself. This makes the pressure bump outside of the planetary orbit an interesting candidate for subsequent planet formation.

## 5.3. Particle filtering by proto-Jupiters

The pebble accretion scenario does not only allow for fast accretion of planetary cores at large distances, it also gives potential solutions to (a) the great dichotomy between the terrestrial planets and the gas giants (Morbidelli et al. 2015), (b) the inward motion of the water ice line as the protoplanetary disc evolves in time and crosses the orbit of the Earth in less than 1 Myr (Morbidelli et al. 2016), and (c) explain the difference between the non-carbonaceous and carbonaceous meteorites through different isolated reservoirs (Kruijer et al. 2017). The solution to all these problems could be related to the growth of the Jupiter core and to the amount of pebbles that can drift past it after reaching pebble-isolation mass.

In these scenarios, the Jupiter core forms in the cold part of the protoplanetary disc ($r > r_{ice}$), where the pebbles are large. This makes the accretion very efficient because larger pebbles can be accreted more efficiently, allowing Jupiter to grow faster than the bodies in the terrestrial region (Morbidelli et al. 2015). Additionally, as soon as Jupiter reaches its pebble-isolation mass, the inward flux of large pebbles ($\tau_f > 10^{-1}$) is stopped. The small pebbles drifting through are accreted very inefficiently (unless they grow again through coagulation), slowing down the growth of the bodies in the terrestrial region.

Morbidelli et al. (2016) reported that the inward flux of icy pebbles was stopped by a proto-Jupiter that had reached pebble-isolation mass, thus fossilizing the water ice line at $\approx 3$ AU because the bulk of the material was stopped outside of the proto-Jupiter. However, a small fraction of water ice is needed to explain ordinary chondrites (and to a lesser extent, even enstatite chondrites) as they show evidence for some water alteration. The amount of water available had to be much lower than expected from solar proportion, however. This shows that these meteorites formed in a cold environment, but the availability of icy grains was somehow strongly reduced. The passage of small grains ($<10 \mu$m) through the barrier at the Jupiter pressure bump coupled with the inefficient accretion of such small grains explains these observations and is in agreement with the blocking efficiency of planets at pebble-isolation mass found in this study.

Kruijer et al. (2017) showed through meteoritic evidence that the reservoir between non-carbonaceous meteorites and carbonaceous meteorites was spatially separated in the protoplanetary disc around the young Sun at about $\approx 1$ Myr. This separation can be achieved by a growing planet that stops the inward flux of particles, which corresponds to the core of Jupiter. Kruijer et al. (2017) gave the mass of the Jupiter core as $\approx 20 M_E$. This is in agreement with the pebble-isolation mass we found here, but the exact mass at which a growing Jupiter generated a pressure bump outside of its orbit depends on the disc properties (viscosity, aspect ratio, and pressure gradient).
## 5.4. Ice giant formation

As soon as the planet reaches its pebble-isolation mass, the envelope of the planets is no longer heated by infalling pebbles, and gas accretion can start (Lambrechts et al. 2014). This initial contraction of the gaseous envelope depends on the cooling for the envelope and with that on the opacity inside it (Piso et al. 2015). At low temperatures (below 1000 K), the opacity is dominated by the dust grains, where a larger amount of dust grains increases the opacity and thus prolongs the contraction of the planetary envelope. However, in the pebble accretion scenario, the planet blocks the influx of new pebbles, and only very tiny dust grains can reach the planet.

As soon as the main flux of pebbles onto the planet is stopped, the core stops to grow, but very small grains might still enter the planetary atmosphere and thus keep the opacity high, prolonging gas envelope contraction (Lambrechts & Lega 2017). However, even in the case of no diffusion, particles with $\tau_f < 0.005$ can reach the planet (Fig. 8) and keep the opacity high. Only for larger planets can the pebble flux sufficiently be reduced. Additionally, this depends on the viscosity of the protoplanetary disc because of the diffusion of particles through the pressure bump, where higher viscosities allow a more efficient diffusion and the planet has to reach higher masses to block pebbles with the same Stokes number compared to planets in low-viscosity environments. This initial growth stage might be very important for the growth of ice giants, preventing them from immediately entering into rapid gas accretion and thus explain why Uranus and Neptune have envelopes of 10–15% of their mass without entering into runaway gas accretion.

Lambrechts et al. (2014) found that Uranus and Neptune may never have reached pebble-isolation mass and that this prevented them for accreting gas. Our findings indicate that these planets might have reached the pebble-isolation mass without diffusion limit (to stop efficient growth of the core), but the influx of small particles prevented an efficient cooling of the atmosphere and thus runaway gas accretion. However, this could only have occurred if the viscosity of the protoplanetary disc was very low ($\alpha \approx 10^{-4}$) because otherwise the pebble-isolation mass is higher than the mass of Uranus and Neptune for typically expected disc aspect ratios ($H/r > 0.04$) in the outer disc at late disc evolution stages. We note that the outer disc structure is dominated by stellar irradiation, so that it is independent of the disc viscosity and $H/r$ becomes smaller through the reduced stellar irradiation as the system ages (Bitsch et al. 2015a).

## 6. Summary

We have conducted 3D hydrodynamical simulations to measure the pebble-isolation mass as a function of the disc structure and turbulence strength. In particular, we investigated the dependence on the disc aspect ratio $H/r$, the disc viscosity (parametrised through $\alpha$), the radial pressure gradient $\partial \ln P/\partial \ln r$, and the particle size described by the Stokes number $\tau_f$. Our findings generally agree with the results presented in Lambrechts et al. (2014), who inferred the cubic dependence on the disc aspect ratio $H/r$, but we refined the pebble-isolation mass to more disc parameters and also confirmed our results in fully radiative discs with heating and cooling. In Eq. (26) we provide the pebble-isolation mass as a function of our investigated disc parameters, which is useful for planet formation simulations involving pebble accretion (Bitsch et al. 2015b; Levison et al. 2015; Chambers 2016; Matsumura et al. 2017).

Our findings result in a pebble-isolation mass that is up to a factor of 2–3 higher than found in Lambrechts et al. (2014) in the high-viscosity case ($\alpha \sim 0.01$) and a factor of $\approx 1.5$ higher than in 2D simulations, see Appendix B. A higher viscosity and a steeper radial pressure gradient both result in a higher pebble-isolation mass. For very low viscosities, our simulations match the results of Lambrechts et al. (2014). This results in higher core masses in planet formation simulations compared to previous simulations (Sect. 4) in discs with $\alpha = 0.0054$. Discs with higher viscosity thus better match the heavy element content of the ice giants in our own solar system, because even as the disc evolves and $H/r$ decreases, the pebble isolation mass stays high enough so that the ice giants did not reach pebble isolation mass. Additionally, discs with higher viscosity can more easily match overall the heavy element content of giant exoplanets (Thorngren et al. 2016).

We also investigated the penetration of particles of various Stokes number $\tau_f$ through the pressure bump by radial advection with the gas and through turbulent diffusion. In the absence of turbulent diffusion, a planet that has reached pebble-isolation mass can readily block pebbles with $\tau_f > 0.005$, while a mass higher by a factor two is necessary to block pebbles with Stokes numbers as low as $\tau_f \sim 10^{-4}$. Including turbulent diffusion of particles due to viscosity changes this picture. Depending on viscosity and the particle size, the generated pressure bump needs to be stronger. To block particles of $\tau_f = 0.05$, a typical size in drift-limited pebble growth models (Birnstiel et al. 2012), the planetary mass has to be increased by up to a factor of $\approx 2$ compared to the pebble-isolation mass without turbulent diffusion (Eq. (10)) for high-viscosity discs. In low-viscosity discs, blocking of particles with $\tau_f = 0.05$ requires a much smaller increase of the planetary mass than in the pebble-isolation mass without diffusion (Fig. 11).

*Acknowledgements.* B.B. was supported by the Knut and Alice Wallenberg Foundation (grant 2012.0150). A.J. thanks the Knut and Alice Wallenberg Foundation (grants 2012.0150, 2014.0017, 2014.0048), the Swedish Research Council (grant 2014-5775) and the European Research Council (ERC Consolidator Grant 724687-PLANETESYS) for their financial support. The computations were done on the "Mesocentre SIGAMM" machine, hosted by the Observatoire de la Côte d'Azur. We thank P. Pinilla for discussions of the diffusion of dust particles. We thank the referee John Chambers for the comments that helped to improve this work.

## References

Ataiee, S., Baruteau, C., Alibert, Y., & Benz, W. 2018, A&A, in press, DOI: 10.1051/0004-6361/201732026

Auffinger, J., & Laibe, G. 2018, MNRAS, 473, 796

Bai, X. N. 2017, ApJ, 845, 75

Bai, X. N., & Stone, J. M. 2010, ApJ, 722, L220

Baillié, K., Charnoz, S., & Pantin, É 2015, A&A, 577, A65

Baruteau, C., & Masset, F. 2008, ApJ, 672, 1054

Baruteau, C., Crida, A., Paardekooper, S. J., et al. 2014, in Protostars and Planets VI (University of Arizona Press), 667

Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2011, A&A, 525, A11

Birnstiel, T., Klahr, H., & Ercolano, B. 2012, A&A, 539, A148

Birnstiel, T., Andrews, S. M., Pinilla, P., & Kama, M. 2015, ApJ, 813, L14

Bitsch, B., & Johansen, A. 2016, A&A, 590, A101

Bitsch, B., Crida, A., Libert, A.-S., & Lega, E. 2013a, A&A, 555, A124

Bitsch, B., Crida, A., Morbidelli, A., Kley, W., & Dobbs-Dixon, I. 2013b, A&A, 549, A124

Bitsch, B., Morbidelli, A., Lega, E., & Crida, A. 2014, A&A, 564, A135

Bitsch, B., Johansen, A., Lambrechts, M., & Morbidelli, A. 2015a, A&A, 575, A28

Bitsch, B., Lambrechts, M., & Johansen, A. 2015b, A&A, 582, A112

Bitsch, B., Lambrechts, M., & Johansen, A. 2018, A&A, 609, C2
Brauer, F., Dullemond, C., & Henning, T. 2008, A&A, 480, 859

Carrera, D., Johansen, A., & Davies, M. B. 2015, A&A, 579, A43

Chambers, J. E. 2016, ApJ, 825, 63

Crida, A., & Bitsch, B. 2017, Icarus, 285, 145

Crida, A., & Morbidelli, A. 2007, MNRAS, 377, 1324

Crida, A., Morbidelli, A., & Masset, F. 2006, Icarus, 181, 587

de Juan Ovelar, M., Pinilla, P., Min, M., Dominik, C., & Birnstiel, T. 2016, MNRAS, 459, L85

Frehlik, R., & Murray-Clay, R. A. 2017, AJ, 154, 98

Goodman, J., & Rafikov, R. R. 2001, ApJ, 552, 793

Gundlach, B., & Blum, J. 2015, ApJ, 798, 34

Güttler, C., Blum, J., Zsom, A., Ormel, C. W., & Dullemond, C. P. 2010, A&A, 513, A56

Ida, S., Guillot, T., & Morbidelli, A. 2016, A&A, 591, A72

Johansen, A., & Lacerda, P. 2010, MNRAS, 404, 475

Johansen, A., & Lambrechts, M. 2017, Annu. Rev. Earth Planet. Sci., 45, 359

Johansen, A., Blum, J., Tanaka, H., et al. 2014, in Protostars and Planets VI (University of Arizona Press), 547

Johansen, A., Mac Low, M. M., Lacerda, P., & Bizzarro, M. 2015, Sci. Adv., 1, 1500109

Klahr, H. H., & Bodenheimer, P. 2003, ApJ, 582, 869

Kley, W., Bitsch, B., & Klahr, H. 2009, A&A, 506, 971

Kruijer, T. S., Kleine, T., Burkhardt, C., & Budde, G. 2017, in 48th Lunar and Planetary Science Conference, 1386

Lambrechts, M., & Johansen, A. 2012, A&A, 544, A32

Lambrechts, M., & Johansen, A. 2014, A&A, 572, A107

Lambrechts, M. & Lega, E. 2017, A&A, 606, A146

Lambrechts, M., Johansen, A., & Morbidelli, A. 2014, A&A, 572, A35

Lega, E., Crida, A., Bitsch, B., & Morbidelli, A. 2014, MNRAS, 440, 683

Levison, H. F., Thommes, E., & Duncan, M. J. 2010, AJ, 139, 1297

Levison, H. F., Kretke, K., & Duncan, M. J. 2015, Nature, 524,322

Machida, M. N., Kokubo, E., Inutsuka, S. I., & Matsumoto, T. 2010, MNRAS, 405, 1227

Matsumura, S., Brasser, R., & Ida, S. 2017, A&A, 607, A67

Morbidelli, A. & Nesvorny, D. 2012, A&A, 546, A18

Morbidelli, A., Szulágyi, J., Crida, A., et al. 2014, Icarus, 232, 266

Morbidelli, A., Lambrechts, M., Jacobson, S. A., & Bitsch, B. 2015, Icarus, 258, 418

Morbidelli, A., Bitsch, B., Crida, A., et al. 2016, Icarus, 267, 368

Nelson, R. P., Gressel, O., & Umurhan, O. M. 2013, MNRAS, 435, 2610

Ormel, C. W. & Klahr, H. H. 2010, A&A, 520, A43

Paardekooper, S.-J., & Mellema, G. 2006, A&A, 453, 1129

Paardekooper, S.-J. & Mellema, G. 2006, A&A, 459, L17

Paardekooper, S. J., Baruteau, C., & Kley, W. 2011, MNRAS, 410, 293

Papaloizou, J. C. B., & Larwood, J. D. 2000, MNRAS, 315, 823

Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, 545, A81

Pinilla, P., Klarmann, L., Birnstiel, T., et al. 2016, A&A, 585, A35

Piso, A. M. A., & Youdin, A. 2014, ApJ, 786, 21

Piso, A. M. A., Youdin, A., & Murray-Clay, R. A. 2015, ApJ, 800, 82

Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62

Raettig, N., Klahr, H. H., & Lyra, W. 2015, ApJ, 804, 35

Rafikov, R. R. 2002, ApJ, 572, 566

Rice, W. K. M., Armitage, P. J., Wood, K., & Lodato, G. 2006, MNRAS, 373, 1619

Ros, K., & Johansen, A. 2013, A&A, 552, A137

Schoonenberg, D., & Ormel, C. W. 2017, A&A, 602, A21

Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337

Simon, J., Armitage, P. J., Li, R., & Youdin, A. 2015, ApJ, 822, 55

Stoll, M., & Kley, W. 2016, A&A, 594, A57

Takeuchi, T. & Lin, D. N. C. 2002, ApJ, 581, 1344

Tanaka, H., Takeuchi, T., & Ward, W. R. 2002, ApJ, 565, 1257

Thorngren, D., Fortney, J., Murray-Clay, R. A., & Lopez, E. 2016, ApJ, 831, 64

Turner, N., Fromang, S., Gammie, C., et al. 2014, in Protostars and Planets VI, (University of Arizona Press), 411

Venturini, J., & Helled, R. 2017, ApJ, 848, 95

Visser, R. G., & Ormel, C. W. 2016, A&A, 586, A66

Weidenschilling, S. J. 1977, MNRAS, 180, 57

Youdin, A., & Lithwick, Y. 2007, Icarus, 192, 588

Zsom, A., Ormel, C. W., Güttler, C., Blum, J., & Dullemond, C. P. 2010, A&A, 513, A57
# Appendix A: Azimuthal disc structure

We present here the 2D surface density structure of a 25 $M_E$ planet embedded in a disc with $\alpha = 0.001$ and $H/r = 0.05$ (top panel of Fig. A.1) as well as the $\eta$ value in the 2D configuration. The data of these simulations correspond to the purple line in Fig. 1, indicating that the pebble-isolation mass has been reached when calculating the azimuthally averaged $\eta$ profile. However, as can be seen in the bottom panel of Fig. A.1, a negative $\eta$ value and with it a blocking of inward-drifting pebbles can be achieved at all azimuthal values when the pebble-isolation mass is reached. In principle, there is a small region in parameter space that might allow pebbles to “tunnel” through the pressure barrier, but the planet and with it the spiral waves orbit the central star at a frequency of $\Omega_p$, whereas the pebbles orbit with a frequency $\Omega(r)$, which results in the trapping of the pebbles inside the pressure bump.

<figure><img src="image_18.png" /><figcaption>**Fig. A.1.** Surface density (top) and $\eta$ value (bottom) for discs with $\alpha = 0.001$, $H/r = 0.05$ with an embedded 25 $M_E$ planet. The planet is located at $r = 1$ and $\phi = -1$. A planet with this mass opens a partial gap in the disc density and has reached the pebble-isolation mass (Fig. 1). The $\eta$ value reaches negative values outside of the planetary orbit over the whole azimuthal range, indicating that taking the azimuthally average $\eta$ value to compute the pebble-isolation mass is correct.</figcaption></figure>

# Appendix B: Comparison to 2D simulations

The pressure in isothermal 3D simulations is related to the gas volume density $P_{3D} = c_s^2 \rho_g$, while for isothermal 2D simulations, the pressure is related to the gas surface density with $P_{2D} = c_s^2 \Sigma_g$. However, the gas surface density and the gas volume density are related through

$$ \rho_g = \frac{\Sigma_g}{\sqrt{2\pi H_g}}, \qquad (B.1) $$

where $H_g$ is the gas scale height. For power-law discs in 3D with $\Sigma_g \propto r^s$ and with $H/r = \text{const.}$, this implies $\rho_g \propto r^{s-1}$ and $c_s^2 \propto r^{-1}$. In 2D discs, the surface density needs to have $s > 1$ to invert the pressure, while in 3D $s > 2$ is needed, making it harder to open a pressure bump in the protoplanetary disc.

In Fig. B.1 we show the $\eta$ value as function of orbital distance for 2D discs with $\alpha = 0.001$, $H/r = 0.05$, and different planetary masses, the same as in Fig. 1. Clearly, a much lower planetary mass is needed to generate a pressure bump outside of the planetary orbit. In 2D simulations, planetary masses of about a factor of 1.5 less are needed to open a pressure bump in the protoplanetary disc than in 3D simulations.

<figure><img src="image_19.png" /><figcaption>**Fig. B.1.** $\eta$ parameter as function of orbital distance from the planet for $\alpha = 0.001$, $H/r = 0.05$ and different planetary masses in 2D discs. The location of the planet is always at $r = 1$. A negative $\eta$ parameter indicates the pressure bump. Here a mass of about $\approx 10$ Earth masses is needed to generate the pressure bump, which is about a factor of 1.5–2 lower than in the 3D case.</figcaption></figure>