Diversity in the outcome of dust radial drift in protoplanetary discs
^{1}
UMIFCA, CNRS/INSU France (UMI 3386), and Departamento de
Astronomía, Universidad de Chile, Casilla 36D Santiago, Chile
email:
christophe.pinte@obs.ujfgrenoble.fr
^{2}
Univ. Grenoble Alpes, IPAG, 38000
Grenoble, France
CNRS, IPAG, 38000
Grenoble,
France
^{3}
Monash Centre for Astrophysics and School of Mathematical
Sciences, Monash University, Clayton, Vic
3800,
Australia
^{4}
School of Physics and Astronomy, University of Saint Andrews,
North Haugh, St
Andrews, Fife
KY16 9SS,
UK
Received:
11
October
2012
Accepted:
18
February
2014
The growth of dust particles into planet embryos needs to circumvent the “radialdrift barrier”, i.e. the accretion of dust particles onto the central star by radial migration. The outcome of the dust radial migration is governed by simple criteria between the dusttogas ratio and the exponents p and q of the surface density and temperature power laws. The transfer of radiation provides an additional constraint between these quantities because the disc thermal structure is fixed by the dust spatial distribution. To assess which discs are primarily affected by the radialdrift barrier, we used the radiative transfer code MCFOST to compute the temperature structure of a wide range of disc models, stressing the particular effects of grain size distributions and vertical settling. We find that the outcome of the dust migration process is very sensitive to the physical conditions within the disc. For high dusttogas ratios (≳0.01) and/or flattened disc structures (H/R ≲ 0.05), growing dust grains can efficiently decouple from the gas, leading to a high concentration of grains at a critical radius of a few AU. Decoupling of grains from gas can occur at a large fraction (>0.1) of the initial radius of the particle, for a dusttogas ratio greater than ≈0.05. Dust grains that experience migration without significant growth (millimetre and centimetresized) are efficiently accreted for discs with flat surface density profiles (p < 0.7) while they always remain in the disc if the surface density is steep enough (p > 1.2). Between (0.7 < p < 1.2), both behaviours may occur depending on the exact density and temperature structures of the disc. Both the presence of large grains and vertical settling tend to favour the accretion of nongrowing dust grains onto the central object, but it slows down the migration of growing dust grains. If the disc has evolved into a selfshadowed structure, the required dusttogas ratio for dust grains to stop their migration at large radius become much smaller, of the order of 0.01. All the disc configurations are found to have favourable temperature profiles over most of the disc to retain their planetesimals.
Key words: circumstellar matter / protoplanetary disks / stars: formation / radiative transfer / methods: analytical / methods: numerical
© ESO, 2014
1. Introduction
The evolution of dust particles embedded in circumstellar discs constitutes the first step towards planet formation (Blum & Wurm 2008; Chiang & Youdin 2010). Micrometresized dust grains collide and grow by coagulation, forming larger pebbles (millimetre and centimetresized, Beckwith et al. 2000; Dominik et al. 2007 and references therein) that may ultimately give birth to kilometresized planetesimals. Dust grains also settle towards the disc midplane and migrate inwards as a result of the conjugate actions of stellar gravity, gas drag, and radial pressure gradient (e.g. Weidenschilling 1977; Nakagawa et al. 1986; Fromang & Papaloizou 2006; Laibe et al. 2012, herafter LGM12). This results in a stratified disc structure with a population of small (micrometresized) grains remaining close to the surface and larger grains located deeper in the disc (e.g. Dullemond & Dominik 2004; Duchêne et al. 2004; D’Alessio et al. 2006; Pinte et al. 2007, 2008). However, the interplay between these two mechanisms remains unclear. For instance, how dust grains overcome the “radialdrift barrier” – i.e. the accretion of dust grains onto the central star on timescales shorter than the disc’s lifetime. Several mechanisms have been invoked to overcome this barrier, such as grain growth (e.g. Laibe et al. 2008, hereafter LGFM08; Brauer et al. 2008; Okuzumi 2009; Birnstiel et al. 2009; Zsom et al. 2010; Windmark et al. 2012a; Garaud et al. 2013) or trapping in pressure maxima (e.g. Pinilla et al. 2012; Gibbons et al. 2012).
Another mechanism has been studied by Stepinski & Valageas (1997); Youdin & Shu (2002), and Youdin & Chiang (2004) and revived by LGM12, Laibe et al. (2014), and Laibe (2014, hereafter L14) for growing dust grains. As the grains move inward, their radial motion is affected by the increasing drag force, the increasing pressure gradient, and eventually a larger size due to grain growth. The combination of these effects can lead to a variety of outcomes. In particular, for specific values of the exponents of the gas surface density (Σ(r) ∝ r^{−p}) and midplane temperature (T(r) ∝ r^{−q}), the time required for grains to reach smaller and smaller radii increases drastically^{1}. This causes dust grains to decouple from the gas phase and “pileup” in the dense, inner regions of the disc: the migration timescale becomes so long that dust particles are virtually stopped and never end up accreted onto the central object.
Growing grains may also break through the drift barrier because of the transition from Epstein to Stokes drag at high gas densities (Birnstiel et al. 2010, see their Fig. 11).
Although the kinematics of the radial motion of the grains depends on the grain size, the outcome of the radial drift motion essentially depends on the values of p and q. LGM12 and L14 derived three analytic criteria that predict whether dust grains are ultimately accreted onto the central star or not. In order for dust to pile up, the criteria are

(i)
q < 1 and Λ(r) > 1 − q for growing grains, where (1)measures the relative efficiency between growth and migration, and ϵ_{0} is the disc’s initial dusttogas ratio;

(ii)
−p + q + 1/2 < 0 for nongrowing grains in the Epstein drag regime (which corresponds to dust grains of size smaller than ≈10 cm for density conditions encountered in protoplanetary discs); and

