Trapping dust particles in the outer regions of protoplanetary disks
P. Pinilla, T. Birnstiel, L. Ricci, C. P. Dullemond, A. L. Uribe, L. Testi, A. Natta
AAstronomy & Astrophysics manuscript no. Pinilla c (cid:13)
ESO 2018June 18, 2018
Trapping dust particles in the outer regions of protoplanetary disks
Pinilla P. , Birnstiel T. , Ricci L. , Dullemond C. P. , Uribe A. L. Testi L. , and Natta A. Zentrum f¨ur Astronomie der Universit¨at Heidelberg, Institut f¨ur Theoretische Astrophysik, Albert-Ueberle-Str. 2, 69120 Heidelberg,Germanye-mail: [email protected] Max-Planck-Institut f¨ur Astronomie, K¨onigstuhl 17, 69117 Heidelberg, Germany University Observatory Munich, Scheinerstr. 1, D-81679 M¨unchen, Germany Excellence Cluster Universe, Boltzmannstr. 2, D-85748 Garching, Germany Department of Astronomy, California Institute of Technology, MC 249-17, Pasadena, CA 91125, USA European Southern Observatory, Karl-Schwarzschild-Strasse 2, 85748 Garching, Germany INAF - Osservatorio Astrofisico di Arcetri, Largo Fermi 5, I-50125 Firenze, ItalyReceived October 5 2011 / Accepted December 5 2011
ABSTRACT
Context.
Dust particles are observed at mm sizes in outer regions of the disk, although theoretically, radial drift does not allow dustparticles to form pebbles.
Aims.
In order to explain grain growth to mm sized particles and their retention in outer regions of protoplanetary disks, as it isobserved at sub-mm and mm wavelengths, we investigate if strong inhomogeneities in the gas density profiles can slow down excessiveradial drift and can help dust particles to grow.
Methods.
We use coagulation / fragmentation and disk-structure models, to simulate the evolution of dust in a bumpy surface densityprofile which we mimic with a sinusoidal disturbance. For di ff erent values of the amplitude and length scale of the bumps, weinvestigate the ability of this model to produce and retain large particles on million years time scales. In addition, we introduced acomparison between the pressure inhomogeneities considered in this work and the pressure profiles that come from magnetorotationalinstability. Using the Common Astronomy Software Applications ALMA simulator, we study if there are observational signatures ofthese pressure inhomogeneities that can be seen with ALMA. Results.
We present the favorable conditions to trap dust particles and the corresponding calculations predicting the spectral slopein the mm-wavelength range, to compare with current observations. Finally we present simulated images using di ff erent antennaconfigurations of ALMA at di ff erent frequencies, to show that the ring structures will be detectable at the distances of the TaurusAuriga or Ophiucus star forming regions. Key words. accretion, accretion disk – circumstellar matter –stars: premain-sequence-protoplanetary disk–planet formation.
1. Introduction
The study of planet formation is an important field in astronomywith an increasing research since the middle of the twentiethcentury, however there are still countless unanswered questions.One of these questions is the observed grain growth to mmsized particles in the outer disk regions (Beckwith & Sargent(1991); Wilner et al. (2000, 2005); Testi et al. (2001, 2003);Andrews & Williams (2005); Rodmann et al. (2006); Natta etal. (2007); Isella et al. (2009); Lommen et al. (2009); Ricciet al. (2010a, 2011); Guilloteau et al. (2011)) that suggests amechanism operating that prevents the rapid inward drift (Klahr& Henning (1997), Brauer et al. (2007), Johansen et al. (2007)).Di ff erent e ff orts are aimed to explain theoretically the growthfrom small dust particles to planetesimals, which have led to thedevelopment of di ff erent numerical models, e.g. Nakagawa etal. (1981), Dullemond & Dominik (2005), Brauer et al. (2008),Zsom & Dullemond (2008), Okuzumi (2009), Birnstiel et al.(2010a). Due to the fact that circumstellar disks exhibit a widerange of temperatures, they radiate from micron wavelengthsto millimeter wavelengths, which is why they can be observedwith infrared and radio telescopes. With the construction ofdi ff erent kinds of these telescopes, e.g. Spitzer, Herschel, SMA,EVLA or ALMA, astronomers can observe with more details the material inside accretion disks around young stars. Theparallel development of theory and observations have allowedastrophysicists to study the di ff erent stages of planet formation,making this topic one of the most active fields in astronomytoday.In the first stage of planet formation, the growth fromsub-micron sized particles to larger objects is a complex processthat contains many physical challenges. In the case of smoothdisk with a radial gas pressure profile that is monotonicallydecreasing, the dust particles drift inwards owing to the factthat the gas moves with sub-keplerian velocity due to the gaspressure gradient. Before a large object can be formed, the radialdrift causes dust pebbles to move towards the star. Moreover, thehigh relative velocities due to turbulence and radial drift causethe solid particles to reach velocities that lead to fragmentationcollisions which do not allow dust particles to form largerbodies (Weidenschilling 1977; Brauer et al. 2008; Birnstiel etal. 2010b). The combination of these two problems is called“meter-size barrier” because on timescales as short as 100 years,a one meter sized object at 1 AU moves towards the star dueto the radial drift, preventing that any larger object could beformed.The observations in the inner regions of the disk, whereplanets like Earth should be formed, are very di ffi cult because a r X i v : . [ a s t r o - ph . E P ] D ec inilla P. et al.: Trapping dust particles in the outer regions of protoplanetary disks Table 1.
Parameters of the model
Parameter Values A { .
1; 0 .
3; 0 .
5; 0 . } f { .
3; 0 .
7; 1 .
0; 3 . } α − R (cid:63) [ R (cid:12) ] 2 . T (cid:63) [ K ] 4300 M disk [ M (cid:12) ] 0 . ρ s [ g / cm ] 1 . v f [ m / s ] 10 these regions are so small on the sky that few telescopes canspatially resolve them. Also, these regions are optically thick.However, what amounts to the meter-size barrier in the innerfew AU is a “millimeter-size barrier” in the outer regions of thedisk. These outer regions ( (cid:38)
50 AU) are much easier to spatiallyresolve and are optically thin. Moreover, one can use millimeterobservations, which probe precisely the grain size range of themillimeter-size. Therefore, the study of dust growth in the outerdisk regions may teach us something about the formation ofplanets in the inner disk regions.Observations of protoplanetary disks at sub-millimeter andmm wavelengths show that the disks remain dust-rich duringseveral million years with large particles in the outer regions(Natta et al. 2007; Ricci et al. 2010a). However, it is still unclearhow to prevent the inward drift and how to explain theoreticallythat mm-sized particles are observed in the outer regionsof the disk. Di ff erent mechanisms of planetesimal formationhave been proposed to resist the rapid inwards drift like:gravitational instabilities (Youdin & Shu 2002), the presence ofzonal flows(Johansen et al. 2009, 2011; Uribe et al. 2011) ordead zones of viscously accumulated gas which form vortices(Varni`ere & Tagger 2006). With the model presented here,we want to imitate mechanisms that allow to have long-livedpressure inhomogeneities in protoplanetary disks, by artificiallyadding pressure bumps onto a smooth density profile.To confront the millimeter-size barrier, it is necessary tostop the radial drift considering a radial gas pressure profilethat is not monotonically decreasing with radius. Instead, wetake a pressure profile with local maxima adding a sinusoidalperturbation of the density profile. These perturbations influencedirectly the pressure, following a simple equation of state forthe pressure in the disk. Depending on the size of the particle,the dust grains will be nearly perfectly trapped in the pressurepeaks, because a positive pressure gradient can cause those dustparticles move outwards. On the other hand, turbulence can mixpart of the dust particles out of the bumps, so that overall theremay still be some net radial inward drift. More importantly, dustfragmentation may convert part of the large particles into micronsize dust particles, which are less easily trapped and thus driftmore readily inward.In the work of Birnstiel et al. (2010b), they compared theobserved fluxes and mm spectral indices from Taurus (Ricciet al. 2010a) and Ophiucus (Ricci et al. 2010b) star-formingregions with predicted fluxes and spectral indices at mmwavelengths. They neglected the radial drift, forcing the dustparticles to stay in the outer disk regions. They aimed to keepthe spectral index at low values, which implies that the dustparticles could acquire millimeter sizes (Beckwith & Sargent1991). However, they found over-predictions of the fluxes. Asan extension of their work, the purpose of this paper is to modelthe combination of three processes: the radial drift, the radialturbulent mixing and the dust coagulation / fragmentation cycle in a bumpy surface density profile. Our principal aim is to findout how the presence of pressure bumps can help explain theretainment of dust pebbles in the outer regions of protoplanetarydisks, while still allowing for moderate drift and thus obtaning abetter match with the observed fluxes and mm spectral indeces.In addition, we show simulated images using di ff erent antennaconfigurations of the complete stage of ALMA, to study if itis possible to detect these kind of inhomogeneities with futureALMA observations.This paper is ordered as follows: Sect. 2 will describethe coagulation / fragmentation model and the sinusoidalperturbation that we take for the initial condition of the gassurface density. Section 3 will describe the results of thesesimulations, the comparison between existing mm observationsof young forming disks and the results from our model. Wediscuss if the type of structures generated by our model canbe detectable with future ALMA observations. In Sect. 4, weexplore the relation of our model with predictions of currentsimulations of the magnetorotational instability (MRI) (Balbus& Hawley 1991). Finally, Sect. 5 will summarize our results andthe conclusions of this work.
2. Dust Evolution Model
We make use of the model presented in Birnstiel et al. (2010a)to calculate the evolution of dust surface density in a gaseousdisk, radial drift, and turbulent mixing. The dust size distributionevolves due to grain growth, cratering and fragmentation.Relative velocities due to Brownian motion, turbulence, radialand azimuthal drift as well as vertical settling are taken intoaccount.In this work, we do not consider the viscous evolution of thegas disk because the aim is to investigate how dust evolution isinfluenced by stationary perturbations of an otherwise smoothgas surface density. The e ff ects of time dependent perturbationand the evolution of the gas disk will be the subject of futurework. For a comprehensive description of the numerical codewe refer to Birnstiel et al. (2010a).In order to simulate radial pressure maxima that allow thetrapping of particles, we consider a perturbation of the gassurface density that it is taken for simplicity as a sinusoidalperturbation such that: Σ (cid:48) ( r ) = Σ ( r ) (cid:32) + A cos (cid:34) π rL ( r ) (cid:35)(cid:33) , (1)where the unperturbed gas surface density Σ ( r ) is given by theself similar solution of Lynden-Bell & Pringle (1974): Σ ( r ) = Σ (cid:32) rr c (cid:33) − γ exp − (cid:32) rr c (cid:33) − γ , (2)where r c is the characteristic radius, taken to be 60AU, and γ is the viscosity power index equal to 1, which are themedian values found from high angular resolution imaging in thesub-mm of disks in the Ophiucus star forming regions (Andrewset al. 2010). The wavelength L ( r ) of the sinusoidal perturbationdepends on vertical disk scale-height H ( r ) by a factor f as L ( r ) = f H ( r ) , (3)with H ( r ) = c s Ω − , where the isothermal sound speed c s isdefined as c s = k B T µ m p , (4) r ( AU ) Σ ( g / c m ) A=0.1A=0.3zonal flowsRWI r ( AU ) −30−20−1001020 d l n P / d l n r A=0.1A=0.3zonal flowsRWI
Fig. 1.
Comparison between: The gas surface density (left plot) taken in this work (Eq. 1) for two di ff erent values of the amplitudeand constant width (dashed and dot-dashed lines). The Rossby wave instability (Regaly et al. 2011), and the presence of zonal flowsdue to MHD instabilities (Uribe et al. 2011). Right plot shows the pressure gradient for each of the gas surface density profiles.and the Keplerian angular velocity Ω is Ω = (cid:114) GM (cid:63) r , (5)with k B being the Boltzmann constant, m p the mass of the protonand µ is the mean molecular mass, which in proton mass units istaken as µ = .
3. For an ideal gas, the pressure is defined as P ( r , z ) = c s ρ ( r , z ) , (6)where ρ ( r , z ) is the gas density, such that Σ (cid:48) ( r ) = (cid:82) ∞−∞ ρ ( r , z ) dz .With the surface density described by Eq. 1, we can havepressure bumps such that the wavelength is increasing withradius. These bumps are static, which may not be entirelyrealistic. However, these can be a good approximation oflong-lived, azimuthally extended pressure bumps, that can bee.g. the result of MHD e ff ects (see Johansen et al. 2009;Dzyurkevich et al. 2010). The influence of time-dependentpressure perturbations (e.g. Laughlin et al. 2004; Ogihara et al.2007) on the dust growth process will be the topic of futurework. Left plot of Fig. 1 (dashed lines) shows the behavior ofthe perturbed surface density for two values of the amplitudeand fixed value of the width. The right plot of Fig. 1 shows thecorresponding pressure gradient.The very fine dust particles move as the gas because they arewell coupled to the gas since the stopping time is very short.In the presence of a drag force, the stopping time is definedas the time that a particle, with a certain momentum, needs tobe aligned to the gas velocity. However, when the particles arelarge enough and they are not forced to move as the gas, theyexperience a head wind, because of the sub-Keplerian velocityof the gas and therefore they lose angular momentum and move inwards. In this case, the resulting drift velocity of the particlesis given by Weidenschilling (1977): u drift = − + St ∂ r P ρ Ω . (7)Comparing Eq. 7 with the expression for the drift velocitygiven by Birnstiel et al. (2010a) (Eq. 19), we can notice thatthe drag term from the radial movement of the gas is not takenhere since we are assuming a stationary gas surface density.The Stokes number denoted by St, describes the coupling of theparticle to the gas. The Stokes number is defined as the ratiobetween the eddy turn-over time (1 / Ω ) and the stopping time.For larger bodies, the Stokes number is much greater than one,which implies that the particles are not a ff ected by the gas drag,consequently they move on Keplerian orbits. When the particlesare very small, St (cid:28)
1, they are strongly coupled to the gas.Since the gas is orbiting at sub-Keplerian velocity because itspressure support, there is a relative velocity between the dustparticles and the gas. The Stokes number equal to unity is acritical value where the particles are still influenced by the gasdrag but they are not completely coupled to the gas, instead theyare marginally coupled, and move at speeds between Keplerianand the sub-Keplerian gas velocity.The retainment of dust particles due to the presence ofpressure bumps depend on the size of the particles. Since verysmall particles, with St (cid:28)
1, are tightly coupled to the gas, theydo not drift inwards. However, radial drift becomes importantwhen the size of the particles increases and it is strongest whenSt = λ mfp and the particle size,denoted by a , satisfies that λ mfp / a ≥ /
9, the Stokes number isgiven by:St = a ρ s Σ g π , (8) where ρ s is the solid density of the dust grains, that is taken tobe constant (see Table 1). In this case particles are small enoughto be in this regime. Parameterizing the radial variation of thesound speed via c s ∝ r − q / , (9)where for a typical disk, the temperature is assumed to bea power law such that T ∝ r − q , which is an approximationof the temperature profile taken for this model. Therefore, thewavelength of the perturbation L ( r ) scales as: L ( r ) = f H ( r ) ∝ f r ( − q + / . (10)The pressure bumps have the same amplitude A andwavelength L ( r ) than the density, because the over pressure isinduced adding inhomogeneities to the gas surface density andparameterizing the temperature on the midplane by a powerlaw (Nakamoto & Nakagawa 1994). The model taken here canartificially imitate e.g. the case of zonal flows in protoplanetarydisks, where over densities create pressure bumps. Zonal flowsare formed due to MRI, which depend on the degree ofionization of the disk, i.e. on the temperature of the disk andother factors as the exposure to cosmic and stellar rays. MRIappears to be the most probable source of turbulence, and if theturbulence is not uniform, there can be excitation of long-livedpressure fluctuations in the radial direction. For instance,Johansen et al. (2009) performed shearing box simulations ofMRI turbulent disk that show large scale radial variations inMaxwell stresses of 10%. Dzyurkevich et al. (2010) presented3D global non-ideal MHD simulations including a dead zonethat induces pressure maxima of 20 − ff erent factors:First, it is important to analyze the necessary conditions tohave local outwards movement of the dust particles. Second, wecompare our assumptions with current studies on zonal flows(Johansen et al. 2009; Uribe et al. 2011). And third, we onlywork with small-enough amplitude disturbances such that thedisk has angular momentum per unit mass increasing outwards,which means it is Rayleigh stable .The Rayleigh criterion establishes that for a rotating fluidsystem, in the absence of viscosity, the arrangement of angularmomentum per unit mass is stable if and only if it increasesmonotonically outward (Chandrasekhar 1961), which implies: ∂∂ r ( rv φ ) > . (11)Since turbulence is necessary to have angular momentumtransport, instabilities may occur if a magnetic field is present,and in that cases the angular velocity decreases with radius(MRI). For a typical α − turbulent disk, the MRI time scale ismuch greater than the time that the disk needs to recover theRayleigh stability, this implies that the disk should remain quasistable at all time (Yang & Menou 2010). Any perturbation inthe gas surface density, that is Rayleigh unstable will almost -4 -3 -2 -1 a ( c m ) −7−6−5−4−3−2−101 l o g σ [ g / c m ] r ( AU ) -4 -3 -2 -1 a ( c m ) − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] Fig. 2.
Vertically integrated dust density distribution at 1 Myrfor A = . f = A = . f = Σ (cid:48) of Eq. 1 (see Eq. 8). The blueline represents the maximum size of the particles before theyreach fragmentation velocities (fragmentation barrier accordingto Eq.15).instantly be smeared out su ffi ciently to make it Rayleigh stableagain, thereby lowering its amplitude. This happens on a timescale much shorter than what MRI could ever counteract. Theangular velocity of Eq. 11 is given by: v φ = v k + r ρ ∂ P ∂ r = v k (1 − η ) (12)with η = − r Ω ρ ∂ P ∂ r . (13)The Rayleigh stability of the disk depends on the amplitudeand the width of the bumps. In this case, we want to study theinfluence of the amplitude of the perturbation, so for this analysiswe constrain the value of the wavelength of the perturbation, f equals to unity, such that it stays consistent with the valuesexpected from predictions of zonal flows models by Uribe et al.(2011) (see Fig. 1). t (Myrs) -2 -1 M d ( t ) / M d ( i n i t i a l ) A = 0.1 & r = [50, 100]AUα =10 −3 A = 0.1 & r = [50, 100]AUα =10 −3 f=0.2f=0.3f=0.7f=1.0f=3.0 t (Myrs) -2 -1 M d ( t ) / M d ( i n i t i a l ) A = 0.1 & r = [50, 100]AUα =10 −4 A = 0.1 & r = [50, 100]AUα =10 −4 f=0.2f=0.3f=0.7f=1.0f=3.0 Fig. 3.
Ratio between final and the initial dust mass between 50 AU and 100AU, at di ff erent times of evolution. Taking a constantvalue of the amplitude A = . ff erent values of the width of the perturbed density (Eq. 1). For α = − left plot and α = − right plotTaking the perturbed gas density of Eq. 1 and f = .
0, it ispossible to find the upper limit of the amplitude to satisfy the Eq.11, i.e. the condition to remain Rayleigh stable at all time. Thiscalculation lets the maximum value of the amplitude A to be atmost ∼
35% of the unperturbed density.Various investigations have been focused on the possibilityof Rayleigh instabilities when disks have sharp profiles for theradial density (see e.g. Papaloizou & Lin (1984), Li et al. (2001)and Yang & Menou (2010)). These kind of profiles can existwhen the temperature in the midplane is not su ffi cient to ionizethe gas (Gammie 1996) and as a result the turbulence parameter α is reduced. These regions are known as “dead zones” and theseare candidates to be places where planet embryos can be formed.In the boundaries of these regions, it was shown by Varni`ere &Tagger (2006) that it is possible to have a huge vortex with alocal bump in the gas surface density. Lovelace et al. (1999)demonstrated that these perturbations create an accumulationof gas that leads to the disk to be unstable to Rossby waveinstability (RWI). As a comparison of the amplitudes generatedby RWI vortices and the amplitudes of our perturbations, the leftplot of Fig. 1 also shows the azimuthally average gas surfacedensity of a large-scale anticyclonic vortex by Regaly et al.(2011). In those cases, the equilibrium of the disk is a ff ectedsuch that the disk may become Rayleigh unstable. Since we arefocused on perturbations that allow to the disks stay Rayleighstable at all evolution time, we do not consider for our perturbeddensity such kind of amplitudes.On the other hand, since the drift velocity is given by Eq. 7,to prevent the inwards drift, the value of η from Eq. 13 mustbe negative, this implies that the pressure gradient has to bepositive for some regions of the disk. Doing this calculation forthe condition that η <
0, we get that the values of the amplitude A have to be at least equal to about 10%. In right plot of Fig. 1 wesee that with an amplitude of 10% the pressure gradient barelyreaches positive values in the inner regions of the disk ( (cid:46) . (cid:46) A (cid:46) .
35, when the width of the perturbationis taken to be f =
1. Taking into account the growth-fragmentation cycle and theexistence of pressure bumps, the radial drift e ffi ciency can bereduced if the bumps have favorable values for the amplitudeand the length scale. When the particles grow by coagulation,they reach a certain size with velocities high enough to causefragmentation (fragmentation barrier). The two main sources ofrelative velocities are radial drift and turbulence. In the bumpsthe radial drift can be zero if the pressure gradient is high enoughas well as azimuthal relative velocities; but there are still relativevelocities due to the turbulence. Therefore, it is necessary to havea condition, such that the particles do not reach the thresholdwhere they fragment. The maximum turbulent relative velocitybetween particles, with St ∼
1, is given by Ormel & Cuzzi(2007), ∆ u max (cid:39) α St c s , (14)for St (cid:46) . ff by a factor of 2. Therefore, to breakthrough the mm size barrier, we must have that ∆ u max has tobe smaller that the fragmentation velocities of the particles v f .Recent collision experiments using silicates (G¨uttler et al. 2010)and numerical simulations (Zsom et al. 2010) show that thereis an intermediate regime between fragmentation and sticking,where particles should bounce. In this work, we do not take intoaccount this regime since there are still many open questionsin this field. For example, Wada et al. (2011) suggest thatthere is no bouncing regime for ice particles, which may bepresent in the outer regions of the disks (see Schegerer & Wolf2010). Laboratory experiments and theoretical work suggest thattypical values for fragmentation velocities are of the order of fewm s − for silicate dust (see e. g. Blum & Wurm 2008). Outsidethe snow-line, the presence of ices a ff ects the material properties,making the fragmentation velocities increase (Sch¨afer et al.2007; Wada et al. 2009). Since in this work we assume a radialrange from 1AU to 300AU, the fragmentation threshold velocityis taken as v f =
10m s − . All the parameters used for this modelare summarized in Table 1. r ( AU ) -2 -1 a ( c m ) St =1α =10 −3 α =10 −2 Fig. 4.
Particle size corresponding to a Stokes number of unityfor A = . f = . ff erent values of the turbulent parameter α . The dashedline corresponds to r = (cid:46)
1, taking the size at whichthe turbulent relative velocities ∆ u max are as high as thefragmentation velocity v f , we can find the maximum value ofthe grain size, which is a max (cid:39) Σ g παρ s v f c s , (15)this a max is valid only for St (cid:46) n ( r , z , a ) is the number of particlesper cubic centimeter per gram interval in particle mass, whichdepends on the grain mass, distance to the star r and height abovethe mid-plane z , such that ρ ( r , z ) = (cid:90) ∞ n ( r , z , a ) · mdm , (16)is the total dust volume density. The quantity n ( r , z , a ) canchange due to grain growth and distribution of masses viafragmentation. The vertically integrated dust surface densitydistribution per logarithmic bin defined as σ ( r , a ) = (cid:90) ∞−∞ n ( r , z , a ) · m · adz (17)and the total dust surface density is then Σ d ( r ) = (cid:90) ∞ σ ( r , a ) d ln a . (18)
3. Results
The simulations have been done with a disk of mass 0 . M (cid:12) ,with a surface density described by Eq. 1 from 1 . α is taken to be 10 − , unless other value is specified.Figure 2 shows the vertically integrated dust density distribution,taking into account: coagulation, radial mixing, radial drift andfragmentation, after 1 Myr of the evolution of the protoplanetarydisk. The solid white line shows the particle size correspondingto a Stokes number of unity. From Eq. 8 we can see that whenSt =
1, the particle size a is proportional to the gas surface density Σ (cid:48) , then the solid line reflects the shape of the surface density.The blue line of Fig. 2 represents the fragmentation barrier,which illustrates the maximum size of the particles before theyreach velocities higher than the fragmentation velocity (see Eq.15). Hence, particles above the fragmentation barrier shouldfragment down to smaller particles, which again contribute tothe growth process.Both plots of Fig. 2 have the same amplitude of thesinusoidal perturbation A = .
1. The factor f which describesthe width of the perturbation, is taken to be f = f = A = . ff ected by radial drift and turbulence such that the particlesdo not grow over mm size in the outer regions.Second, taking a greater value of the factor f , at the sameamplitude, implies that the retention of particles is even weaker.This is because with a wider perturbation is it harder to havepositive pressure gradient. It is expected that for a smallervalue of f , the pressure gradient is higher since the surfaceprofile should be steeper and therefore a dust trapping muchmore e ffi cient. However, the di ff usion timescale τ ν dependsquadratically on the length scale (cid:96) . Therefore, when the wigglesof the perturbation are taken with a smaller wavelength,the di ff usion times become much shorter, implying that theturbulence mixes the dust particles out of the bumps faster, evenwhen the pressure gradient is higher for narrow wiggles. Moreprecisely, τ ν ∝ (cid:96) ν − where the viscosity is defined as ν = α c s h .For this reason, we notice in Fig. 3-left plot that the trappingis more e ff ective taking values of the width less than one, butonly by a very small margin. As a result, the ratio betweenthe final and the initial mass of the dust for r ∈ [50 , . A = ff ective when thewavelength of the perturbation is taken shorter.Only when the di ff usion timescales become larger or equalto the drift timescales for a given pressure profile, a turbulenceparameter and a Stokes number, the dust particles can be retainedinside the bumps and therefore they can grow. From Eq. 7 we candeduce that the drift time scales as τ drift ∝ (cid:96) ( ∂ (cid:96) P ) − , where insidethe bumps (cid:96) is the width of the perturbations (which dependsdirectly on f ). As a consequence, after an equilibrium state isreached, drift and di ff usion timescales are both proportional tothe square of the width. Accordingly, for a given value of f , thee ff ect of turbulent mixing and radial drift o ff set.We can notice in Fig. 3 that for f = { } there is a smalle ff ect in the retention of particles for α = − (left plot) and animportant e ff ect for α = − (right plot). This implies that forthese values of f and A , the pressure gradient becomes positiveenough in the outer regions ( r ∈ [50 , A and f are the same for both cases, the pressuregradient is exactly the same. However, for α = − , the small -4 -3 -2 -1 a ( c m ) t = 0.5 Myr − 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 t = 0.5 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 1.0 Myr -4 -3 -2 -1 t = 1.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 3.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 t = 3.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] r ( AU ) -4 -3 -2 -1 a ( c m ) t = 5.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] r ( AU ) -4 -3 -2 -1 t = 5.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 0.5 Myr -4 -3 -2 -1 t = 0.5 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 1.0 Myr -4 -3 -2 -1 t = 1.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 3.0 Myr − 5− 4− 2− 101 10 -4 -3 -2 -1 t = 3.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] r ( AU ) -4 -3 -2 -1 a ( c m ) t = 5.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] r ( AU ) -4 -3 -2 -1 t = 5.0 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 0.5 Myr -4 -3 -2 -1 t = 0.5 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 1.0 Myr -4 -3 -2 -1 t = 1.0 Myr -4 -3 -2 -1 a ( c m ) t = 3.0 Myr -4 -3 -2 -1 t = 3.0 Myr r ( AU ) -4 -3 -2 -1 a ( c m ) t = 5.0 Myr r ( AU ) -4 -3 -2 -1 t = 5.0 Myr -4 -3 -2 -1 a ( c m ) t = 0.5 Myr -4 -3 -2 -1 t = 0.5 Myr − 7− 6− 5− 4− 3− 2− 101 l og σ [ g / c m ] -4 -3 -2 -1 a ( c m ) t = 1.0 Myr -4 -3 -2 -1 t = 1.0 Myr -4 -3 -2 -1 a ( c m ) t = 3.0 Myr -4 -3 -2 -1 t = 3.0 Myr r ( AU ) -4 -3 -2 -1 a ( c m ) t = 5.0 Myr r ( AU ) -4 -3 -2 -1 t = 5.0 Myr Fig. 5.
Vertically integrated dust density distribution with fixed value of length scale as f =
1, for A = . A = . ff erent times 0 . . . . e ffi ciency that becomes visible reducing f , vanishes after twoMyr, because turbulent mixing and radial drift cancel each other.Nevertheless, when α is reduced one order of magnitude(right panel of Fig. 3), the di ff usion timescales are now longer.In this case, we have that the drift timescales are shorter thanthe di ff usion timescales, hence the ratio between the final andthe initial dust mass increase in average for each f . When f issmall enough to have positive pressure gradient ( f = { } ),outward drift wins over turbulent mixing, and as a result there isa visible e ff ectiveness in the trapping of particles. However, thereis almost no di ff erence between f = f = ff ect between radial drift and turbulence when f varies does not happen when the amplitude of the perturbationchanges.We fixed the value of the width of the perturbation to unity,because this value is consistent with current model predictionsof zonal flows. The comparison between our assumption of thedensity inhomogeneities and the work from Uribe et al. (2011)is discussed in Sect. 4. In addition, for longer values of thewidth, we should have higher values of the amplitudes in orderto have a positive pressure gradient. In that case, however thedisk becomes easily Rayleigh unstable when the amplitude inincreased. These are the reasons why we fix the value of thewidth to unity and not higher.Simulations of MRI-active disks suggest that the typicalvalues for the turbulence parameter α are in the range of 10 − − − (Johansen & Klahr 2005; Dzyurkevich et al. 2010). In thiswork, we focus on the results for α = − , because with a largervalue of the turbulence the viscous time scales become shortercompared with the growth time scales of the dust, making theparticles mix out of the bumps and then drift radially inwardsbefore any mm-sizes are reached. In addition, if α is taken oneorder of magnitude higher, the fragmentation barrier is lowerby about one order of magnitude in grain size, implying thatparticles do not grow over mm-sizes in the outer regions of thedisk (see Eq. 14). Fig. 4 shows the location of the fragmentationbarrier for the case of A = . f = . ff erentvalues of α . We can notice that the maximum value of the grainsize for the case of α = − is of the order of few mm in theouter regions of the disk r > α = − thegrains even reach cm-sizes.Figure 5 compares the surface density distribution for twodi ff erent values of the amplitude of the perturbation A = . A = .
3, at di ff erent times of evolution. Taking A = . A = . r (cid:38) ff erent times of evolution. It is important to noticethat the line of the fragmentation barrier (15) is calculatedtaking into account only turbulent relative velocities, since radialand azimuthal turbulence relative velocities were assumingzero at the peaks of the bumps. Particles with St > ff ected by gasturbulence, then relative velocities due to turbulence are lower,which implies that they can grow over the fragmentation barrier.Moreover, we can see in the right plot of Fig. 1, that the pressuregradient for A = . r (cid:38) ffi cacy of radial drift. Therefore, the total relative velocitiesfor r (cid:38) >
1) to grow over the fragmentation barrier. Only the particleswith St < r (cid:46) t (Myrs) -2 -1 M d ( t ) / M d ( i n i t i a l ) f = 1.0 & r = [50, 100]AUf = 1.0 & r = [50, 100]AU A=0.7A=0.5A=0.3A=0.1A=0.0
Fig. 6.
Ratio between final and the initial dust mass between50 AU and 100 AU, at di ff erent times of evolution. Taking aconstant value of the width f = . ff erent values of theamplitude of the perturbed density (Eq. 1).particles continuously grow to mm-size particles by coagulationbecause collision velocities due to turbulence are lower that thetaken fragmentation velocity v f .As we mentioned before, the e ffi ciency of the dust trappingwill depend on the amplitude of the pressure bumps. It isexpected that for higher amplitude there is more trapping ofparticles, since the pressure gradient is also higher and positive(see right plot of Fig. 1). Taking the perturbed density of Eq. 1,we can see in Fig. 6 that between 50 AU and 100 AU from thestar, the amount of dust grows significantly from A = . A = .
3. From A = . A = .
5, there is still a considerablygrowth, but the rate of growth is slower. From A = . A = . f = . . (cid:46) A (cid:46) .
35, which meansthat for those values of the amplitude, these kind perturbationscan be explained via MRI without any Rayleigh instabilitypresent at any evolution time. The amplitude of A = . A = . ff erent times of the simulation, for two values for theperturbation amplitude and wavelength without the gas inwardmotion. For A = . f = . ff ective. Consequently, due toturbulence the dust particles collide, fragment and become evensmaller, so they mix and the retention of those small particles,with St <
1, becomes more di ffi cult. Thus the radial drift is notreduced and particles move inwards. Therefore particles with -5 -4 -3 -2 -1 d u s t - g a s r a t i o A=0.1 & f=1.0 -5 -4 -3 -2 -1 A=0.1 & f=0.7 r (AU) -5 -4 -3 -2 -1 d u s t - g a s r a t i o A=0.3 & f=1.0 r (AU) -5 -4 -3 -2 -1 A=0.3 & f=0.7
Fig. 7.
Dust to gas mass ratio at the disk for di ff erent times of disk evolution and the parameters summarizing in Table 1: A = . f = A = . f = A = . f = . A = . f = . (cid:46) (cid:28)
1, and are well coupled to the gas. Hencethe dust-gas ratio initially decreases quickly and then becomesalmost constant with time, which implies that after several Myronly the very small dust particles remain. Top-right plot of Fig. 7shows that taking the same amplitude but a smaller wavelength,the dust to gas radio has the same behavior. This confirms thatwhen α turbulence is constant, a decrease of the wavelengthimplies shorter di ff usion time scales. Therefore the trapping isnot more e ff ective even if the pressure gradient is higher fornarrow bumps.Conversely, due to the strong over-pressures at A = . r < ∼ − to ∼ − . This oscillating behavior, even after 5 Myrof dust evolution, is possible thanks to the fact that the particlesare retained in the bumps and grow enough to make the dust-gasratio higher inside the bumps. Only around ∼ r < f = . f = . ff ect over the dust-to gas ratio as wasexplained before. In this Section we compare the models predictions of thedisk fluxes at millimeter wavelenghts with observational dataobtained for young disks in Class II Young Stellar Objects(YSOs).To do this, we calculate the time-dependent flux for thedisk models described in Sect. 2. For the dust emissivity, weadopted the same dust model as in Ricci et al. (2010a,b, 2011)and Birnstiel et al. (2010b), i.e. spherical composite porousgrains made of silicates, carbonaceous materials and waterice, with relative abundances from Semenov et al. (2003). Ateach stellocentric radius in the disk, the wavelength-dependentdust emissivity is calculated considering the grain size number A =0.1 f =1 A =0.3 f =1 Fig. 8.
Comparison of the observed fluxes at mm-wavelengthsof young disks in Taurus (red dots; from Ricci et al. (2010a)and Ricci (priv comm)), Ophiucus (blue dots; from Ricci et al.(2010b)), and Orion Nebula Cluster (green dots; from Ricci etal. (2011)) star forming regions with the predictions of the diskmodels at di ff erent times of the disk evolution (star symbols).Disk ages are indicated by numbers, in Myr, above the starsymbols. The predicted ∼ ∼ f = A = . A = . ∼ / to account for the di ff erent distancesestimated for the Orion Nebula Cluster ( ∼ ∼ n ( r , z , a ) derived from the dust evolution models at thatradius, as described in Sect. 2.The opacity, in the millimeter wavelength regime, can beapproximated by a power law (Miyake & Nakagawa 1993),which means that the flux can be approximated to F ν ∝ ν α mm ,where α mm is known as the spectral index. The spectral indexgives us information about the size distribution of the dust inthe disk. Figure 8 shows the time-dependent predicted fluxesat ∼ F ) and spectral index between ∼ α − ) for a disk model with f = A = . A = . Fig. 9.
Disk image at 2 Myr and observing wavelength of 0.45mm, the amplitude of the perturbation is A = . f = { , , , } thecorresponding rms value (see Table 2).As detailed in Birnstiel et al. (2010b) the F vs α − plot reflects some of the main properties of the dust populationin the disk outer regions, which dominate the integrated flux atthese long wavelengths. In particular the 1mm-flux density isproportional to the total dust mass contained in the outer disk.The mm-spectral index is instead related to the sizes of grains:values lower than about 3 are due to emission from grains largerthan about ∼ f = A = . R (cid:38)
50 AU), and most of them areinitially retained in those regions ( (cid:46) . f = A = . ffi cient in retaining mm-sized particles in the outerdisk. The 1 mm-flux density significantly decreases because ofthe loss of dust from the outer regions, especially the mm-sizedgrains which are e ffi cient emitters at these wavelengths. Giventhat the spectral index is a proxy for the grain size, that is alsoa ff ected by radial drift: its value increases with time becauseof the gradual loss of mm-sized pebbles in the outer disk. For Fig. 10.
Comparison between the simulated images for anobserving wavelength of 1 mm and 2 Myr of evolution, usingfull antenna configuration of ALMA for two di ff erent values ofthe amplitude of the perturbation: A = . A = . { , , , } the correspondingrms value (Table 2 ).this case, the under-predicted fluxes are not consistent withobservational data for disk ages (cid:38) A = . ff erentresults. In this case the trapping of particles in the pressurebumps is e ffi cient enough to retain most of the large pebblesformed in the outer disk (see bottom plot of Fig.8). Since radialdrift is much less e ffi cient in this case, the predicted 1mm-fluxdensity is less a ff ected than in the A = . ff at a value of ∼ A = . ff radial drift, therefore they restricted particles tostay artificially in the disk outer regions. Specifically, for a diskwith the same unperturbed disk structure presented here, theyfound a larger 1mm-flux density than what we obtained in the A = . Fig. 12.
Spectral index α − of the model data (top plot) andthe spectral index taking two simulated images at 0 . . The Atacama Large Millimeter / sub-millimeter Array (ALMA)will provide an increase in sensitivity and resolution to observein more detail the structure and evolution of protoplanetarydisks. With a minimum beam diameter of ∼ ff er a resolution down to 2 AU for disks observedin Orion and sub-AU for disk in Taurus-Auriga (Cossins et al.2010). Using the Common Astronomy Software Applications(CASA) ALMA simulator (version 3.2.0), we run simulationsto produce realistic ALMA observations of our model usingALMA array of 50 antennas 12m-each.The selection of observing mode to obtain the images, waschosen to have simultaneously the most favorable values for theresolution and sensitivity that should be available with ALMA.The spatial resolution depends on the observing frequency andthe maximum baseline of the array. We do not take the largestarray because for very large baselines, the sensitivity couldbe not enough for the regions that we will need to observe.Therefore, we used di ff erent antenna arrays depending on theobserving frequency to get the best possible resolution with Table 2.
Atmospheric conditions, total flux and rms for the simulated observations at 140pc and at di ff erent observing wavelengths.The pwv value takes into account the expected conditions for ALMA. The simulated images are using the full ALMA, but theantenna configuration is chosen in order to have the best conciliation between resolution and sensitivity. Amplitude Wavelength Atmospheric conditions Total Flux rms(mm) pwv (mm) τ (Jy) (Jy)0.45 0 . .
60 6 . × − . × − A = . .
40 7 . × − . × − . .
20 2 . × − . × − . .
03 4 . × − . × − A = . .
20 1 . × − . × − Fig. 11.
Disk simulated images with parameters of Table 1, A = . f = { , , , } the corresponding rms value (see Table 2)enough sensitivity. The sensitivity depends on the number ofantennas, the bandwidth (which is taken as ∆ ν = . . . −
750 GHz). It is important tonote that the simulated images take into account the atmosphericconditions and the expected receiver noise based on technicalinformation of the antennas, but the residual noise after datacalibrations and its uncertainties are not considered. We can note(see Fig. 9-bottom plot) that with one of the full configurations(max. baseline ∼ A = .
3, then the e ff ects will be observable with ALMA.Both images of the figure have been compute with the completeantenna configuration of ALMA for an observing wavelengthof 1 mm and 2 Myr of the disk evolution. We can see that for A = . ∼ − ff erentobserving wavelength using di ff erent antenna configurationsof ALMA. The antenna configuration is chosen by CASAdepending on the expected resolution. The best image isobtained at 100 GHz and a maximum baseline of 16 km (mostextended ALMA configuration), where it is possible to detectclearly the most external ring structure and some internal ringstructures. Nevertheless, with more compact configurations atdi ff erent frequencies is still possible to detect some structuresfrom the presence of the pressure bumps which allow theformation of mm-sized particles. However, it is important to takeinto account that the simulated images of Figs. 9, 10 11 and12 are considering a perfect data calibration after observations,for long baselines and high frequencies the calibration e ff ectsbecome more important.Taking the ratio of the images at two di ff erent wavelengths,we will have the values of the spectral index α − mm , whichindicates the location of mm-sized grains. With the fullconfiguration of ALMA and a maximum value of the baselineof 12 km for both observing frequencies, some regions withlarge particles are distinguished which are the regions with lowspectral index α − (cid:46)
3, as it was explained in Sec. 3.2. InFig. 12 is the spectral index of the model data (top plot) andthe spectral index taking two simulated images at 0 . .
4. Approach to zonal flows predictions
In this paper we have so far assumed ad-hoc models of pressurebumps. But which processes may cause such long lived bumpsin protoplanetary disks? In this section we study zonal flows as apossible explanation of the origin of long lived pressure bumps.One possible cause of pressure bumps originates from MRIturbulence. Hawley et al. (1995) and Brandenburg et al. (1995)were the first attempts to simulate nonlinear evolution of MRI inaccretions disk, taking a box as representation of a small part ofthe disk. More recent simulations have been done taking higherresolution (see e.g. Johansen et al. 2009) and global setup (Flocket al. 2011).In magnetorotational instability, “zonal flows” are excitedas a result of the energy transportation from the MRI unstablemedium scales, to the largest scales, causing an inverse cascadeof magnetic energy, and creating a large-scale variation in theMaxwell stress (Johansen et al. 2009).Di ff erent three-dimensional MHD simulations show thatwith the presence of zonal flows, pressure bumps can appear ifthere are drops in magnetic pressure through the disk Johansenet al. (2009), Dzyurkevich et al. (2010) and Uribe et al. (2011).Nevertheless, in recent simulations with higher resolution byFlock et al. (2011), pressure bumps are not formed. There is stillnot a final answer about how and why these pressure bumps canor cannot be created via zonal flows.Another alternative for the origin of pressure bumps dueto MRI is the change of the degree of ionization. The diskbecomes MRI active if the degree of ionization is su ffi cient,for the magnetic field to be well couple to the gas. Variabledegrees of ionization in the disk could cause local changes in
20 30 40 50 60 70 800.900.951.001.051.10 Σ m a x / Σ m i n
20 30 40 50 60 70 80 r ( AU ) −0.020−0.015−0.010−0.0050.0000.0050.010 ( v φ − v K ) / v K Fig. 13.
Top plot: Ratio between the surface density at twodi ff erent azimuthal angles of of the disk from zonal flowssimulation by (Uribe et al. 2011). The azimuthal angles arechosen such that for a specific radius, the amplitude of thepressure bump has a maximum Σ max and a minimum Σ min .Bottom plot: Azimuthal velocity with respect to Keplerianvelocity for the azimuthal and time-averaged surface density ofthe midplane.the magnetic stress, which could induce structures in the densityand pressure.The question we now wish to answer is: are the pressurebumps caused by zonal flows of MRI-turbulence strong enoughto trap the dust like in the models of Sect. 3? To find this out, wetake the three-dimensional MHD simulations of MRI-turbulentprotoplanetary disks by Uribe et al. (2011). These models havea resolution of ( N r , N θ , N φ ) = (256 , , ff ordable.Hence, the strategy is to find first a quasi-steady state of thegas surface density from MRI evolution, in which structures inthe pressure survive the entire simulation (around 1000 innerdisk orbits). Afterwards, to do the coagulation / fragmentationsimulation of the dust in 1D taking the gas surface density for r ( AU ) -4 -3 -2 -1 a ( c m ) t = 2 Myr r ( AU ) -4 -3 -2 -1 a ( c m ) t = 4 Myr − 7− 6− 5− 4− 3− 2− 10 l og σ [ g / c m ] Fig. 14.
Vertically integrated dust density distribution after 2 Myr and 4 Myr of dust evolution; taking the azimuthal and time-averageprofile of the gas surface density in the midplane from the MHD simulations (see Uribe et al. 2011, figure 3) without a planet. Thesolid white line shows the particle size corresponding to a Stokes number of unity, which has the same shape of the density profile.a specific azimuthal angle in the midplane from the results fromMHD simulations.Top plot of Fig. 13 shows the ratio between the time averagedsurface density at two di ff erent azimuthal angles where theamplitude of a pressure bump is maximum and minimum fora specific radius. We can notice, that the variations on theazimuthal angle are very uniform, around ∼ ff used on turbulent di ff usion timescales.In the future, the continuous generation and evolution of thesestructures should be implemented alongside the dust evolution.However, for lack of a better model of this time-dependence atthis stage, we assume these structures to be static.Since these MHD simulations use a radial domain where r ∈ [1; 10], we rescale this grid logarithmically, so the gassurface density is taken from 10AU to 100AU, and also thesurface density has been scaled such that the total disk massis 0 . M (cid:12) , see Fig. 1-left plot (solid-line). Comparing thegas surface density obtained from MHD simulations with theassumed perturbed density Σ (cid:48) (Eq. 1), we see that the amplitudeof the surface density perturbation from zonal flows is around25% and comparable with the amplitude of 30% of Σ (cid:48) . Thewidths of the bumps from Uribe et al. (2011) are not uniform,but our assumption of f =
5. Conclusions
Theoretical models of dust evolution in protoplanetary disksshow that the growth from sub-micron sized particles to largerobjects is prevented basically by two phenomena: radial driftand fragmentation. Nevertheless, infrared and radio observationsshow that millimeter sized particles can survive under thosecircumstances in the outer regions of disks. Therefore, varioustheoretical e ff orts have been focused on explaining the survivalof those bodies.Taking into account strong inhomogeneities expected tobe in the gas density profile e.g. zonal flows, and usingthe coagulation / fragmentation and disk-structure models byBirnstiel et al. (2010a), we have investigated how the presence ofpressure bumps can cause the reduction of radial drift, allowingthe existence of millimeter sized grains in agreement withobservations. In this work, we assumed a sinusoidal functionfor the gas surface density to simulate pressure bumps. Theamplitude and wavelength disturbances are chosen consideringthe necessary conditions to have outward angular momentumtransport in an α -turbulent type disk, outward radial drift of dustand reasonable values compare to predictions from the recentwork of zonal flows (Uribe et al. 2011). The results presented here suggest that the presence ofpressure bumps with a width of the order of the disk scale-heightand an amplitude of 30% of the gas surface density of thedisk, provide the necessary physical conditions for the survivalof larger grains in a disk with properties summarized in Table1. Comparisons between the observed fluxes of the Taurus,Ophiucus and Orion Nebula Cluster star forming regions withthe results of the models ratify that the e ff ect of the radialdrift is reduced allowing particles to grow. Figure 8 shows howmodels with these kind of disturbances reproduce much bettermm-observations than models with full or without radial drift.In addition, we presented a comparison between the bumpydensity profile assumed in this work and 3D MHD models ofzonal flows that can cause long lived bumps in protoplanetarydisks. We showed that the pressure bumps cause by zonal flowsof (Uribe et al. 2011) are in agreement with the amplitudes andwavelengths used in this work. Therefore, taking those bumps,the survival of dust particles is possible in the outer regions aftersome Myr.The simulated images using CASA ALMA simulator(version 3.2.0) show that, with di ff erent antenna configurationof the final ALMA stage, the ring structures, due to the presenceof the pressure bumps, should be detectable. Future ALMAobservations will have an important impact for understandingthe first stages of planet formation and it will be very importantto investigate if the grain growth and retetion can be explainedwith the presence of these kind of inhomogeneities in the gasdensity profile. Acknowledgements.
We acknowledge Francesco Trotta for his help with thecode that we used in this work to derive the mm-fluxes. We would like to thankthe referee, Wladimir Lyra, for his useful suggestions. This work was supportedin part through the 3rd funding line of German Excellence Initiative. T. Birnstielacknowledges ESO O ffi ce for Science which provided funding for the visitsin Garching. L. Ricci acknowledges the PhD fellowship of the InternationalMax-Planck-Research School. A. L. Uribe acknowledges the CPU time forrunning the simulations in the Bluegene / P supercomputer and the THEO clusterat the Rechenzentrum Garching (RZG) of the Max Planck Society. Finally,L. Testi acknowledges ASI contract to the INAF-Osservatorio Astrofisico DiArcetri.
References
Andrews, S. M., & Williams, J. P. 2005, ApJ, 631, 1134Andrews, S. M., Wilner, D. J., Hughes, A. M., Qi, C., & Dullemond, C. P. 2010,ApJ, 723, 1241Balbus, S. A., & Hawley, J. F. 1991, ApJ, 376, 214Beckwith, S. V. W., & Sargent, A. I. 1991, ApJ, 381, 250Bertout, C., Robichon, N., & Arenou, F. 1999, A&A, 352, 574Birnstiel, T., Dullemond, C. P., & Brauer, F. 2010a, A&A, 513, A79Birnstiel, T., et al. 2010b, A&A, 516, L14Blum, J., & Wurm, G. 2008, ARA&A, 46, 21Brandenburg, A., Nordlund, A., Stein, R. F., & Torkelsson, U. 1995, ApJ, 446,741Brauer, F., Dullemond, C. P., Johansen, A., Henning, T., Klahr, H., & Natta, A.2007, A&A, 469, 1169Brauer, F., Dullemond, C. P., & Henning, T. 2008, A&A, 480, 859Chandrasekhar, S. 1961, International Series of Monographs on Physics, Oxford:Clarendon, 1961,Cossins, P., Lodato, G., & Testi, L. 2010, MNRAS, 407, 181Dullemond, C. P., & Dominik, C. 2005, A&A, 434, 971Dzyurkevich, N., Flock, M., Turner, N. J., Klahr, H., & Henning, T. 2010, A&A,515, A70Flock, M., Dzyurkevich, N., Klahr, H., Turner, N. J., & Henning, T. 2011, ApJ,735, 122Gammie, C. F. 1996, ApJ, 457, 355Guilloteau, S., Dutrey, A., Pi´etu, V., & Boehler, Y. 2011, A&A, 529, A105G¨uttler, C., Blum, J., Zsom, A., Ormel, C. W., & Dullemond, C. P. 2010, A&A,513, A56Johansen, A., & Klahr, H. 2005, ApJ, 634, 1353 Johansen, A., Oishi, J. S., Mac Low, M.-M., Klahr, H., Henning, T., & Youdin,A. 2007, Nature, 448, 1022Johansen, A., Youdin, A., & Klahr, H. 2009, ApJ, 697, 1269Johansen, A., Klahr, H., & Henning, T. 2011, A&A, 529, A62Klahr, H. H., & Henning, T. 1997, Icarus, 128, 213Hawley, J. F., Gammie, C. F., & Balbus, S. A. 1995, ApJ, 440, 742Isella, A., Carpenter, J. M., & Sargent, A. I. 2009, ApJ, 701, 260Laughlin, G., Steinacker, A., & Adams, F. C. 2004, ApJ, 608, 489Li, H., Colgate, S. A., Wendro ff , B., & Liska, R. 2001, ApJ, 551, 874Lommen, D., Maddison, S. T., Wright, C. M., van Dishoeck, E. F., Wilner, D. J.,& Bourke, T. L. 2009, A&A, 495, 869Lovelace, R. V. E., Li, H., Colgate, S. A., & Nelson, A. F. 1999, ApJ, 513, 805Lynden-Bell, D., & Pringle, J. E. 1974, MNRAS, 168, 603Menten, K. M., Reid, M. J., Forbrich, J., & Brunthaler, A. 2007, A&A, 474, 515Miyake, K., & Nakagawa, Y. 1993, Icarus, 106, 20Nakagawa, Y., Nakazawa, K., & Hayashi, C. 1981, Icarus, 45, 517Nakamoto, T., & Nakagawa, Y. 1994, ApJ, 421, 640Natta, A., Testi, L., Calvet, N., Henning, T., Waters, R., & Wilner, D. 2007,Protostars and Planets V, 767Ogihara, M., Ida, S., & Morbidelli, A. 2007, Icarus, 188, 522Okuzumi, S. 2009, ApJ, 698, 1122Ormel, C. W., & Cuzzi, J. N. 2007, A&A, 466, 413Papaloizou, J., & Lin, D. N. C. 1984, ApJ, 285, 818Regaly, Z., Juhasz, A., Sandor, Z., & Dullemond, C. P. 2011, arXiv:1109.6177Ricci, L., Testi, L., Natta, A., Neri, R., Cabrit, S., & Herczeg, G. J. 2010a, A&A,512, A15Ricci, L., Testi, L., Natta, A., & Brooks, K. J. 2010b, A&A, 521, A66Ricci, L., Mann, R. K., Testi, L., Williams, J. P., Isella, A., Robberto, M., Natta,A., & Brooks, K. J. 2011, A&A, 525, A81Rodmann, J., Henning, T., Chandler, C. J., Mundy, L. G., & Wilner, D. J. 2006,A&A, 446, 211Sch¨afer, C., Speith, R., & Kley, W. 2007, A&A, 470, 733Semenov, D., Henning, T., Helling, C., Ilgner, M., & Sedlmayr, E. 2003, A&A,410, 611Schegerer, A. A., & Wolf, S. 2010, A&A, 517, A87Testi, L., Natta, A., Shepherd, D. S., & Wilner, D. J. 2001, ApJ, 554, 1087Testi, L., Natta, A., Shepherd, D. S., & Wilner, D. J. 2003, A&A, 403, 323Uribe, A. L., Klahr, H., Flock, M., & Henning, T. 2011, ApJ, 736, 85Varni`ere, P., & Tagger, M. 2006, A&A, 446, L13Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2009, ApJ, 702,1490Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2011, ApJ, 737,36Weidenschilling, S. J. 1977, MNRAS, 180, 57Wilking, B. A., Gagn´e, M., & Allen, L. E. 2008, Handbook of Star FormingRegions, Volume II, 351Wilner, D. J., Ho, P. T. P., Kastner, J. H., & Rodr´ıguez, L. F. 2000, ApJ, 534,L101Wilner, D. J., D’Alessio, P., Calvet, N., Claussen, M. J., & Hartmann, L. 2005,ApJ, 626, L109Yang, C.-C., & Menou, K. 2010, MNRAS, 402, 2436Youdin, A. N., & Shu, F. H. 2002, ApJ, 580, 494Zsom, A., & Dullemond, C. P. 2008, A&A, 489, 931Zsom, A., Ormel, C. W., G¨uttler, C., Blum, J., & Dullemond, C. P. 2010, A&A,513, A57, B., & Liska, R. 2001, ApJ, 551, 874Lommen, D., Maddison, S. T., Wright, C. M., van Dishoeck, E. F., Wilner, D. J.,& Bourke, T. L. 2009, A&A, 495, 869Lovelace, R. V. E., Li, H., Colgate, S. A., & Nelson, A. F. 1999, ApJ, 513, 805Lynden-Bell, D., & Pringle, J. E. 1974, MNRAS, 168, 603Menten, K. M., Reid, M. J., Forbrich, J., & Brunthaler, A. 2007, A&A, 474, 515Miyake, K., & Nakagawa, Y. 1993, Icarus, 106, 20Nakagawa, Y., Nakazawa, K., & Hayashi, C. 1981, Icarus, 45, 517Nakamoto, T., & Nakagawa, Y. 1994, ApJ, 421, 640Natta, A., Testi, L., Calvet, N., Henning, T., Waters, R., & Wilner, D. 2007,Protostars and Planets V, 767Ogihara, M., Ida, S., & Morbidelli, A. 2007, Icarus, 188, 522Okuzumi, S. 2009, ApJ, 698, 1122Ormel, C. W., & Cuzzi, J. N. 2007, A&A, 466, 413Papaloizou, J., & Lin, D. N. C. 1984, ApJ, 285, 818Regaly, Z., Juhasz, A., Sandor, Z., & Dullemond, C. P. 2011, arXiv:1109.6177Ricci, L., Testi, L., Natta, A., Neri, R., Cabrit, S., & Herczeg, G. J. 2010a, A&A,512, A15Ricci, L., Testi, L., Natta, A., & Brooks, K. J. 2010b, A&A, 521, A66Ricci, L., Mann, R. K., Testi, L., Williams, J. P., Isella, A., Robberto, M., Natta,A., & Brooks, K. J. 2011, A&A, 525, A81Rodmann, J., Henning, T., Chandler, C. J., Mundy, L. G., & Wilner, D. J. 2006,A&A, 446, 211Sch¨afer, C., Speith, R., & Kley, W. 2007, A&A, 470, 733Semenov, D., Henning, T., Helling, C., Ilgner, M., & Sedlmayr, E. 2003, A&A,410, 611Schegerer, A. A., & Wolf, S. 2010, A&A, 517, A87Testi, L., Natta, A., Shepherd, D. S., & Wilner, D. J. 2001, ApJ, 554, 1087Testi, L., Natta, A., Shepherd, D. S., & Wilner, D. J. 2003, A&A, 403, 323Uribe, A. L., Klahr, H., Flock, M., & Henning, T. 2011, ApJ, 736, 85Varni`ere, P., & Tagger, M. 2006, A&A, 446, L13Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2009, ApJ, 702,1490Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2011, ApJ, 737,36Weidenschilling, S. J. 1977, MNRAS, 180, 57Wilking, B. A., Gagn´e, M., & Allen, L. E. 2008, Handbook of Star FormingRegions, Volume II, 351Wilner, D. J., Ho, P. T. P., Kastner, J. H., & Rodr´ıguez, L. F. 2000, ApJ, 534,L101Wilner, D. J., D’Alessio, P., Calvet, N., Claussen, M. J., & Hartmann, L. 2005,ApJ, 626, L109Yang, C.-C., & Menou, K. 2010, MNRAS, 402, 2436Youdin, A. N., & Shu, F. H. 2002, ApJ, 580, 494Zsom, A., & Dullemond, C. P. 2008, A&A, 489, 931Zsom, A., Ormel, C. W., G¨uttler, C., Blum, J., & Dullemond, C. P. 2010, A&A,513, A57