Issue 
A&A
Volume 615, July 2018



Article Number  A110  
Number of page(s)  11  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/201732026  
Published online  23 July 2018 
How much does turbulence change the pebble isolation mass for planet formation?
^{1}
University of Bern, Physics Institute, Space Research and Planetary Sciences,
Sidlerstrasse 5,
3012
Bern, Switzerland
email: sareh.ataiee@space.unibe.ch
^{2}
IRAP, Université de Toulouse, CNRS, UPS,
Toulouse, France
Received:
2
October
2017
Accepted:
24
March
2018
Context. When a planet becomes massive enough, it gradually carves a partial gap around its orbit in the protoplanetary disk. A pressure maximum can be formed outside the gap where solids that are loosely coupled to the gas, typically in the pebble size range, can be trapped. The minimum planet mass for building such a trap, which is called the pebble isolation mass (PIM), is important for two reasons: it marks the end of planetary growth by pebble accretion, and the trapped dust forms a ring that may be observed with millimetre observations.
Aims. We study the effect of disk turbulence on the PIM and find its dependence on the gas turbulent viscosity, aspect ratio, and particles Stokes number.
Methods. By means of 2D gas hydrodynamical simulations, we found the minimum planet mass to form a radial pressure maximum beyond the orbit of the planet, which is the necessary condition to trap pebbles. We then carried out 2D gas plus dust hydrodynamical simulations to examine how dust turbulent diffusion impacts particles trapping at the pressure maximum. We finally provide a semianalytical calculation of the PIM based on comparing the radial drift velocity of solids and the root mean square turbulent velocity fluctuations around the pressure maximum.
Results. From our results of gas simulations, we provide an expression for the PIM vs. disk aspect ratio and turbulent viscosity. Our gas plus dust simulations show that the effective PIM can be nearly an order of magnitude larger in highviscosity disks because turbulence diffuse particles out of the pressure maximum. This is quantified by our semianalytical calculation, which gives an explicit dependence of the PIM with Stokes number of particles.
Conclusions. Disk turbulence can significantly alter the PIM, depending on the level of turbulence in regions of planet formation.
Key words: protoplanetary disks / planetdisk interactions / submillimeter: planetary systems
© ESO 2018
1 Introduction
Pebbles – solids in the mmcm size range – are important ingredients for planet formation in protoplanetary disks as they can efficiently speed up the formation of planetary cores (e.g. Ormel & Klahr 2010; Johansen & Lacerda 2010; Lambrechts & Johansen 2012). Solids couple more or less to the disk gas depending on their size, internal density, and on gas density and temperature. The accretion of wellcoupled solids on a planetary core is inefficient as they are fully dragged by the gas across the orbit of the planet. The accretion of decoupled solids is also inefficient because of their high velocity relative to the planet. But marginally coupled solids can be efficiently accreted as the drag force around the planetary core notably reduces their velocity relative to the planet (Ormel & Klahr 2010). This is the case for solids in the pebble size range in typical regions of planet formation.
Another interesting aspect of pebbles is their potential observational signatures in the (sub)mm continuum emission of protoplanetary disks, in particular when they are trapped at pressure maxima. These pebbles can form bright and dark emission rings when there is a radial pressure bump (e.g. Andrews et al. 2016) or lopsided emission structures if there is an azimuthal pressure bump such as a vortex in the disk gas (e.g. Fuente et al. 2017). Observing the cold dust emission of disks can give constraints about the physical process behind the dust trapping in these pressure bumps such as the possible existence of a planet (Regály et al. 2012; Akiyama et al. 2016).
The pebble isolation mass (PIM) is defined as the minimum planetary mass that prevents pebbles from drifting to and being accreted by a planet in a protoplanetary disk. As the planet grows, the spiral wakes that it generates in the protoplanetary disk progressively become shocks closer to the planet. This process leads to the opening of a partial gap in the gas disk about the orbital radius of the planet (Rafikov 2002). A pressure maximum builds up at the outer edge of the gap of the planet, where pebbles should be trapped (e.g. Haghighipour & Boss 2003; Pinilla et al. 2012). This trapping may have two major implications: the growth of the planet by pebble accretion ceases and the trapped pebbles can produce ringlike structures in the (sub)mm continuum emission. Said differently, the PIM can be seen as the maximum mass of a planetary core for insitu pebblemediated planet formation, and as the minimum planetary mass that may produce observable ringlike features in observations of the radio continuum.
An accurate determination of the PIM, particularly its functional dependence on disk parameters, is especially important for global models of planet formation and evolution. Lambrechts et al. (2014) investigated the dependency of the PIM on the scale height of the disk and location of the planet, but not on turbulent viscosity of the disk. However, the formation of even a partial gap in the disk gas around the orbit of the planet is sensitive to the turbulent viscosity of the disk, as shown for instance in the gapopening criterion formulated in Crida et al. (2006). Furthermore, the presence of a pressure maximum is a necessary condition for trapping solids that are marginally coupled to the gas beyond the orbit of the planet, but it is not sufficient. Turbulent diffusion also acts on solids and may kick them out of the equilibrium location defined by the pressure maximum. This depends on how strong turbulent diffusion is compared to the radial drift due to gas drag. Because turbulent diffusion depends both on the turbulent viscosity of the disk and the size of the solids via their stopping time (Youdin & Lithwick 2007), the PIM should therefore also depend on solids size. We note that very small particles, which are well coupled to the gas, are not trapped at the pressure maximum (e.g. Zhu et al. 2012) unless perhaps the planet is so massive that it keeps the gas outside the orbit of the planet from flowing across the gap.
The aim of this study is to examine how the turbulent viscosity of the disk impacts the PIM with a combination of hydrodynamical simulations and semianalytical calculations. The manuscript is organized as follows: in Sect. 2, we describe our physical model and numerical setup, and we present the results of simulations. Our gas plus dust simulations show that the inclusion of dust turbulent diffusion can increase the PIM by nearly an order of magnitude at high gas turbulent viscosities. In Sect. 3, we derive a semianalytical expression for the PIMthat is valid for a broad range of disk scale heights, viscosities, and pebble sizes, and we discuss its accuracy. This expression can easily be used in population synthesis models of planet formation and evolution, or for the interpretation of radio disk observations. We finally summarize our results in Sect. 4.
2 Hydrodynamical simulations
The PIM is a priori a function of the disk’s scale height H, its turbulentkinematic viscosity ν, and the particles stopping time t_{s}. These quantities can be expressed as dimensionless parameters that are the disk’s aspect ratio h, α turbulent viscosity, and Stokes number St via (1)
where r denotes the distance from the central star and c_{s} is the sound speed. The Stokes number in the disk midplane is defined as the ratio between the particles stopping time t_{s} (the time particles need to adjust their motion to the gas because of drag forces), and the eddy turnover time t_{eddy}, that is the correlation timescale of turbulent fluctuations, which is typically of the order of the orbital timescale (see e.g. Youdin & Lithwick 2007). For particle sizes much smaller than the molecular meanfree path (Epstein regime), which is the case in our semianalytic calculation and our 2D simulations, the stopping time is t_{s} = πρ_{s}s∕2ΣΩ with ρ_{s} the particles internal density, s their physical radius, Σ and Ω the surface density and the angular velocity of the gas at the location of the particles (Weidenschilling 1977). Adopting t_{eddy} = Ω^{−1} results in St = πsρ_{s}∕2Σ.
To examine how disk turbulence impacts the PIM, we adopted a twostep strategy. In the first step, we carried out series of 2D gas hydrodynamical simulations of diskplanet interactions for a range of values for h and α. More specifically, for each pairof {h, α}, we found the minimum planet mass for which the gas azimuthal velocity (v_{φ}) can reach or slightly exceed the Keplerian velocity (v_{K}) beyond the orbital radius of the planet due to the local deposition of angular momentum by the outer wake of the planet. This is the necessary condition for pebbles to be trapped beyond the orbit of the planet. We provide a fitting formula for this PIM as a function of h and α. In this step, the turbulent viscosity modifies the gas velocities through the gas momentum equation and results in different PIM for different viscosities. That means that in this step we investigated the indirect effect of turbulent viscosity on particle trapping. In the second step, we carried out a few 2D gas plus dust hydrodynamical simulations of diskplanet interactions with a range of dust sizes to investigate how much the turbulent diffusion of dust alters the PIM as inferred from gasonly simulations. In this step, we examined the direct influence of viscous turbulent diffusion on the pebble trapping in addition to the indirect effect that is measured in the first step. These gas plus dust simulations are carried out for a specific value of h, but for a range of α values. After describing the physical model and numerical setup of our simulations in Sect. 2.1, we present the results of the gas simulations in Sect. 2.2 and those of the gas plus dust simulations in Sect. 2.3.
2.1 Physical model and numerical setup
We carried out our gas and gas plus dust hydrodynamical simulations with the code Dusty FARGOADSG. It is an extended version of the public code FARGOADSG (Masset 2000; Baruteau & Masset 2008a,b), which includes dust (Baruteau & Zhu 2016).
Gas – The gas continuity and momentum equations are solved in a 2D grid with polar coordinates {r, φ}. We use 666 grid cells uniformly spaced in radius between 0.5a and 2.5a (a denotes the semimajor axis of the planet), and 600 cells evenly spaced in azimuth between 0 and 2π. The initial surface density profile of the gas is with M_{⋆} the mass of the central star. This corresponds to a gas surface density ~900 g cm^{−2} at 1 au for a solarmass star. The ToomreQ parameter never falls below ≈25 throughout our series of simulations, and for this reason gas selfgravity is discarded. Furthermore, for simplicity we adopt a locally isothermal equation of state with the gas pressure given by , where the sound speed (or temperature) are fixed in time. The disk aspect ratio, h = c_{s}∕v_{K}, is assumed to be uniform, which translates to a temperature profile T ∝ r^{−1}. The aspect ratio h is varied from 0.03 to 0.06 by steps of 5 × 10^{−3}. The α turbulent viscosities taken in this work are 5 × 10^{−4}, 10^{−3}, 2 × 10^{−3}, 5 × 10^{−3}, and 10^{−2}. This range of h and α values are meant to reflect the range of temperatures and rates of turbulent transport expected in the inner regions (r ≲ 20 AU) of protoplanetary disks (e.g. Simon et al. 2015).
Dust – When dust is included in the simulations, it is modelled as Lagrangian test particles that feel the gravity of the star and planet, as well as the gas drag. We do not consider the dust backreaction on the gas in this study. This effect might be important in longterm studies as we briefly discuss in Sect. 4. Turbulence is modelled by applying stochastic kicks to the particles position following the method in Charnoz et al. (2011), in which the spatial dependence of the turbulent diffusion of the dust is discarded, however. Specifically, the radial kicks have mean ⟨δr⟩ = D_{d} dt ∂_{r}Σ∕Σ and standard deviation , and the azimuthal kicks have mean ⟨δφ⟩ = D_{d} dt ∂_{φ}Σ∕(r^{2}Σ) and standard deviation σ_{φ} = σ_{r}∕r. For the dust’s turbulent diffusion D_{d}, we use from Youdin & Lithwick (2007). When the particles stopping time becomes shorter than the hydrodynamical timestep set by the gas CFL condition in the simulations, the socalled short friction time approximation is used to update particles velocities (Johansen & Klahr 2005). In this work, we use 50 000 dust particles with a size distribution n(s) ∝ s^{−1} for sizes s between 1 mm and 300 mm. The particles internal density is set to 2 g cm^{−3}. These particle sizes and internal density correspond to initial Stokes numbers St ∈ [3 × 10^{−4} − 0.2]. The Stokes number of the particles changes as they drift or diffuse in the disk inversely proportionally to the gas surface density at the particles location. The size range of the particles is chosen such that, first, there are enough particles with St > 0.05 that are trapped at the edge of the gap of the planet and, second, there are still particles with St < 0.05 that are not trapped at the gap edge and drift to the inner disk. The St ≈ 0.05 threshold is meant to be a slightly lower bound for the typical range of Stokes numbers for which pebble accretion is most efficient, which is St ∈ [0.1 − 1] (see e.g. Baruteau et al. 2016). Said differently, since this study deals with the PIM, we do not investigate the trapping of particles with a smaller Stokes number, as they should not contribute significantly to pebble accretion (e.g. Lambrechts & Johansen 2014). However, one can estimate the PIM for smaller Stokes numbers from the semianalytical formula presented in Sect. 3.
Planet – The planet is held on a fixed circular orbit in all our simulations, and its gravitational potential is smoothed over a softening parameter ϵ = 0.4H(a). Its mass is gradually increased over the first five orbits in the simulations to let the gas adjust to the insertion of the planet in the computational grid. To avoid reflections of the planet wakes near to the radial edges of the grid, socalled wavekilling zones are used in which gas fields are damped towards their initial radial profile (de ValBorro et al. 2006). All calculations are carried out in a frame corotating with the planet.
Units and running times – Results of simulations are expressed in the following units: the mass unit is the mass of the central star (M_{⋆}), the length unit is the planet’s (fixed) semimajor axis (a), and the time unit is . We denote by q the planettostar mass ratio and by q_{⊕} the EarthtoSun mass ratio (q_{⊕}≈ 3 × 10^{−6}). The simulations with only gas, which are presented in Sect. 2.2, are run over t_{run} = 5000, 4000, 3000, 2000 or 1000 planet orbits for our models with lowest to highest viscosity. These running times are chosen such that the disk profiles reach a steady state. Simulations with gas plus dust, which are presented in Sect. 2.3, are restarts of the gasonly simulations, where dust particles are inserted between r = 1.4a and 1.41a uniformly in radius and azimuth.
2.2 Results of gasonly simulations: minimum planet mass to form a local pressure bump beyond the planet’s orbit
Our gas hydrodynamical simulations were carried out to find the minimum planet mass for which a radial pressure maximum forms beyond the orbit of the planet as a function of disk aspect ratio h and α turbulent viscosity. Asalready mentioned at the beginning of Sect. 2, this was carried out by looking for the minimum planettostar mass ratio (q_{g}) for which the gas azimuthal velocity slightly exceeds the Keplerian velocity beyond the planet’s orbit, for each pair of {h, α}; the subscript g in q_{g} indicates that this is the PIM inferred from the gasonly simulations without considering the effect of turbulent diffusion of the dust. Equivalently, since with η = −h^{2}∂logP∕∂logr, another way of looking for q_{g} is to find the planettostar mass ratio from which η cancels out beyond the orbit of the planet. Figure 1 shows the azimuthallyaveraged radial profile of η for α = 5 × 10^{−4}, obtained when the planet mass reaches the PIM for different values of the disk’s aspect ratio h. As can be seen in the figure, η > 0 everywhere in the disk except at the location of the radial pressure maximum, where η ≈ 0.
The values of q_{g} obtained for all of our models are shown in Fig. 2. The quantity is displayed as a function of α to highlight a linear relationship between both quantities. We find the best agreement for the following expression: (2)
Although this relationship shows that the PIM has a much stronger dependency on h than α, the gas turbulent viscosity can notably affect the PIM since α can take a large range of values in regions of planet formation. For example, in a disk with h = 0.04 near to the location of the planet, the PIM increases from ~6 to ~20 Earth masses when α increases from 5 × 10^{−4} to 10^{−2}.
Equation (2) can be regarded as the condition for the opening of a partial gap in the gas around the orbit of the planet. This certainly differs from the gapopening criterion formulated in Crida et al. (2006), which assumes planets carve a gap with a gas density drop of about 90%. The minimum planettostar mass ratio (q_{gap}) that satisfies the gapopening criterion of Crida et al. (2006), (3)
is given by Eq. (10) in Baruteau et al. (2014), which is written as (4)
We find thatq_{g}∕q_{gap} varies from about 3 to 12% for our range values of h and α. This point is further accentuated in Fig. 3, which shows that the gas surface density depletion in our models is between 10 and 20%. From this we can take away that a planet that builds up a 20% drop in the gas surface density about its orbit has reached the PIM.
We point out that Eq. (2) implies small PIMs for lowviscosity disks. According to Eq. (2), for α = 0, the PIM ranges from ~1.5 Earth masses for h = 0.03 to ~12.5 Earth masses for h = 0.06. This is consistent with the results of Duffell & MacFadyen (2013) that basically every nonmigrating planet is able to carve a partial gap or dip about its orbit in an inviscid disk.
Fig. 1
Illustration of the determination of the PIM via gasonly hydrodynamical simulations, where we find the minimum planettostar mass ratio (q_{g}) for the formation of a pressure maximum beyond the orbit of the planet. The quantity η = −h^{2} ∂logP∕∂logr (averaged over azimuth) is shown as a function of orbital radius r (normalized to the planet’s semimajor axis a). Each curve is for a pair {q_{g}, h} at α = 5 × 10^{−4}. The shaded black area shows where η ≈ 0 for the practical determination of q_{g}. 
Fig. 2
Calculation of the PIM via gasonly hydrodynamical simulations (q_{g}, PIM normalized to the star mass). Results are shown for the range of disk aspect ratios h and α turbulent viscosities considered in this work. The blue dashed curve shows the best fit for vs. α, given by Eq. (2). 
Fig. 3
Gas surface density profile in all the gasonly hydrodynamical simulations where the planet takes its PIM. Line transparency varies with the α turbulent viscosity, and colour varies with the disk aspect ratio (see legend). 
Fig. 4
Results of gas plus dust hydrodynamical simulations without dust turbulent diffusion, which show that a planet with a planettostar mass ratio q = q_{g} (the normalized PIM that is inferred from our gasonly simulations) does trap solid particles with St > 0.05 at the pressure maximum beyond its orbit. 
2.3 Results of gas plus dust simulations: importance of turbulent diffusion on dust trapping at the pressure bump
In the previous section, we used gasonly hydrodynamical simulations to obtain the minimum planettostar mass ratio (q_{g}) for a planet to form a pressure maximum beyond its orbit, which is a necessary condition for trapping pebbles. To ensure that planets with q = q_{g} trap pebbles in the absence of dust turbulent diffusion, we performed several gas plus dust simulations without dust turbulent diffusion for h = 0.04 and all our α values. Figure 4 demonstrates the trapping of particles with St ≳ 0.05 at the pressure maximum beyond the orbit of the planet, as expected. We note that the trapping location increases with decreasing the alpha turbulent viscosity, which is due to planets carving wider gaps in less viscous disks (e.g. Crida et al. 2006). The particles with St ≪ 0.05 are well coupled to the gas and have not had enough time to drift significantly over the simulation time.
While the formation of a radial pressure bump is a necessary condition to trap pebbles beyond the orbit of the planet, it is not sufficient. We also needed to make sure that the turbulent diffusion of the dust does not diffuse particles out of the pressure bump. Because of our large parameter space, finding the PIM via gas plus dust hydrodynamical simulations for every single set of {h, α, St} would be computationally expensive. We thus chose to focus on one particular aspect ratio, h = 0.04, and find the minimum planettostar mass ratio that effectively traps dust particles with St > 0.05 at the pressure bump. Those particles have a size s > 10 mm for all our α turbulent viscosities. As already argued in Sect. 2.1, the St = 0.05 threshold can be seen as a lower bound for the typical range of particle Stokes numbers for which pebble accretion is efficient.
For each value of α, we start from q = q_{g} (the PIM as inferred from the gasonly simulations presented in Sect. 2.2) and increase the planet mass until dust trapping for the aforementioned sizes is observed for the whole running time of our gas + dust simulations. The corresponding planettostar mass ratio is denoted by q_{d}, where the dsubscript highlights that it is the mass ratio of PIM to the central star as inferred from gas plus dust simulations. Dust trapping can be quantified by counting the number of particles (N) with size s > 10 mm that remainbeyond the orbital radius of the planet over time. We denote by N_{0} the number of such particles at the beginning of the gas + dust simulations (at time t = t_{run}). If the ratio N∕N_{0} stays close to unity for the whole running time of the gas plus dust simulations, the planet has then reached its actual PIM.
The results of our gas plus dust simulations are summarized in Fig. 5. The first two rows of panels show the location of the dust particles overplotted on the relative perturbation of the gas surface density. The first row indicates q = q_{g} and the second row indicates q = q_{d}. The alpha turbulent viscosity decreases from left to right throughout the panels. The comparison between these first two rows of panels makes it clear that for q = q_{g}, most of the large particles have diffused out of the pressure bump and have mainly drifted inwards towards the star. This is further highlighted in the third row of panels, which represents the time evolution of N∕N_{0}. Clearly for the smallest viscosities, q_{d} and q_{g} are in close agreement, but the discrepancy significantly increases with increasing viscosity. For instance, for α =10^{−2}, the PIM obtained from the gas plus dust simulations is about 8.5 times larger than when obtained from the gasonly simulations. It is actually about half a Jupiter mass in the former case, and about 20 Earth masses in the latter case for a solar mass central star.
We note that the above values of q_{d} are for the trapping of St > 0.05 particles at the pressure bump. A different (in particular smaller) threshold value should result in a different value for q_{d}. Instead of finding the Stokes or sizedependence of q_{d} via gas plus dust simulations, we now turn to the use of analytical calculations, which are presented in the next section.
Fig. 5
Results of the gas plus dust hydrodynamical simulations with dust turbulent diffusion for h = 0.04. The alpha turbulent viscosity decreases from left to right. First and second rows show the relative perturbation of the gas surface density (in black and white) with dots showing the location of dust particles (the colour of the dots indicates particles size). The first row shows results for q = q_{g}, that is when the mass of the planet equals the PIM inferred from the gasonly simulations (q_{⊕} is the Earth to Sun mass ratio). The second row represents q = q_{d}, that is when the mass of the planet equals the PIM obtained from the gas plus dust simulations. In each panel the white dashed circle indicates the orbital radius of the planet. The third row of panels represent the time evolution of N∕N_{0}, which is the ratio of particles number with size >10 mm that remain outside of the orbital radius of the planet. 
3 Semianalytical formula for the pebble isolation mass
In the previous section we showed that dust turbulent diffusion could significantly alter the PIM that is inferred from the necessary condition to form a radial pressure bump beyond the orbit of the planet. As dust turbulent diffusion depends on the Stokes number (in addition to the gas turbulent viscosity), the PIM should explicitly depend on the Stokes number (in addition to the gas turbulent viscosity and aspect ratio). To examine how the PIM depends on the Stokes number, we turn to analytical and semianalytical calculations. The idea is that particles get trapped at the pressure bump outside th orbit of the planet if the radial drift velocity (v_{d rift}) of the particles inside of the pressure bump, which is positive, becomes larger than the root mean square turbulent velocity fluctuations (v_{t urb}). When such a condition is satisfied, the particles that are kicked into the gap by the turbulence drift back towards the pressure maximum. We use Eq. (37) of Youdin & Lithwick (2007) to calculate the radial profile of v_{t urb}, and their Eqs. (57) and (59) for that of v_{drift}. These equations are written as
Equating these two velocities gives (7)
This is the minimum (positive) pressure gradient required in the gap to keep the particles with Stokes number St from diffusing inwards through the gap. Equation (7) shows that for St ≫ 1, ∂ logP∕∂logr tends to , independently of St, while it tends to for St ≪ 1. We thus expect the PIM to change behaviour when varying St around St = 1, as shown below.
To find the planet mass that produces the pressure gradient satisfying Eq. (7) in the gap, we need to know how the planet modifies the gas pressure profile, or actually the gas density profile, in a locally isothermal disk model (where the planet does not alter the gas temperature). We use the density profile given by Eqs. (12)–(16) of Duffell (2015) which is, overall, in fairly good agreement with the results of diskplanet simulations for partial gapopening planets; we return to this point later in this section. This density profile can be written as Σ (r) = Σ_{0}(r) C_{Σ}(r), where Σ _{0} denotes the initial (unperturbed) density profile, and C_{Σ} is a function of q, h, α, and r through an integral function that depends on s = −∂logΣ_{0}∕∂logr and f = ∂logh∕∂logr. Assuming a steadystate viscous disk model, where the mass accretion rate Ṁ = 3πνΣ_{0} is uniform, we have s = 1∕2 + 2f (since α is taken uniform), which implies that (8)
Anticipatingthat C_{Σ}(r) and the PIM have a weak dependence on the disk’s flaring index f, as is shown below, we take f = 0 and s = 1∕2, as in our hydrodynamical simulations. Because C_{Σ}(r) does not have an easily tractable form, we use an iterative method to find the right planettostar mass ratio q such that Eqs. (7) and (8) are satisfied, that is such that (9)
For each pair {h, α}, we first find the planettostar mass ratio q_{g,A} such that a radial pressure maximum builds up beyond the orbit of the planet (∂ logC_{Σ}∕∂logr ≈ 3∕2 from Eq. (8)). Then, for each Stokes number St, we increase the planettostar mass ratio from q = q_{g,A} until Eq. (9) is satisfied in the gap. The planettostar mass ratio that we obtain, which we denote by q_{d,A}, corresponds to the PIMtostar mass ratio for the effective trapping of particles with Stokes number ≥ St beyond the planet gap. The subscript A in q_{g,A} and q_{d,A} indicates that these quantities are calculated using our semianalytical method.
The calculation of q_{g,A} with the analytic density profile of Duffell (2015) is well approximated by (10)
By comparison with Eq. (2), we see that for each pair {h, α}, q_{g,A} is smaller than its value in our hydrodynamical simulations by a factor ≈1.5. This comparison is illustrated in Fig. 6.
Figure 7 shows the results of our semianalytical calculation of the ratio q_{d,A}∕q_{g,A}. For the range of Stokes numbers, aspect ratios, and alpha turbulent viscosities that we explored, we obtain good agreement with the following expression: (11)
which, from Eq. (7), can be recast more simply as (12)
As shown by the lower panel in Fig. 7, the agreement between Eq. (11) and our results of semianalytical calculations is within ≈10%. In particular, we note that for a given pair {h, α}, q_{d,A}∕q_{g,A} decreases with increasing St for St ≪ 1, and progressively becomes independent of St for St ≳ 1, in agreement with the qualitative behaviour of ∂logP∕∂logr given by Eq. (7) and discussed above. We finally stress that it is not entirely surprising that q_{d,A}∕q_{g,A} features the pressure gradient given by Eq. (7), since it is the parameter that sets the effective trapping of particles with a given Stokes number in the presence of dust turbulent diffusion.
Figure 8 compares the PIM inferred from the gas plus dust hydrodynamical simulations presented in Sect. 2.3 with the PIM obtained with the above semianalytical calculation, and expressed via Eq. (11); we denote by q_{d} its value normalized to the mass of the central star. Results are for h = 0.04 and St ≥ 0.05 particles. We see that q_{d} and q_{d,A} are in good agreement overall, although q_{d} tends to be larger than q_{d,A} at large viscosities. For example, for α = 10^{−2}, q_{d} is about 2.9 times larger than q_{d,A}.
The difference between q_{d} and q_{d,A} can be due to the discrepancy between the analytical gap density profile of Duffell (2015) and that in our hydrodynamical simulations. To investigate this, we show in Fig. 9 the radial profiles of the gas pressure (upper panel) and the logarithmic radial pressure gradient (lower panel) obtained from our simulations (solid curves) and from (Duffell (2015); dashed curves). Curves are shown for various alpha turbulent viscosities, but for the same planet mass and disk aspect ratio. For α ∈ [5 × 10^{−4} − 2 × 10^{−3}], there is very good agreement in both the maximum value of ∂logP∕∂logr and its radial location, which explains why q_{d}∕q_{d,A} ≈ 1 for these values of α. For α = 5 × 10^{−3} and 10^{−2}, however, the agreement is not so good, the maximum in ∂logP∕∂logr being much larger and located closer to the planet in the analytical profile than in the simulations.Consequently, a smaller planet mass is needed to satisfy the condition v_{drift} ≃ v_{turb} than the PIMs obtained in the simulations, which explains why q_{d}∕q_{d,A} > 1 for these values of α. Based on this comparison, we thus suggest that Eq. (11) should be considered as a lower limit for the PIM in highviscosity disks, and a good estimate for the PIM in lowviscosity disks.
The expressions of q_{g,A} and q_{d,A} in Eqs. (10) and (11) are for a steadystate viscous disk model with flaring index f = 0. As stated above, the background pressure profile scales as r^{−3∕2} for a steadystate viscous disk model, whatever the specific choice of density and temperature profiles. However, the perturbed gap density profile (C_{Σ} in Eq. (8)) depends on the gradients of the background density and temperature profiles in a way that does not feature the background pressure gradient. The PIM should therefore also depend on the background temperature gradient (or, equivalently, on the surface density gradient, since both quantities are related in a steadystate viscous disk model). We checked with our semianalytical calculation that this dependence remains weak. This is illustrated in Fig. 10, which shows that q_{d,A} slightly decreases with increasing flaring index. This decrease is by ≈7% for f = 0.5, which corresponds to a uniform background temperature. We also verified this weak dependence by running gas plus dust simulations for a radiative disk model with background gas surface density in r^{−15∕14} and temperature in r^{−3∕7} (e.g. Bitsch et al. 2014), for which we find basically the same PIM as with our fiducial disk profiles.
Fig. 6
PIM (normalized to the Earth mass) without effects of dust turbulent diffusion. A comparison between our results of gasonly hydrodynamical simulations (solid curves, Eq. (2)) with the prediction of the analytical gap model (dashed curves, Eq. (10) of Duffell (2015);. 
Fig. 7
Semianalytical calculation showing the dependence of the PIM on dust turbulent diffusion, using the method in Sect. 3. Upper panel: ratio between the normalized PIM with dust turbulent diffusion (q_{d,A}) and that without (q_{g,A}), as a functionof the particles Stokes number. Both q_{g,A} and q_{d,A} are calculated using the analytical gap density profile of Duffell (2015). Lower panel: relative difference between the normalized PIM obtained with our semianalytical calculation (q_{d,a}) and that obtained with our fitting formula (q_{d,A fit}, Eq. (11)). 
Fig. 8
Ratio of the normalized PIM for St ≥ 0.05 particles obtained in our gas plus dust hydrodynamical simulations (q_{d}) and by the semianalytical calculation given by Eq. (11) (q_{d,A}). Results are for a disk aspect ratio h = 0.04 and for all our α values. 
Fig. 9
Radial profiles of the gas pressure (upper panel) and of the logarithmic radial pressure gradient (lower panel) from our hydrodynamical simulations (solid curves) and from the analytic gap model (dashed curves) of Duffell (2015). 
Fig. 10
Dependence of the PIM on the disk’s flaring index f, as obtained from our semianalytical calculation for a disk aspect ratio h = 0.04, turbulent viscosity α = 10^{−3} and for particles Stokes number St = 0.1. 
4 Summary
In this work we have examined how disk turbulence affects the PIM. Our findings can be summarized as follows:
 1.
First, we carried out gas hydrodynamical simulations of planetdisk interactions to find the minimum planet mass tobuild up a radial pressure bump outside the orbit of the planet, as a function of disk aspect ratio h and αturbulent viscosity.The formation of such a pressure bump is a necessary but not sufficient condition for trapping pebbles beyond the orbit of the planet. By denoting q_{g} the ratio between the PIM and the mass of the central star, our results of simulations are best reproduced by the followingexpression: . Planets that reach the PIM deplete the surface density in the disk gas around their orbit by 10 to 20%. Our semianalytical calculation based on the gap density profile of Duffell (2015) gives a PIMtostar mass ratio q_{g,A} such that . The difference between q_{g} and q_{g.A} resides in a slightly steeper gap pressure profile in Duffell (2015) than in our hydrodynamical simulations.
 2.
In addition to the presence of a radial pressure maximum, particle trapping requires that the drift velocityremains larger than the root mean square turbulent velocity fluctuations in the gap. This sets a minimum pressure gradient in the gap, which is given by Eq. (7). Our semianalytical calculation shows that with dust turbulent diffusion, the effective trapping of particles with Stokes number ≥St gives riseto a correction factor by which q_{g,A} should be multiplied to obtain the actual PIMtostar mass ratio (q_{d,A}). This correction factor is given by Eq. (11). Our final expression for the PIMtostar mass ratio is written
 3.
We also carried out gas plus dust hydrodynamical simulations with a range of particle sizes and we find that the PIM given by Eq. (13) is in good agreement with our results of simulations for lowviscosity disks. This PIM takes, however, smaller values than in the simulations for highviscosity disks (by up to a factor ~3 for α = 10^{−2}).
On the simulation side, this work is based on 2D hydrodynamical simulations of diskplanet interactions for a nonmigrating planet. Such simulations require the use of a softening length in the planet’s gravitational potential, which can be seen as a crude way to mimic the effects of a finite vertical thickness. Morbidelli & Nesvorny (2012), who performed 2D simulations with a softening length slightly larger than ours (ϵ = 0.6H instead of 0.4H), found that for α = 6 × 10^{−3} and h = 0.056 the PIM is about 50 Earth masses, while our Eq. (2) gives ~ 42 Earth masses. Another comparison can be made with Lambrechts et al. (2014), whose 3D simulations resulted in a PIM of ~20 Earth masses for α = 6 × 10^{−3} and h = 0.05, while we find ~30 Earth masses. This comparison indicates that our choice of softening length gives overall good agreement with 3D simulations. Furthermore, if the migration of the planet is considered, a possible flux of pebbles coming from the disk inside the orbit of the planet may increase the PIM, depending on how the migration timescale and the pebbles radial drift timescale compare. This particular point requires a dedicated study.
A study of the PIM based on 3D gasonly hydrodynamical simulations along with 2D integrations of particle trajectories in the disk midplane has been recently undergone by Bitsch et al. (2018). A detailed comparison of our results with those of Bitsch et al. (2018) is found in Appendix A.
On the analytical side, our semianalytical calculation of the PIM discards dustgas interactions in the spirals and the direct gravitational interaction of the planet on the dust particles. Both of these effects could alter the predicted PIMs.
The stability of the dust ring that forms by piling up and growth of the particles should also be investigated in future works. If the lifetime of the ring is short compared to the age of the disk, the observed ring in millimetre observationswould need more massive planets to sustain a stronger pressure bump. On the other hand, the streaming instability might lead to a periodic growth of the planet by frequent destruction of the outer pressure bump. Therefore, studying the longterm evolution of dust trapping at the pressure maximum of planets with masses about the PIM and in lowviscosity disks (such that the low diffusion allows dusttogas ratios to reach unity) might probably need to take into account effects of dust drag onto the gas. Dust fragmentation and growth, which also depend on the disk turbulent diffusion, could also alter the PIM. These effects should be examined in future studies.
Acknowledgements
This work has been carried out within the frame of the National Centre for Competence in Research PlanetS supported by the Swiss National Science Foundation. We are grateful to the referee for her/his constructivereport that helped clarify this paper. We would also like to thank Paola Pinilla for her useful comments.
Appendix A Comparison with Bitsch et al. (2018)
Very recently, Bitsch et al. (2018) have studied the PIM dependency on the disk aspect ratio h, the α turbulent viscosity, and the initial pressure gradient ∂logP_{0}∕∂logr using gasonly 3D hydrodynamical simulations. Upon integrating the trajectory of particles with different Stokes numbers in the midplane of the disk, these authors have shown that their PIM could indeed trap particles with St ≥ 0.01 at the outer edge of the gap of the planet in the absence of dust turbulent diffusion. Using analytical estimations, they have obtained a correction factor that should be applied to their PIM expression to account for the effects of dust turbulent diffusion. Because the PIM expression obtained in Bitsch et al. (2018) can be seen as the 3D counterpart to that of the present study, we compare here both PIM formulae.
Fig. A.1
PIM (in Earth masses) without effects of dust turbulent diffusion. A comparison between the results of our 2D hydrodynamical simulations (solid curves) and the 3D simulations (dashed curves) of Bitsch et al. (2018); for Σ ∝ r^{−0.5} and uniformaspect ratios. 
To ease the comparison, we express the formula of Bitsch et al. (2018) with our notations. We denote by q_{g3d} their PIMtostar mass ratio without dust turbulent diffusion, and by q_{d3d} that obtained when accounting for dust turbulent diffusion. Eqs. (10) and (11) of Bitsch et al. (2018) can be recast as (A.1)
with α_{3} = 10^{−3}. In calculating ∂ logP∕∂logr in Eq. (11) of Bitsch et al. (2018), we use with .
In Fig. A.1, we compare q_{g} as obtained in our 2D simulations with q_{g3d} from the 3Dsimulations of Bitsch et al. (2018). Both of these show overall good agreement and similar behaviours when varying h and α. Our PIM values are a factor 1.5 to 2 times smaller than in the 3D simulations of Bitsch et al. (2018). As shown in Bitsch et al. (2018), this can be understood by a more difficult gap formation in a 3D disk model than in 2D. Bitsch et al. (2018) have also carried out a few 2D simulations, which are directly comparable to ours. For a background surface density profile in r^{−1∕2}, h = 0.05, and α =10^{−3}, their Fig. B1 shows a PIM of about 17M_{⊕}, while we obtain 14M_{⊕} in our 2D simulations. While this is an overall good agreement, the 20% relative difference can be due to different grid resolutions (we use about twice as many grid cells in the radial direction), or perhaps to different simulation times (which are not specified in Bitsch et al. 2018). As we have experienced during preliminary tests, the gap profile can evolve on long timescales (e.g. a few thousand planet orbits for α =10^{−3}), and running times that are too short lead to overestimating the PIM.
Fig. A.2
PIM (in Earth masses) with effects of dust turbulent diffusion. A comparison between the results of our semianalytical calculation (solid curves) and those of (Bitsch et al. 2018; dashed curves) for Σ ∝ r^{−0.5}, uniform aspect ratios and two values of the turbulent viscosity of the disk: α = 10^{−2} (upper panel) and α = 5 × 10^{−4} (lower panel). The red stars indicate the PIM obtained from our 2D gas plus dust simulations (for h = 0.04). 
To account for the effects of dust turbulent diffusion on the PIM, Bitsch et al. (2018) have used a similar strategy as ours, in the sense that they also look for the planet mass in their simulations from which the root mean square turbulent velocity fluctuations of the particles (v_{turb}) become equal to their radial drift velocity (v_{drift}) within the gap. For the radial drift velocity, they have used (A.3)
which coincides with Eq. (6) in the limit St ≪ 1. For the root mean square turbulent velocity fluctuations, instead of Eq. (5), Bitsch et al. (2018) have used (A.4)
with H_{d} a radial length scale for particle trapping in the pressure maximum, which is taken equal to the gas pressure scale height H, and with the dust’s turbulent diffusion D_{d} taken equal to the gas kinematic viscosity ν. We note that our expression for D_{d} tends to ν only in the limit St ≪ 1. Also, we stress that in the same limit St ≪ 1, Eq. (A.4) gives v_{turb} ≈ αc_{s}, while our Eq. (5) gives . Eq. (5) thus leads to (much) larger turbulent velocity fluctuations than Eq. (A.4). Using Eqs. (A.3) and (A.4), the pressure gradient in the gap such that v_{turb} = v_{drift} is (A.5)
which is a factor α^{−1∕2} smaller than our expression given by Eq. (7) in the limit St ≪ 1, due to the different expressions for v_{turb} in that limit.
Using Eq. (A.5), Bitsch et al. (2018) have found that the PIM with dust turbulent diffusion is written (A.6)
which, using Eq. (A.5), can be recast as (A.7)
Interestingly, Eq. (A.7) formally resembles our expression given by Eq. (12), except mainly for the scaling with the pressuregradient (which in our case comes to the 0.7 power). However, because we have different expressions for the pressure gradient ∂logP∕∂logr, our PIM values should differ a priori.
Figure A.2 compares q_{d3d} (see Eqs. (A.6), (A.1), and (A.2)) with our expression for q_{d} (see Eq. (13)) for α = 10^{−2} (upper panel) and α = 5 × 10^{−4} (lower panel). For the regime of Stokes numbers that we are interested in, our PIM with dust turbulent diffusion is overall in good agreement with Bitsch et al. (2018) for α = 10^{−2}. However, for α = 5 × 10^{−4}, we predicted smaller PIMs than Bitsch et al. (2018), for which the PIM displays no dependency on Stokes number for that viscosity. The differences are due to the combined effects of (i) different prescriptions for calculating the particle velocities (particularly the turbulent velocity fluctuations, as shown above), (ii) different dependencies of the PIM with the pressure gradient (compare Eqs. (A.7) and (12)), and also (iii) a larger q_{g} in Bitsch et al. (2018).
A last point of comparison with Bitsch et al. (2018) is the effect of varying the initial density and temperature gradients onthe PIM. While we have specialized to steadystate viscous disk models where both gradients are related, Bitsch et al. (2018) have varied both gradients independently. In Fig. 10, we show with our semianalytical method that the PIM slightly decreases with increasing the background temperature gradient ∂ logT∕∂logr. Bitsch et al. (2018) have found with their simulations that the PIM slightly decreases with increasing the initial pressure gradient, which looks consistent with our findings, except that this dependence is found with increasing the initial density gradient for a fixed temperature profile (f = 0, see their Fig. 3), whereas no dependence of the PIM is found when varying the initial temperature gradient for a fixed surface density profile (s = 1∕2, see their Sect. 2.5). The reason for this is not totally clear, and may have to do with the steadystate assumption in our disk models. As stated before, starting with a nonsteadystate disk model implies larger simulation times for a proper estimate of the PIM.
References
 Akiyama, E., Hasegawa, Y., Hayashi, M., & Iguchi, S. 2016, ApJ, 818, 158 [NASA ADS] [CrossRef] [Google Scholar]
 Andrews, S. M., Wilner, D. J., Zhu, Z., et al. 2016, ApJ, 820, L40 [NASA ADS] [CrossRef] [Google Scholar]
 Baruteau, C., & Masset, F. 2008a, ApJ, 672, 1054 [NASA ADS] [CrossRef] [Google Scholar]
 Baruteau, C., & Masset, F. 2008b, ApJ, 678, 483 [NASA ADS] [CrossRef] [Google Scholar]
 Baruteau, C., & Zhu, Z. 2016, MNRAS, 458, 3927 [NASA ADS] [CrossRef] [Google Scholar]
 Baruteau, C., Crida, A., Paardekooper, S.J., et al. 2014, Protostars Planets VI, 667 [Google Scholar]
 Baruteau, C., Bai, X., Mordasini, C., & Mollière, P. 2016, Space Sci. Rev., 205, 77 [NASA ADS] [CrossRef] [Google Scholar]
 Bitsch, B., Morbidelli, A., Lega, E., & Crida, A. 2014, A&A, 564, A135 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bitsch, B., Morbidelli, A., Johansen, A., et al. 2018, A&A, 612, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Charnoz, S., Fouchet, L., Aleon, J., & Moreira, M. 2011, ApJ, 737, 33 [NASA ADS] [CrossRef] [Google Scholar]
 Crida, A., Morbidelli, A., & Masset, F. 2006, Icarus, 181, 587 [Google Scholar]
 de ValBorro, M., Edgar, R. G., Artymowicz, P., et al. 2006, MNRAS, 370, 529 [NASA ADS] [CrossRef] [Google Scholar]
 Duffell, P. C. 2015, ApJ, 807, L11 [NASA ADS] [CrossRef] [Google Scholar]
 Duffell, P. C., & MacFadyen, A. I. 2013, ApJ, 769, 41 [NASA ADS] [CrossRef] [Google Scholar]
 Fuente, A., Baruteau, C., Neri, R., et al. 2017, ApJ, 846, L3 [NASA ADS] [CrossRef] [Google Scholar]
 Haghighipour, N. & Boss, A. P. 2003, ApJ, 583, 996 [NASA ADS] [CrossRef] [Google Scholar]
 Johansen, A., & Klahr, H. 2005, ApJ, 634, 1353 [NASA ADS] [CrossRef] [Google Scholar]
 Johansen, A., & Lacerda, P. 2010, MNRAS, 404, 475 [NASA ADS] [Google Scholar]
 Lambrechts, M., & Johansen, A. 2012, A&A, 544, A32 [Google Scholar]
 Lambrechts, M., & Johansen, A. 2014, A&A, 572, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lambrechts, M., Johansen, A., & Morbidelli, A. 2014, A&A, 572, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Masset, F. 2000, A&AS, 141, 165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Morbidelli, A.,& Nesvorny, D. 2012, A&A, 546, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ormel, C. W., & Klahr, H. H. 2010, A&A, 520, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, 545, A81 [Google Scholar]
 Rafikov, R. R. 2002, ApJ, 572, 566 [NASA ADS] [CrossRef] [Google Scholar]
 Regály,Z., Juhász, A., Sándor, Z., & Dullemond, C. P. 2012, MNRAS, 419, 1701 [NASA ADS] [CrossRef] [Google Scholar]
 Simon, J. B., Lesur, G., Kunz, M. W., & Armitage, P. J. 2015, MNRAS, 454, 1117 [NASA ADS] [CrossRef] [Google Scholar]
 Weidenschilling, S. J. 1977, MNRAS, 180, 57 [NASA ADS] [CrossRef] [Google Scholar]
 Youdin, A. N., & Lithwick, Y. 2007, Icarus, 192, 588 [NASA ADS] [CrossRef] [Google Scholar]
 Zhu, Z., Nelson, R. P., Dong, R., Espaillat, C., & Hartmann, L. 2012, ApJ, 755, 6 [NASA ADS] [CrossRef] [Google Scholar]
All Figures
Fig. 1
Illustration of the determination of the PIM via gasonly hydrodynamical simulations, where we find the minimum planettostar mass ratio (q_{g}) for the formation of a pressure maximum beyond the orbit of the planet. The quantity η = −h^{2} ∂logP∕∂logr (averaged over azimuth) is shown as a function of orbital radius r (normalized to the planet’s semimajor axis a). Each curve is for a pair {q_{g}, h} at α = 5 × 10^{−4}. The shaded black area shows where η ≈ 0 for the practical determination of q_{g}. 

In the text 
Fig. 2
Calculation of the PIM via gasonly hydrodynamical simulations (q_{g}, PIM normalized to the star mass). Results are shown for the range of disk aspect ratios h and α turbulent viscosities considered in this work. The blue dashed curve shows the best fit for vs. α, given by Eq. (2). 

In the text 
Fig. 3
Gas surface density profile in all the gasonly hydrodynamical simulations where the planet takes its PIM. Line transparency varies with the α turbulent viscosity, and colour varies with the disk aspect ratio (see legend). 

In the text 
Fig. 4
Results of gas plus dust hydrodynamical simulations without dust turbulent diffusion, which show that a planet with a planettostar mass ratio q = q_{g} (the normalized PIM that is inferred from our gasonly simulations) does trap solid particles with St > 0.05 at the pressure maximum beyond its orbit. 

In the text 
Fig. 5
Results of the gas plus dust hydrodynamical simulations with dust turbulent diffusion for h = 0.04. The alpha turbulent viscosity decreases from left to right. First and second rows show the relative perturbation of the gas surface density (in black and white) with dots showing the location of dust particles (the colour of the dots indicates particles size). The first row shows results for q = q_{g}, that is when the mass of the planet equals the PIM inferred from the gasonly simulations (q_{⊕} is the Earth to Sun mass ratio). The second row represents q = q_{d}, that is when the mass of the planet equals the PIM obtained from the gas plus dust simulations. In each panel the white dashed circle indicates the orbital radius of the planet. The third row of panels represent the time evolution of N∕N_{0}, which is the ratio of particles number with size >10 mm that remain outside of the orbital radius of the planet. 

In the text 
Fig. 6
PIM (normalized to the Earth mass) without effects of dust turbulent diffusion. A comparison between our results of gasonly hydrodynamical simulations (solid curves, Eq. (2)) with the prediction of the analytical gap model (dashed curves, Eq. (10) of Duffell (2015);. 

In the text 
Fig. 7
Semianalytical calculation showing the dependence of the PIM on dust turbulent diffusion, using the method in Sect. 3. Upper panel: ratio between the normalized PIM with dust turbulent diffusion (q_{d,A}) and that without (q_{g,A}), as a functionof the particles Stokes number. Both q_{g,A} and q_{d,A} are calculated using the analytical gap density profile of Duffell (2015). Lower panel: relative difference between the normalized PIM obtained with our semianalytical calculation (q_{d,a}) and that obtained with our fitting formula (q_{d,A fit}, Eq. (11)). 

In the text 
Fig. 8
Ratio of the normalized PIM for St ≥ 0.05 particles obtained in our gas plus dust hydrodynamical simulations (q_{d}) and by the semianalytical calculation given by Eq. (11) (q_{d,A}). Results are for a disk aspect ratio h = 0.04 and for all our α values. 

In the text 
Fig. 9
Radial profiles of the gas pressure (upper panel) and of the logarithmic radial pressure gradient (lower panel) from our hydrodynamical simulations (solid curves) and from the analytic gap model (dashed curves) of Duffell (2015). 

In the text 
Fig. 10
Dependence of the PIM on the disk’s flaring index f, as obtained from our semianalytical calculation for a disk aspect ratio h = 0.04, turbulent viscosity α = 10^{−3} and for particles Stokes number St = 0.1. 

In the text 
Fig. A.1
PIM (in Earth masses) without effects of dust turbulent diffusion. A comparison between the results of our 2D hydrodynamical simulations (solid curves) and the 3D simulations (dashed curves) of Bitsch et al. (2018); for Σ ∝ r^{−0.5} and uniformaspect ratios. 

In the text 
Fig. A.2
PIM (in Earth masses) with effects of dust turbulent diffusion. A comparison between the results of our semianalytical calculation (solid curves) and those of (Bitsch et al. 2018; dashed curves) for Σ ∝ r^{−0.5}, uniform aspect ratios and two values of the turbulent viscosity of the disk: α = 10^{−2} (upper panel) and α = 5 × 10^{−4} (lower panel). The red stars indicate the PIM obtained from our 2D gas plus dust simulations (for h = 0.04). 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.