(iii)
q < 2 / 3 for planetesimals (>1 m) in the Stokes drag regime. In the following, we use the terms “growing grains”, “nongrowing grains”, and “planetesimals” to refer to these three cases.
In this paper, we present the first results of the effect of radiative transfer on the radialdrift barrier in both the Epstein and Stokes regimes. Specifically, we determine i) which parts of the (p, q) diagram are populated with physical models; ii) whether these discs are strongly affected by the radialdrift barrier process; and iii) the key parameters that affect the stability of the dust particles with respect to radial migration.
2. Methodology
In the following, we discuss the outcome of grains that are large enough (i.e. > 10 μm) to be significantly affected by radial migration. These grains first settle toward the disc midplane before they migrate inward. LGM12’s criteria were originally derived for the migration of dust grains in the disc midplane but still remain valid if the grains are located close to the midplane, i.e. at lower altitudes than the local scale height (see Appendix A).
The criteria defined in LGM12 and L14 apply for p and q corresponding to the gas phase. Here we only consider young gasrich and radially optically thick discs (corresponding to the classical T Tauri stage). We assume that the density in the disc midplane is high enough to ensure thermal equilibrium between gas and dust and we take both T_{gas} = T_{dust} and q_{gas} = q_{dust}, and only use the notation q in the following. Disc temperature structures are computed for given density structures, i.e. for given values of p_{dust}. We assume that in the initial state of the discs, gas and dust grains are radially well mixed, implying p_{gas} = p_{dust}. We discuss deviations from this case in Sect. 3.5.
2.1. Model description
We assume an axisymmetric flared density structure with a Gaussian vertical profile ρ(r,z) = ρ_{0}(r) exp( − z^{2}/ 2 h(r)^{2}) where r is the distance from the star in the disc midplane and z the altitude above the midplane. We use powerlaw distributions for the dust surface density Σ(r) = Σ_{0} (r/r_{0})^{−pdust} and the scale height h(r) = h_{0} (r/r_{0})^{β} where h_{0} is the scale height at radius r_{0} = 100 AU and β the disc flaring exponent^{2}. We fix the density structure rather than assume hydrostatic equilibrium in order to assess the effect of each model parameter individually. In particular, the flaring parameter and scale height are free parameters that are not selfconsistently derived from the temperature calculation. The disc extends from an inner cylindrical radius R_{in} = 0.1 AU to an outer limit R_{out} = 300 AU. Since we are only interested in the disc’s temperature gradient and not its absolute value, we use a fixed central object (T_{eff} = 4000 K, R = 2 R_{sun}) and a fixed dust disc mass of M_{dust} = 10^{4} M_{sun}. The disc mass directly influences the local stopping time (or equivalently, the Stokes number) and thus the grains’ kinematics. However, the Stokes number is only a scaling factor in the equations of motion of the grains, and its value does not affect whether the grains will pile up or decouple at a finite radius. The key parameters are the dependencies of the Stokes number and the Stokes number’s rate with respect to both the grain size and distance to the central star, which are independent of the disc mass/density (LGM12; L14). In particular, the decoupling radius for growing dust grains is independent of the disc mass.
Dust grains are modelled as homogeneous and spherical particles with sizes distributed according to the power law dn(a) ∝ a^{3.5} da between a_{min} (fixed to 0.01 μm) and a_{max}. Absorption and scattering opacities, scattering phase functions, and Mueller matrices are calculated using the Mie theory.
Range of values explored for the various parameters in the grid of models computed with MCFOST.
The thickness of the dust layer results from the competition between gravity and turbulent stirring and depends on both the grain size and the turbulent viscosity coefficient α_{SS} (Shakura & Sunyaev 1973). Dust settling is implemented using a diffusionadvection equation with a variable diffusion coefficient, which has been found this to be a reasonable approximation to global MHD simulations of stratified and turbulent discs (Fromang & Nelson 2009, Eq. (26)). The degree of dust settling is changed by varying α_{SS}. For comparison, models without any settling are also included.
2.2. Computing the stability criteria
To ascertain the importance of radiative transfer on the outcome of the radial migration of dust grains, we computed a grid of disc models with the MCFOST radiative transfer code (Pinte et al. 2006, 2009) varying both the dust optical properties and the disc density distribution and geometry (Table 1). This includes the p_{dust} exponent, which is treated here as an input parameter and is not evolved with time.
MCFOST solves the continuum transfer with a Monte Carlo method. In short, the code stochastically propagates photon packets in 3D through the disc: the transport of packets is governed by successive scattering, absorption, and reemission events, which are determined by the local dust properties. The temperature structure is calculated using the immediate reemission algorithm of Bjorkman & Wood (2001) combined with a continuous deposit of energy to estimate the mean intensity (Lucy 1999). The resulting midplane temperature profiles in the outer parts of the discs (r ≳ 10 AU) are fitted by a power law T(r) ∝ r^{−q}. Radiative transfer calculations have shown that this is a good approximation of the temperature profile in the outer disc (see for instance the radial temperature profiles presented in D’Alessio et al. 1998; Dullemond 2002; Pascucci et al. 2004; Pinte et al. 2009). We additionally compute the local slope of the midplane temperature as q_{loc} = −dlog (T(r))/dlog (r) in order to determine in which regions of the discs the temperature profile can be approximated by a power law and the simple analytical criteria can be applied, and to assess how the migration will be affected when the temperature profiles deviate strongly from the power law.
2.3. Model limitations
This work is based on analytical calculations of the dust evolution in discs, which necessarily implies some simplifying assumptions:

The equations are mathematically valid only when thetemperature profile is a power law. We found that this is the case forr > 10 AU in most of the disc configurations. When the temperature profile differs from a power law, the general trend remains, but the predictions become more qualitative.

The dusttogas ratio is assumed to be constant when integrating the analytical expressions. But because dust grains migrate, this is not the case in discs.

The underlying growth model remains simple and does not consider more complex effects that have been shown to affect the evolution of the dust population in young discs (see Dullemond 2013, for a recent review): for instance, bouncing or fragmentation barrier, charging barrier or porosity evolution (Okuzumi et al. 2012).

The radial gas flow, which may be important in the very central parts of the disc was not considered here. It may also act to remove the dust particles from the disc, preventing them from being accreted onto planetesimals, as it may transport particles that are coupled to the gas inwards, even when they are not drifting.
Additionally, we do not consider the case of large dust grains that are observed at large distances from the star (see for instance Hughes et al. 2008; Andrews et al. 2009, 2011; Isella et al. 2009, who measured critical radii^{3} ranging from ≈15 to ≈200 AU). The presence of millimetre sized grains at distance much greater than a few 10 AU requires additional mechanisms, such as local pressure maxima to trap the dust grains (e.g. Pinilla et al. 2012).
Fig. 1 Temperature exponent q as a function of the surface density exponent p_{dust}. Large bullets show the median values, smaller bullets the 1st and 2nd quartiles, and vertical bars the complete distribution of models; i.e., the horizontal ticks on each bar represent the minimum and maximum values reached in the whole set of models. Models with maximum grain sizes ranging from 1 μm to 10 μm are shown in red, while models with a_{max} ranging from 100 μm to 10 mm are shown in blue. The different models are slightly shifted along the p axis for clarity. The plane is divided by the stability criteria for dust grains: q = 1 with growth, −p_{gas} + q + 1/2 = 0 (assuming p_{gas} = p_{dust}, black dashed line) without growth and q = 2 /3 for planetesimals (black dotdashed lines). The thick gray line indicates the location of steadystate, constant α viscous disc models p_{gas} = 3/2 − q. 

Open with DEXTER 
3. Results and discussion
Fig. 2 Stability criteria Λ(r) + q − 1 > 0 for the growing grains as a function of the surface density exponent for r = 10 AU (left panel) and r = 100 AU (right panel) and an initial dusttogas ratio ϵ_{0} = 0.01 in both cases. Models with scale height of 5, 10, and 15 AU at a radius of 100 AU are shown in blue, black, and red respectively. The vertical axes are different in both panels. The second part of the stability criteria for growing grains: q < 1 is always satisfied in the disc models we have explored, see Fig. 1. 

Open with DEXTER 
Fig. 3 Relative decoupling radius for growing dust grains (i.e. decoupling radius divided by initial radius of the grains in the disc) as a function of Λ for various values of p (left) and q (right). The red thick solid line shows a typical Classical T Tauri Star disc with p = 1, q = 0.5. 

Open with DEXTER 
The values of q obtained from radiative transfer modelling for the various disc models are presented in Fig. 1.
3.1. Stability of growing dust grains relative to radial migration
Initially small, interstellarlike dust grains (≲1 μm) experience simultaneous migration and growth in the outer disc, until they reach a size of at least a few millimetres, when the sticking probabilities during collisions are strongly reduced (“bouncing barrier”, see Sect. 3.2). For these growing dust grains, the criterion q < 1 is satisfied for every disc’s configuration (see Fig. 1). If this first requirement were not fulfilled, growing grains would not pile up in discs. Figure 2 shows the values of Λ(r) + q − 1 for an initial gasdustratio of ϵ_{0} = 0.01. The pileup of growing grains mainly depends on the value of ϵ_{0}: higher values of ϵ_{0} favour the pileup, which becomes systematic for ϵ_{0} ≳ 0.02; i.e., Λ(r) is always larger than q − 1. The case ϵ_{0} = 0.01 results in values of Λ + q − 1 close to 0, i.e. at the limit of stability. In that case, the stability of growing dust grains depends mainly on the local scale height (as previously discussed by Brauer et al. 2008, see their Eq. (43)): smaller scale heights tend to favour dust pileup. For instance, scale heights smaller than 10 AU at a reference radius of 100 AU result in a pileup for grains starting at 100 AU. Λ(r) increases as r decreases and the conditions for pileup become less stringent as smaller radii, see for instance the righthand panel of Fig. 2 at r = 10 AU, where most of the disc configurations will be favourable for dust grain pileup.
The previously described criteria only tell us whether pileup will occur. For this pileup to help the formation of planetary cores, it has to occur at a distance that is sufficiently far from the star that the dust grains do not reach the very central regions where dust sublimation or the gas accretion flow will remove the dust particles before they can be incorporated in planetesimals. The pileup of growing dust grains, when it happens, is extremely efficient (the migration time goes as an exponential of a power law with respect to the grain’s position), and the grains literally stop their migration. The stopping radius of the grains depends on the relative efficiency of the growth and migrations processes and the value of Λ. Figure 3 shows the decoupling radius R_{d} as a fraction of the initial radius and as a function of the Λ parameter for the various values of p and q obtained in our grid of radiative transfer models. We estimated R_{d} numerically, by solving simultaneously the equations for grains migration and growth using the (YL; Vnew) model of L14. We calculated R_{d} as the grain’s radius in the asymptotic limit. We used t = 10^{4} years here and checked that the asymptotic limit was reached. The decoupling radius is essentially not sensitive to the value of the surface density exponent p. However, a noticeable dispersion is obtained when varying the temperature exponent q. The decoupling of the grains is more efficient for discs with smooth thermal profiles. For Λ = 1, the decoupling radius can be located at up to 6% of the initial distance to the central star for low values of q (0.15). For larger q, the corresponding decoupling radius decreases and reach values lower than 1% for q = 0.65. For Λ = 5, the decoupling occurs at values ranging between 18 and 30% of the initial radius. This process was also observed in numerical simulations as presented in LGFM08, explaining why dust grains can be found to overcome the radialdrift barrier, depending on the physical conditions in the disc.
Fig. 4 Initial dusttogas ratio required for growing dust grains to reach a relative decoupling radius R_{dec} = 0.01 (left) and 0.1 (right) for a scale height of 0.1 times the local radius. The ranges for p and q values corresponding to the values we have obtained in the radiative transfer calculations are delimited by the dashed lines (see Fig. 1). Dusttogas ratio between ≈0.01 and 0.03, and between ≈0.05 and 0.1 are sufficient for a grain starting at 100 AU to stop its migration at 1 and 10 AU, respectively. 

Open with DEXTER 
Figure 4 presents the dusttogas ratios that are required for the dust grains to decouple at a final radius equal to 1% (left) and 10% (right) of their initial radius. Dust grains starting a few hundred AU will stop their migration at a distance greater than a few AU as long as the dusttogas ratio is higher than the generally assumed interstellar value 0.01. A slight enhancement compared to this value, for instance ϵ_{°} = 0.03, is enough to obtain a decoupling in almost all the disc configurations we have explored. For a higher dusttogas ratio, the decoupling will occur at even greater distances. With ϵ_{°} between 0.05 and 0.1, most of the disc configurations will show dust grains decoupling at a radius larger than 10% of their initial radius. The required ϵ_{°} for dust grains to pile up at a given radius is smaller for discs with flatter surface density profiles and shallower midplane temperature gradients. For a given surface density, there are large variations in the required ϵ_{°} depending on the temperature structure of the disc. For instance, for p = 0.75, the required values for ϵ_{°} to obtain dust pileup range from 7.5 × 10^{3} to 0.032 for a relative decoupling radius of 0.01 and from 0.043 to 0.085 for a relative decoupling radius of 0.1. Taking the details of the temperature structure and radiative transfer effects into account is therefore critical in dust evolution simulations to accurately predict the outcome of the dust migration processes.
3.2. Stability of nongrowing dust grains and planetesimals relative to radial migration
When dust grains reach centimetric size, they start bouncing instead of sticking during collisions. This reduces the growth efficiency dramatically (Güttler et al. 2010; Zsom et al. 2010; Windmark et al. 2012a), though stochastic processes may assist the formation of larger bodies (Okuzumi et al. 2011; Windmark et al. 2012b; but see also their Corrigendum Windmark et al. 2012c; and Garaud et al. 2013). When the growth becomes negligible, the migration outcome is given by the criterion derived in LGM12: −p_{gas} + q + 1/2 < 0 (displayed in Fig. 1). We note that all discs with shallow initial profiles (p_{dust} ≲ 0.7) will lose a large fraction of their grains by accretion onto the central star. In the absence of additional mechanisms to slow down the migration, and if these grains do not cross regions where the conditions allow them to grow again, these discs will be efficiently depleted in millimetresize grains. For profiles with p_{dust} > 1.2, grains pile up and the discs efficiently retain their grain population in the inner disc. For intermediate surface density profiles (0.7 < p_{dust} < 1.2), either dust migration or dust pileup can happen, depending on the details of the disc structure and grain properties. Interestingly, this includes the standard case p_{gas} = 1, which is the limiting profile of a viscous evolution after several viscous timescales (for a constant α_{SS}). This range of values for p_{dust} is close to the observed surface density in the millimetre regime. For instance, Andrews et al. (2011) found a narrow distribution of surface density gradients for a sample of discs in Ophiuchus: ⟨p_{dust}⟩ ≃ 0.9 ± 0.2, which does not seem to depend on the millimetre flux (in their Table 4, only one disc in 12 is found to have p_{dust} < 0.7), suggesting that most of the discs are located at the interface between the stable and unstable regimes with respect to the migration of their millimetre/centimetre grains.
The highest value of q we obtain remains less than the critical value of 2/3 predicted for planetesimals to be efficiently accreted onto the central star. The thermal profile of the outer disc is therefore favourable to retaining potential planetesimals, regardless of the grain distribution.
For nongrowing grains or planetesimals, the pileup process is smoother than for the growing grain case. There is no actual stopping point since the pileup occurs progressively, slowing down the particles more and more as they reach the disc’s inner region. (The migration time goes as a power law with respect to the grain’s position.) Quantitatively, particles reduce their migration efficiency compared to the case without pileup when they reach a radius that is roughly 10% of their initial distance to the central star (Fig. 3 of LGM12). Then, it takes an order of magnitude more time for the grain to decrease its radial position by another factor 10.
For instance, in a disc with p = 0.4 and q = 0.4 (no pileup) and a disc mass of 0.1 M_{⊙} and gastodust ratio of 100, it takes ≃3.8 × 10^{5} years for a 1 mm grain initially located at 100 AU to reach R = 10 AU, a mere additional ≃3 × 10^{4} years to reach R = 1 AU and an even fewer additional ≃4 × 10^{3} years to reach R = 0.1 AU (i.e. ≃4.1 × 10^{5} years in total). For comparison, in the same disc but with p = 1.6 and q = 0.4 (pileup), it takes ≃3.9 × 10^{5} years for a 1 mm grain initially located at 100 AU to reach R = 10 AU (i.e. roughly the same time as without the pileup) but more additional ≃1.8 × 10^{6} years to reach R = 1 AU and an even more additional ≃9.1 × 10^{6} years to reach R = 0.1 AU (i.e. ≃1.1 × 10^{7} years in total). Figures 6 to 10 of LGM12 shows a parameter study for the time required for a grain to reach the inner regions of typical classical TTauri discs. The global behaviour is not modified if the disc mass is changed and only the respective times are scaled. For a lower disc mass, the corresponding times are reduced in the case of pileup, while they are increased when there is no pileup.
Fig. 5 Same as Fig. 1 but considering only models with a_{max} = 10 mm. Models are separated according to their degree of settling: no settling (black), low settling (α_{SS} ≥ 10^{2}, red), or strong settling (α_{SS} < 10^{2}, blue). 

Open with DEXTER 
3.3. Effect of grain size distribution and vertical settling
In the case of optically thin dust with an opacity law that can be described by a power law κ_{abs} ∝ λ^{−βdust}, the temperature is determined by (2)resulting in a temperature profile . As the effective grain size grows, β_{dust} varies from 2 for small ISMlike grains to 0 for millimetre or larger grains, the resulting q values range from 1/3 to 1/2. Taking the details of radiative transfer into account, e.g. high optical depth and shadowing effects, we obtain values of q ranging from 0.10 to almost 0.65 with the same trend of higher q values for larger grains, on average. (Figure 1 illustrates this effect with 2 subsets of the disc models: a_{max} ≤ 10 μm and and a_{max} ≥ 100 μm.) As a consequence, discs where dust grain have grown will have a steeper temperature profile and will be more prone to losing their nongrowing grain population by accretion onto the central object. In contrast, a possible regeneration of micrometresized grains (Brauer et al. 2008; Güttler et al. 2010; Zsom et al. 2010) will make the temperature profile shallow and stabilise the outcome of the millimetre grains.
Fig. 6 Stability criteria Λ(r) + q − 1 > 0 for growing dust grains when already grown dust grains (left) and vertical settling (right) are present. Colour schemes are the same as for Fig. 1 (left): black: all models; red: a_{max} = 1 μm–10 μm, blue; a_{max} = 100 μm–10 mm and for Fig. 5 (right): no settling (black); α_{SS} ≥ 10^{2} (low settling, red); or α_{SS} < 10^{2} (strong settling, blue). H(100 AU) is defined here as where T(r) is the dust midplane temperature, r = 100 AU, M_{star} = 1 M_{⊙}, and μ = 2.3 is the mean molecular weight. 

Open with DEXTER 
The temperature gradient becomes steeper (larger q) as the degree of dust settling increases (smaller α_{SS}, Fig. 5). Dust settling lowers the altitude of the τ = 1 surface, thus reducing the amount of light intercepted at large radii and steepening the radial temperature profile. This is critical for discs with a surface density slope close to –1.0, since dust settling may in that case be enough for −p + q + 1/2 to become positive resulting in a strongly increased migration rate for mm/cm sized grains.
For growing dust grains, two opposite effects have to be considered. Grain growth and/or dust settling will steepen the temperature profile, but they will also reduce the midplane temperature in the outer disc. As illustrated in Fig. 6 (where results are presented using the local hydrostatic scale height), the latter effect dominates: Λ + q − 1 increases in the presence of already grown dust grains and/or settling. Both the presence of grain growth and vertical settling tend to favour pileup of the grains. In particular, if a_{max} > 100 μm and α_{SS} ≤ 10^{2} (blue models in the right panel of Fig. 6), pileup of dust grains will almost always occur, even with ϵ_{°} = 0.01. (87% of this subset of models have Λ(100 AU) > 1 − q.)
If the dust opacities are significantly lowered in the outer disc due to grain growth or dust settling, the outer disc structured can be strongly modified, changing from a flared geometry to selfshadowed one (e.g. Dullemond & Dominik 2004). Because the outer parts of the disc do not intercept much stellar light directly, the temperature drops significantly, resulting in reduced local sound speed and scale height. Around a T Tauri star, for instance, the midplane temperature can drop from 40 K down to 10 K at 10 AU and from 20 K down to 5 K at at 100 AU when changing the geometry from a flared to a selfshadowed geometry. This translates into a factor 2 on the gas scale height and the resulting dusttogas ratio required for dust grains to pile up is also reduced by 4. In a selfshadowed disc with a dusttogas ratio of 0.01, growing dust grains will thus stop their migration at a radius larger than 10% of their initial radius. Such discs represent excellent sites for forming giant planet cores even in the absence of dust trapping mechanisms. We note, however, that discs are immersed in molecular clouds, whose temperatures range from ≈7 K to 15 K (e.g. Schnee et al. 2005). This is likely to set a lower limit to the temperature reached in the coldest parts of the disc’s midplane and somehow counterbalance the effects of selfshadowing at very large distances from the star.
Fig. 7 Local stability criteria for dust grains (left panel) and planetesimals (central panel) and growing dust grains (central and right panel) as a function of the radius in the disc, assuming various surface density slopes: p_{dust} = 0.5 (blue dashed line), p_{dust} = 1 (red full line), p_{dust} = 1.5 (green dashdotted line). q_{loc} is defined as −dlog (T(r)) /dlog (r). 

Open with DEXTER 
3.4. The radialdrift barrier in the central few astronomical units
The midplane temperature can only be approximated by a power law in the outer regions (where the heating is due to the reprocessing of the stellar light by the disc’s surface). However, a local temperature slope can be computed at any point to derive a local stability criteria, describing the proneness of grains to be accreted across a specific region of the disc. The stability criteria were only formally established for a power law temperature distribution, but their evolution as a function of radius, when the temperature structure differs from a power law, can be used to estimate if the radial drift will be enhanced or slowed down.
Figure 7 shows the stability criteria for discs of different surface density slopes. In the outer regions, q_{loc} remains almost constant (since the temperature structure is similar to a power law), but with a shallow increase with decreasing radius. As a consequence, the conditions for dust stability relative to migration remain mostly unchanged over a wide range of radii (>1 – a few 100 AU). The p ≈ 1 case is of particular interest since −p + q + 1/2 can become positive at a relatively large radius (≈10 AU). Inside this radius, migration of mm/cm sized grains becomes more efficient, and the disc may display a lack of continuum emission at millimetre wavelengths in the central regions. In contrast, the pileup of still growing dust grains becomes more efficient at smaller radius, and the stability criteria starts to be validated for radii smaller than 15 to 30 AU depending on the surface density profile.
The temperature profiles in the very inner disc become extremely steep (q ≫ 1) due to direct heating from the stellar photosphere. As a consequence, all dust grains in the Epstein regime, growing or nongrowing, as well as planetesimals in the Stokes regime, will experience very efficient migration, independently of the value of p. Accretion heating (not included in our modelling) can affect the midplane temperature but only within a few AU of the star (e.g. D’Alessio et al. 1998) where it can result in an enhanced migration rate for the dust grains and planetesimals.
3.5. Long term evolution of the migration process
As the dust migrates inward the dust surface density starts to depart from the initial gas distribution. Neglecting grain growth, Youdin & Shu (2002) show that the drift velocity of grains in the Epstein regime is given by v_{drift}(r) ∝ r^{pgas−q + 1/2}. This differential drift velocity results in dust concentration in the central regions, and steeper dust surface density profiles if −p_{gas} + q + 1/2 < 0 (and flatter profiles if − p_{gas} + q + 1/2 > 0). More complex profiles are obtained if grain growth is taken into account (e.g. LGFM08).
The stability criteria for nongrowing grains and planetesimals are only marginally affected though, because (i) the temperature gradient on large scales depends moderately on p_{dust} (all other parameters being fixed, see Fig. 7 right panel for instance) and (ii) this gradient is mainly set by the distribution of small grains (<1 μm), which intercept most of the stellar radiation in the disc surface and which remain coupled to the gas distribution. We assume here that the small dust population is mostly uninfluenced by the growth process. In the outer disc, where the “plateau of fast migration” occurs, the maximum expected grain size is roughly equal to the optimal size of migration (3)where ρ_{g} and c_{s} are the gas density and sound speed, ρ the intrinsic dust density and Ω_{K} the Keplerian frequency. The parameter s_{opt} scales like the disc’s surface density; i.e., s_{opt} ∝ r^{− pgas}. We computed a second version of the grid presented in Table 1 where the maximum grain size at each radius is set s_{opt}, whereas the mass and density slope for each grain size bin is kept unchanged. As expected and as shown in Fig. B.1,the resulting temperature distribution is indeed hardly modified.
Assuming that the gas surface density remains more or less constant (the viscous timescale of the disc being much longer than the migration timescale), the ability of nongrowing dust grains and planetesimals to overcome the radialdrift barrier will thus also only slightly evolve with time.
In the case of growing dust grains, if the condition for pileup are met, the pile up efficiency will increase with time as migrating dust grains will enter regions where previous dust grains have already piledup and where the local dusttogas ratio is higher (L14). Interestingly, LGFM08 found that dust grains break through the migration barrier at a few 10 AU with p = 3/2, q = 3/4, and ϵ_{°} = 0.01. These simulations were performed with H/r = 0.05, meaning that dust grains should pile up at large distances for ϵ_{°} ≳ 0.04. (This number is derived from the right panel of Fig. 4, scaled by 1/4 to account for the different H/r.) A proper treatment of hydrodynamics effects results in a pileup with a smaller dusttogas ratio illustrating that the criteria we describe are strict conditions but that pileup can occur for less stringent conditions. Dust grains that pile up in the inner disc regions increase the local dusttogas ratio and trigger the pileup of the dust arising from the outer disc (see L14). Moreover, migration becomes less efficient as the dusttogas ratios reaches values close to unity (Nakagawa et al. 1986). These two effects favour the pileup of the grains and were not accounted for in the L14 model. Conversely, pileup will always occur if the criteria (and the assumptions they are based on) are met.
4. Conclusions
Based on analytical calculations of the dust evolution in discs (we highlight the limitations in Sect. 2.3), we have investigated the effects of a complete treatment of radiative transfer on the stability of dust grains and planetesimals relative to radial migration within protoplanetary discs. We find that for a significant fraction of the discs we considered, the dust particles pile up in the disc, potentially providing material for the formation of planetary cores.
Excluding the first few central tenths of AU, all the disc configurations we explored lead to favourable temperature profiles (q < 2/3) for the discs to retain their planetesimals.
The necessary criterion q < 1 for discs to retain their growing grains is also always fulfilled. In most of the disc configurations, the conditions are met for dust grains to decouple at a radius larger than 1% if the initial dusttogas ratio is higher than 0.03. In the case of a flat surface density, lower initial dusttogas ratios, around 0.01, are enough. In those cases, dust grains starting at a few hundred AU will pile up and stop migrating in the regions where telluric planets are formed. If the initial dusttogas ratio is higher, between 0.05 and 0.1, dust grains pile up at a radius larger than 10% of their initial radius and could provide the material for forming the cores of massive planets.
Three classes of discs can be distinguished with respect to the outcome of the radial motion of nongrowing grains. Discs with steep density profiles (p_{dust} > 1.2) will retain their grains, whereas dust migration will be very efficient for flat surface densities (p_{dust} < 0.7). For intermediate cases (0.7 < p_{dust} < 1.2), significant migration of dust grains can occur, but its outcome needs to be studied via casebycase detailed simulations coupling hydrodynamics and radiative transfer. Interestingly, this includes the case p_{dust} ≃ 1, which corresponds to the peak of the distribution of the observed disc profiles.
The presence of large grains and vertical settling steepens the temperature profile and tends to enhance the radial migration of the nongrowing grains, favouring accretion onto the central object. For growing dust grains, however, the reduced temperature in presence of large grains and/or settling will slow down the migration and result in pileup for almost all the disc configurations.
If the settling and/or grain growth increases, the disc structure can become selfshadowed, resulting in much lower gas temperature and sound speed. In this case, the conditions become much more favorable for dust grains piling up and dusttogas ratios over around 0.01 are enough to stop the migration of the grains coming from 100 AU at distances over a few AU.
In the central regions of the disc where the stellar radiation can penetrate directly or accretion heating can contribute, the temperature profile is so steep (q ≫ 1) that both the dust grains, whether they are still growing or not, and planetesimals are rapidly accreted in all cases. This will also be the case if accretion heating becomes significant and results in a steep temperature profile in these central regions.
Understanding the detailed outcome of the radial migration in protoplanetary discs requires numerical simulations to catch the complexity of all the processes at play (LGFM08, Birnstiel et al. 2009, 2010). Ideally, these simulations should be coupling radiative transfer and dynamics, since the conditions for the migration or pileup of the dust grains are strongly affected by the thermal structure of the disc. In this simple study, we performed an exhaustive exploration of the parameter space and, based on analytical criteria, we give the main results of the migration outcome. These results can provide valuable help in tailoring the input parameters of the dust evolution codes used to understand the now growing evidence of dust radial segregation in discs.
Σ(r) and T(r) as observed in discs strongly differ from the socalled minimal mass solar nebula disc models (e.g. Weidenschilling 1977; Hayashi 1981) used in most of theoretical studies.
Acknowledgments
We thank F. Ménard, C. Dougados, J.C. Augereau, J.F. Gonzalez, S. T. Maddison, D. J. Price and J. B. Kajtar for useful discussions. C. Pinte acknowledges funding from the European Commission’s FP7 (contract PERG06GA2009256513) and the Agence Nationale pour la Recherche of France (contract ANR2010JCJC050401). G. Laibe is grateful to the Australian Research Council for funding (contract DP1094585) and acknowledges funding from the European Research Council for the FP7 ERC advanced grant project ECOGAL. Computations were performed at the Service Commun de Calcul Intensif de l’Observatoire de Grenoble (SCCI).
References
 Andrews, S. M., Wilner, D. J., Hughes, A. M., Qi, C., & Dullemond, C. P. 2009, ApJ, 700, 1502 [NASA ADS] [CrossRef] [Google Scholar]
 Andrews, S. M., Wilner, D. J., Hughes, A. M., Qi, C., & Dullemond, C. P. 2011, ApJ, 723, 1241 [NASA ADS] [CrossRef] [Google Scholar]
 Beckwith, S. V. W., Henning, T., & Nakagawa, Y. 2000, Protostars and Planets IV, 533 [Google Scholar]
 Birnstiel, T., Dullemond, C. P., & Brauer, F. 2009, A&A, 503, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Birnstiel, T., Dullemond, C. P., & Brauer, F. 2010, A&A, 513, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bjorkman, J. E., & Wood, K. 2001, ApJ, 554, 615 [NASA ADS] [CrossRef] [Google Scholar]
 Blum, J., & Wurm, G. 2008, ARA&A, 46, 21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brauer, F., Dullemond, C. P., & Henning, T. 2008, A&A, 480, 859 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Chiang, E., & Youdin, A. N. 2010, Ann. Rev. Earth Planet. Sci., 38, 493 [NASA ADS] [CrossRef] [Google Scholar]
 D’Alessio, P., Canto, J., Calvet, N., & Lizano, S. 1998, ApJ, 500, 411 [NASA ADS] [CrossRef] [Google Scholar]
 D’Alessio, P., Calvet, N., Hartmann, L., FrancoHernández, R., & Servín, H. 2006, ApJ, 638, 314 [NASA ADS] [CrossRef] [Google Scholar]
 Dominik, C., Blum, J., Cuzzi, J. N., & Wurm, G. 2007, Protostars and Planets V, 783 [Google Scholar]
 Duchêne, G., McCabe, C., Ghez, A. M., & Macintosh, B. A. 2004, ApJ, 606, 969 [NASA ADS] [CrossRef] [Google Scholar]
 Dullemond, C. P. 2002, A&A, 395, 853 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Dullemond, C. P. 2013, Astron. Nachr., 334, 589 [NASA ADS] [Google Scholar]
 Dullemond, C. P., & Dominik, C. 2004, A&A, 421, 1075 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Fromang, S., & Nelson, R. P. 2009, A&A, 496, 597 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Fromang, S., & Papaloizou, J. 2006, A&A, 452, 751 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Garaud, P., Meru, F., Galvagni, M., & Olczak, C. 2013, ApJ, 764, 146 [NASA ADS] [CrossRef] [Google Scholar]
 Gibbons, P. G., Rice, W. K. M., & Mamatsashvili, G. R. 2012, MNRAS, 426, 1444 [NASA ADS] [CrossRef] [Google Scholar]
 Güttler, C., Blum, J., Zsom, A., Ormel, C. W., & Dullemond, C. P. 2010, A&A, 513, A56 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hayashi, C. 1981, Progr. Theoret. Phys. Suppl., 70, 35 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Hughes, A. M., Wilner, D. J., Qi, C., & Hogerheijde, M. R. 2008, ApJ, 678, 1119 [NASA ADS] [CrossRef] [Google Scholar]
 Isella, A., Carpenter, J. M., & Sargent, A. I. 2009, ApJ, 701, 260 [NASA ADS] [CrossRef] [Google Scholar]
 Laibe, G. 2014, MNRAS, 437, 3037 (L14) [NASA ADS] [CrossRef] [Google Scholar]
 Laibe, G., Gonzalez, J.F., Fouchet, L., & Maddison, S. T. 2008, A&A, 487, 265 (LGFM08) [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Laibe, G., Gonzalez, J.F., & Maddison, S. T. 2012, A&A, 537, A61 (LGM12) [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Laibe, G., Gonzalez, J.F., & Maddison, S. T. 2014, MNRAS, 437, 3025 [NASA ADS] [CrossRef] [Google Scholar]
 Lucy, L. B. 1999, A&A, 344, 282 [NASA ADS] [Google Scholar]
 Nakagawa, Y., Sekiya, M., & Hayashi, C. 1986, Icarus, 67, 375 (NSH86) [NASA ADS] [CrossRef] [Google Scholar]
 Okuzumi, S. 2009, ApJ, 698, 1122 [NASA ADS] [CrossRef] [Google Scholar]
 Okuzumi, S., Tanaka, H., Takeuchi, T., & Sakagami, M.A. 2011, ApJ, 731, 95 [NASA ADS] [CrossRef] [Google Scholar]
 Okuzumi, S., Tanaka, H., Kobayashi, H., & Wada, K. 2012, ApJ, 752, 106 [NASA ADS] [CrossRef] [Google Scholar]
 Pascucci, I., Wolf, S., Steinacker, J., et al. 2004, A&A, 417, 793 [NASA ADS] [CrossRef] [EDP Sciences] [MathSciNet] [Google Scholar]
 Pinilla, P., Birnstiel, T., Ricci, L., et al. 2012, A&A, 538, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pinte, C., Ménard, F., Duchêne, G., & Bastien, P. 2006, A&A, 459, 797 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pinte, C., Fouchet, L., Ménard, F., Gonzalez, J.F., & Duchêne, G. 2007, A&A, 469, 963 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pinte, C., Padgett, D. L., Ménard, F., et al. 2008, A&A, 489, 633 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pinte, C., Harries, T. J., Min, M., et al. 2009, A&A, 498, 967 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schnee, S. L., Ridge, N. A., Goodman, A. A., & Li, J. G. 2005, ApJ, 634, 442 [NASA ADS] [CrossRef] [Google Scholar]
 Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
 Stepinski, T. F., & Valageas, P. 1997, A&A, 319, 1007 [NASA ADS] [Google Scholar]
 Weidenschilling, S. J. 1977, MNRAS, 180, 57 [NASA ADS] [CrossRef] [Google Scholar]
 Windmark, F., Birnstiel, T., Güttler, C., et al. 2012a, A&A, 540, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Windmark, F., Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2012b, A&A, 544, L16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Windmark, F., Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2012c, A&A, 548, C1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Youdin, A. N., & Chiang, E. I. 2004, ApJ, 601, 1109 [NASA ADS] [CrossRef] [Google Scholar]
 Youdin, A. N., & Shu, F. H. 2002, ApJ, 580, 494 [NASA ADS] [CrossRef] [Google Scholar]
 Zsom, A., Ormel, C. W., Güttler, C., Blum, J., & Dullemond, C. P. 2010, A&A, 513, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Appendix A: Validity of the criteria outside the disc midplane
In protoplanetary discs, grains migrate because the local azimuthal velocity of the gas v_{gθ} is slightly slower than the Keplerian velocity (i.e. subKeplerian rotation). Equation (B.21) of LGM12 gives the analytic expression for v_{gθ} with respect to the radial and vertical coordinates (r,z): (A.1)where is the gas pressure in the disc midplane, c_{s}(r,z) is the sound speed and f is a function defined by (A.2)The first term of the righthand side of Eq. (A.1) is the Keplerian contribution. The second term is the pressure gradient term, which depends only on r. LGM12’s criteria originate from this contribution. The third term is the baroclinic term. It contributes to the deviation to the Keplerian velocity only if z ≠ 0. However, the ratio between the baroclinic term and the pressure gradient term is of order . This implies that for particles slightly outside the disc midplane, the baroclinic term has a negligible contribution and all the criteria used above apply.
Appendix B: Models including simple dust migration
Fig. B.1 Same as Fig. 1 but with migration for s>s_{opt}. Taking into account the migration only marginally affects the temperature profile and subsequent migration of the remaining grains. The 2 dotted lines show the envelope of the models without migration (displayed in Fig. 1) for comparison. 

Open with DEXTER 
The results of the grid of models where we assume that a_{max} is a function of the radius and is set to be equal to the optimum size of migration: (B.1)are presented in Fig. B.1.
All Tables
Range of values explored for the various parameters in the grid of models computed with MCFOST.
All Figures
Fig. 1 Temperature exponent q as a function of the surface density exponent p_{dust}. Large bullets show the median values, smaller bullets the 1st and 2nd quartiles, and vertical bars the complete distribution of models; i.e., the horizontal ticks on each bar represent the minimum and maximum values reached in the whole set of models. Models with maximum grain sizes ranging from 1 μm to 10 μm are shown in red, while models with a_{max} ranging from 100 μm to 10 mm are shown in blue. The different models are slightly shifted along the p axis for clarity. The plane is divided by the stability criteria for dust grains: q = 1 with growth, −p_{gas} + q + 1/2 = 0 (assuming p_{gas} = p_{dust}, black dashed line) without growth and q = 2 /3 for planetesimals (black dotdashed lines). The thick gray line indicates the location of steadystate, constant α viscous disc models p_{gas} = 3/2 − q. 

Open with DEXTER  
In the text 
Fig. 2 Stability criteria Λ(r) + q − 1 > 0 for the growing grains as a function of the surface density exponent for r = 10 AU (left panel) and r = 100 AU (right panel) and an initial dusttogas ratio ϵ_{0} = 0.01 in both cases. Models with scale height of 5, 10, and 15 AU at a radius of 100 AU are shown in blue, black, and red respectively. The vertical axes are different in both panels. The second part of the stability criteria for growing grains: q < 1 is always satisfied in the disc models we have explored, see Fig. 1. 

Open with DEXTER  
In the text 
Fig. 3 Relative decoupling radius for growing dust grains (i.e. decoupling radius divided by initial radius of the grains in the disc) as a function of Λ for various values of p (left) and q (right). The red thick solid line shows a typical Classical T Tauri Star disc with p = 1, q = 0.5. 

Open with DEXTER  
In the text 
Fig. 4 Initial dusttogas ratio required for growing dust grains to reach a relative decoupling radius R_{dec} = 0.01 (left) and 0.1 (right) for a scale height of 0.1 times the local radius. The ranges for p and q values corresponding to the values we have obtained in the radiative transfer calculations are delimited by the dashed lines (see Fig. 1). Dusttogas ratio between ≈0.01 and 0.03, and between ≈0.05 and 0.1 are sufficient for a grain starting at 100 AU to stop its migration at 1 and 10 AU, respectively. 

Open with DEXTER  
In the text 
Fig. 5 Same as Fig. 1 but considering only models with a_{max} = 10 mm. Models are separated according to their degree of settling: no settling (black), low settling (α_{SS} ≥ 10^{2}, red), or strong settling (α_{SS} < 10^{2}, blue). 

Open with DEXTER  
In the text 
Fig. 6 Stability criteria Λ(r) + q − 1 > 0 for growing dust grains when already grown dust grains (left) and vertical settling (right) are present. Colour schemes are the same as for Fig. 1 (left): black: all models; red: a_{max} = 1 μm–10 μm, blue; a_{max} = 100 μm–10 mm and for Fig. 5 (right): no settling (black); α_{SS} ≥ 10^{2} (low settling, red); or α_{SS} < 10^{2} (strong settling, blue). H(100 AU) is defined here as where T(r) is the dust midplane temperature, r = 100 AU, M_{star} = 1 M_{⊙}, and μ = 2.3 is the mean molecular weight. 

Open with DEXTER  
In the text 
Fig. 7 Local stability criteria for dust grains (left panel) and planetesimals (central panel) and growing dust grains (central and right panel) as a function of the radius in the disc, assuming various surface density slopes: p_{dust} = 0.5 (blue dashed line), p_{dust} = 1 (red full line), p_{dust} = 1.5 (green dashdotted line). q_{loc} is defined as −dlog (T(r)) /dlog (r). 

Open with DEXTER  
In the text 
Fig. B.1 Same as Fig. 1 but with migration for s>s_{opt}. Taking into account the migration only marginally affects the temperature profile and subsequent migration of the remaining grains. The 2 dotted lines show the envelope of the models without migration (displayed in Fig. 1) for comparison. 

Open with DEXTER  
In the text 