Dust retention in protoplanetary disks

T. Birnstiel; C. P. Dullemond; F. Brauer

doi:10.1051/0004-6361/200912452

Free Access

Issue		A&A Volume 503, Number 1, August III 2009


Page(s)		L5 - L8
Section		Letters
DOI		https://doi.org/10.1051/0004-6361/200912452
Published online		15 July 2009

LETTER TO THE EDITOR

Dust retention in protoplanetary disks

T. Birnstiel - C. P. Dullemond - F. Brauer

Junior Research Group at the Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany

Received 8 May 2009 / Accepted 3 July 2009

Abstract
Context. Protoplanetary disks are observed to remain dust-rich for up to several million years. Theoretical modeling, on the other hand, raises several questions. Firstly, dust coagulation occurs so rapidly, that if the small dust grains are not replenished by collisional fragmentation of dust aggregates, most disks should be observed to be dust poor, which is not the case. Secondly, if dust aggregates grow to sizes of the order of centimeters to meters, they drift so fast inwards, that they are quickly lost.
Aims. We attempt to verify if collisional fragmentation of dust aggregates is effective enough to keep disks ``dusty'' by replenishing the population of small grains and by preventing excessive radial drift.
Methods. With a new and sophisticated implicitly integrated coagulation and fragmentation modeling code, we solve the combined problem of coagulation, fragmentation, turbulent mixing and radial drift and at the same time solve for the 1D viscous gas disk evolution.
Results. We find that for a critical collision velocity of 1 m s^-1, as suggested by laboratory experiments, the fragmentation is so effective, that at all times the dust is in the form of relatively small particles. This means that radial drift is small and that large amounts of small dust particles remain present for a few million years, as observed. For a critical velocity of 10 m s^-1, we find that particles grow about two orders of magnitude larger, which leads again to significant dust loss since larger particles are more strongly affected by radial drift.

Key words: accretion, accretion disks - circumstellar matter - stars: formation - stars: pre-main-sequence - infrared: stars

1 Introduction

The quest for a comprehensive understanding of the formation of planets in the dusty circumstellar disks around many young stars has gained significant impetus. This is partly because of developments in the numerical modeling of these processes, but also because of the high-quality infrared and (sub-)millimeter data that are now available of hundreds of T Tauri and Herbig Ae star disks from observatories such as the Spitzer Space Telescope (e.g. Kessler-Silacci et al. 2006; Furlan et al. 2006), the Very Large Telescope (e.g. van Boekel et al. 2004), the Submillimeter Array (e.g. Andrews & Williams 2007) and the Plateau de Bure Interferometer (e.g. Piétu et al. 2007). The challenge for theoretical astrophysicists is now not only to come up with a plausible model of how to overcome the various problems in our current understanding of how planets form, but also to ensure that new models are consistent with the observational data obtained from these planetary birthplaces.

One of the challenges of theories of planet formation is what is often called the ``meter-size barrier''. The first stages of planet formation are believed to involve the coagulation of dust particles to aggregates of ever larger size, and eventually to the size of asteroids or comets (Weidenschilling 1980,1984,1997; Nakagawa et al. 1981). It was shown by Weidenschilling (1977) that as dust particles grow by coagulation, they acquire increasingly large relative velocities with respect to other particles in their vicinity, as well as a systematic inward drift velocity as originally described by Whipple (1972). The high relative velocities lead to destructive collisions (Blum & Wurm 2008), and thereby limit the growth of the aggregates to some maximum size. The radial drift also leads to a loss of solids toward the star well before large bodies can even form (Brauer et al. 2008a, hereafter Brauer:2008p215).

Dullemond & Dominik (2005) showed that if, hypothetically, coagulation can proceed without limit, i.e., if no fragmentation occurs and no radial drift is present, the growth is so efficient that within 10⁵ years very few small ( $a\la 1$ mm) opacity bearing grains remain in the disk. The disk becomes optically thin and the mid- to far-infrared flux of the disk drops to levels well below that observed in the majority of gas-rich circumstellar disks around pre-main sequence stars. By excluding other mechanisms, Dominik & Dullemond (2008) argue that fragmentation is the most promising process for maintaining a high abundance of fine-grained dust in evolved disks.

However, it has never been demonstrated explicitly that fragmentation is efficient enough to retain the population of small grains at the required levels. Moreover, if we include the effect of radial drift, it remains to be seen whether this radial drift is suppressed sufficiently by keeping the grains small through fragmentation. It is the goal of this Letter to investigate this issue and to determine how efficient fragmentation should be to agree with observations.

2 Model

The model presented here is a combination of a 1D viscous gas disk evolution code and a dust evolution code treating the turbulent mixing, radial drift, coagulation and fragmentation of the dust. A detailed description of the model will follow in a later paper (Birnstiel et al., in prep.). Although the evolution of dust can influence the gas disk evolution, we neglect this effect in our model.

The gas disk model is similar to that described in Hueso & Guillot (2005). The mid-plane temperature of the disk is approximated by following Nakamoto & Nakagawa (1994), taking irradiation by the central star and viscous heating into account. The viscous evolution of the disk surface density $\Sigma_{\rm g}$ is described by (see Lynden-Bell & Pringle 1974),

$\begin{displaymath}\frac{\partial\Sigma_{\rm g} }{\partial t} - \frac{1}{r}\frac... ...}{\partial r}\left( \Sigma_{\rm g} ~ r ~ u_{\rm g}\right) = 0, \end{displaymath}$

(1)

where the gas radial velocity $u_{\rm g}$ is given by

$\begin{displaymath}u_{\rm g} = - \frac{1}{\Sigma_{\rm g}\sqrt{r}} \frac{\partial}{\partial r} \left( \Sigma_{\rm g} \nu_{\rm g} \sqrt{r} \right). \end{displaymath}$

(2)

$\nu_{\rm g}$ is the gas turbulent viscosity,

$\begin{displaymath}\nu_{\rm g} = \alpha ~ c_{\rm s} ~ H_{\rm p}, \end{displaymath}$

(3)

$c_{\rm s}$ denotes the sound speed, $H_{\rm p}$ the pressure scale height, and $\alpha$ is the turbulence parameter (Shakura & Sunyaev 1973).

Grains are affected by a systematic radial inward drift related to headwind caused by the pressure-supported gas (Weidenschilling 1977),

$\begin{displaymath}u_{\rm drift} = \frac{1}{\ensuremath{{\rm St}}\xspace^{-1}+\e... ...ac{ \partial_{r}{P_{\rm g}}}{\rho_ {\rm g} \: \Omega_{\rm k}}, \end{displaymath}$

(4)

where $\ensuremath{{\rm St}}$ is the Stokes number of the particle (a dimensionless representation of the particle size), $P_{\rm g}$ is the gas pressure, $\rho_{\rm g}$ is the gas volume density, and $\Omega_{\rm k}$ is the Kepler frequency. In the Epstein regime, $\ensuremath{{\rm St}}$ is given by

$\begin{displaymath}\ensuremath{{\rm St}}\xspace= \frac{a \rho_{\rm s}}{\Sigma_{\rm g}} \frac{\pi}{2}, \end{displaymath}$

(5)

where $\rho_{\rm s}$ is the solid density of the dust grains.

The second contribution to radial dust velocities is the gas drag due to the radial motion of the gas,

$\begin{displaymath}u_{\rm drag} = \frac{u_{\rm g}}{\ensuremath{{\rm Sc}}\xspace}, \end{displaymath}$

(6)

where $\ensuremath{{\rm Sc}}\xspace=1 + \ensuremath{{\rm St}}\xspace^2$ is the Schmidt number, following Youdin & Lithwick (2007).

The coupling to the gas turbulence leads to turbulent mixing of each dust species with a diffusion constant that is taken to be

$\begin{displaymath}D_{\rm d} = \frac{\nu_{\rm g}}{\ensuremath{{\rm Sc}}\xspace}\cdot \end{displaymath}$

(7)

The dust grain number density n(m,r,z) is a function of mass m, radius r, height above the mid-plane z, and time. If we define

$\begin{displaymath}\sigma(m,r) \equiv \int_{-\infty}^{\infty} n(m,r,z) ~ m^2 \ensuremath{{\rm d}} z, \end{displaymath}$

(8)

the vertically integrated time-evolution of this distribution can now be described by a general two-body process and an advection-diffusion equation,

$\displaystyle \frac{\partial \sigma(m,r)}{\partial t}$	=	$\displaystyle \iint_{0}^{\infty} K(m,m',m'') ~ \sigma(m',r) ~ \sigma(m'',r)~ \ensuremath{{\rm d}} m'~ \ensuremath{{\rm d}} m''$
		$\displaystyle -\frac{1}{r}\frac{\partial}{\partial r} \left[ r \cdot \sigma(m,r) \cdot (u_{\rm drag}+u_{\rm drift}) \right]$
		$\displaystyle +\frac{1}{r}\frac{\partial}{\partial r} \left[ r \cdot D_{\rm d} ... ... \left( \frac{\sigma(m,r)}{\Sigma_{\rm g}}\right) \cdot \Sigma_{\rm g} \right].$	(9)

Here, the right-hand side terms (from top to bottom) correspond to coagulation/fragmentation, advection, and turbulent mixing.

Since the detailed version of the coagulation/fragmentation equation is lengthy Weidenschilling 1980; Nakagawa et al. 1981; Dullemond & Dominik 2005; Brauer:2008p215, we include here only the general integral of a two-body process (first term on the right hand side). For a detailed description of the physics of coagulation/fragmentation, we refer to Brauer:2008p215; the numerical implementation will be discussed in Birnstiel et al. (in prep.).

The physical effects that produce the relative particle velocities considered here are Brownian motion, relative radial motion, vertical settling [see][]Brauer:2008p215, and turbulent motion (see Ormel & Cuzzi 2007).

If particles collide with a relative velocity higher than the critical velocity $\ensuremath{u_{\rm f}}$ , they either fragment into a power-law size distribution of fragments [i.e., $n(m)\propto m^{-1.83}$ , see][]Brauer:2008p215 or the smaller body excavates mass from the larger one (cratering). We assume that the amount of excavated mass equals the mass of the smaller body. The fragmentation velocity $\ensuremath{u_{\rm f}}$ is a free parameter in our model, and unless otherwise noted is taken to be 1 m s^-1, as suggested by experiments (e.g., Blum & Muench 1993) and by theoretical modeling (e.g., Stewart & Leinhardt 2009).

We note that the velocities and therefore also the relative velocities are dependent on the particle size. If we assume that relative particle velocities are dominated by turbulent motion,

$\begin{displaymath}\Delta u_{\rm turb} \propto \sqrt{\alpha\:\ensuremath{{\rm St}}\xspace} \: c_{\rm s}, \end{displaymath}$

(10)

then we can estimate the Stokes number of the largest particles by equating the relative velocity and fragmentation velocity,

$\begin{displaymath}\ensuremath{{\rm St}}\xspace_{\rm max} \simeq \frac{\ensuremath{u_{\rm f}}\xspace^2}{\alpha\:c_{\rm s}^2}\cdot \end{displaymath}$

(11)

With Eq. (5), we can relate this maximum ``dimensionless size'' $\ensuremath{{\rm St}}\xspace_{\rm max}$ to a particle size,

$\begin{displaymath}a_{\rm max} \simeq \frac{2\Sigma_{\rm g} }{\pi \alpha \rho_{\... ...\cdot \frac{\ensuremath{u_{\rm f}}\xspace^2}{c_{\rm s}^2}\cdot \end{displaymath}$

(12)

It should be noted that Eq. (10) is an approximation of the equations derived in Ormel & Cuzzi (2007). The detailed relative velocities of particles depend on the Stokes numbers of both particles. Hence, Eq. (12) is an order of magnitude estimate.

The combined equations for coagulation, fragmentation, and radial motion due to drift, drag, and mixing are solved using the technique of implicit integration [similar to][]Brauer:2008p215. We use the SixPACK F-90 linear algebra package to solve the sparse matrix equation.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{12452fg1.eps} \end{figure}$

Figure 1:

Snapshots of the vertically integrated dust density distribution of model 1 (only coagulation), model 2 (fragmentation at 10 m s^-1) and model 3 (fragmentation at 1 m s^-1). The grain size is given by $a=(3m/4\pi \rho _{\rm s})^{1/3}$ . The dashed line shows the evaporation radius within which no coagulation is calculated. The solid line shows the particle size corresponding to a Stokes number of unity.

Open with DEXTER

3 Results

The simulation results presented here begin with a 0.07 $M_\odot$ disk of power-law surface density $\Sigma(r)\propto r^{-1}$ (from 0.05 to 150 AU) around a 0.5 $M_\odot$ star. Unless otherwise noted, we used a turbulence parameter of 10^-2. The three models differ in the prescription of fragmentation. In model 1, we accounted only for growth without fragmentation. In model 2, particles fragment at $\ensuremath{u_{\rm f}}\xspace=10$ m s^-1, whereas in model 3, the fragmentation speed is taken to be 1 m s^-1. Selected results of all three models are depicted in Fig. 1.

The top row shows the evolution in the surface density distribution of dust because of coagulation, radial mixing, and radial drift in an evolving protoplanetary disk, corresponding to model 1. The dashed line denotes the grain evaporation radius, which moves inwards as the temperature falls (lower surface density caused by accretion produces less viscous heating). The evaporation radius in our model is defined as the radius where the temperature rises above 1500 K assuming that all particles evaporate at this temperature. The solid line denotes the dust grain size, which corresponds to a Stokes number of unity. It can be seen that particles grow to decimeter size and then quickly drift inside the evaporation radius.

The results of simulations that also include grain fragmentation are plotted in the second and third row of Fig. 1. In the innermost regions of model 2, relative particle velocities are high enough (as a result of shorter dynamical times) to fragment particles already at relatively small Stokes number. However, the maximum size of the grain distribution is still strongly affected by radial drift, especially at larger radii. Larger particles are subject to radial drift and therefore drift to smaller radii, where they are pulverized.

In contrast, fragmentation becomes the main limiting factor for particle growth throughout the disk if the critical fragmentation velocity is 1 m s^-1 (model 3). The result is that particles fragment at a small Stokes number. Particles with a Stokes number less than about 10^-3 are still strongly coupled to the gas and not subject to strong radial drift (see Brauer et al. 2007). Therefore, significant amounts of dust can be retained for several million years if fragmentation is included in the calculations.

Figure 2 compares the evolution in the dust-to-gas ratio (solid line) and the 0.55 $\mu$ m optical depth at 10 AU (dashed line) of all three simulations. The optical depth is defined as

$\begin{displaymath}\tau_\nu = \int_{0}^\infty \sigma(m) ~ \kappa_\nu(m) \ensuremath{{\rm d}} m. \end{displaymath}$

(13)

It can be seen that the dust-to-gas ratio declines dramatically in the first two cases. The optical depth in model 1 drops on even shorter timescales than the dust-to-gas ratio because not only the dust loss but also the dust growth causes the optical depth to decrease.

Model 2 can retain the optical depth longer because small particles remain (due to fragmentation), providing enough surface to keep the disk optically thick for about 0.9 Myr. The initial dip in the optical depth is caused by the initial condition: since initially, all mass is in small grains, the optical depth decreases as particles grow to larger sizes until fragmentation sets in, which increases the optical depth until a semi-equilibrium of growth and fragmentation is reached.

If particles are already fragmented at 1 m s^-1 (model 3), the dust-to-gas ratio declines only slightly, after 2 My, only about 4% of the dust being lost (compared to 96% and 94% in models 1 and 2). The optical depth decreases on longer timescales than in model 2 since the dust mass is lost only on the viscous timescale.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{12452fg2.eps} \end{figure}$	Figure 2: Dust-to-gas ratio (solid line, linear scale) and 0.55 $\mu$ m optical depth at 10 AU (dashed line, logarithmic scale) as function of time for model 1 ( top), model 2 ( center) and model 3 ( bottom).
Open with DEXTER

This mechanism of grain retention relies on the particles remaining small due to fragmentation. From Eqs. (11) and (12), it follows that the maximum size or Stokes number depends strongly on the critical velocity $\ensuremath{u_{\rm f}}$ and less strongly on $\alpha$ . If $\ensuremath{u_{\rm f}}$ is 10 m s^-1 instead of 1 m s^-1, the corresponding particle size is 100 times larger meaning that the largest particles are now also experiencing the radial drift barrier and therefore drifting quickly towards the central star.

4 Discussion and conclusion

We have shown that grain fragmentation plays an important role, not only in keeping the disk optically thick (by replenishing small grains) but also in maintaining a relatively high abundance of solids of all sizes. If particles fragment at relatively small sizes, they can stay below the size-regime where they are affected by radial drift.

The simulation results presented here suggest that fragmentation of relatively low critical velocity (a few m s^-1) is needed to retain the dust required for planet formation, although the formation of planetesimals probably requires periods (and/or regions) of quiescence such as dead zones or pressure bumps (see Brauer et al. 2008b; Kretke & Lin 2007). These locally confined quiescent environments could decrease both the relative and radial velocities of particles, thus allowing particles possibly to either grow to larger sizes or to accumulate until gravitational instabilities become important (see Lyra et al. 2009; Johansen et al. 2007). These effects are not yet included in this model and will be the subject of future research.

Assuming that the proposed mechanism of dust retention is effective in protoplanetary disks, we can constrain the critical fragmentation velocity to be smaller than 10 m s^-1 since velocities $\ga$ 10 m s^-1 cannot explain the observed lifetimes of disks even in the case of small $\alpha$ .

Acknowledgements

We like to thank Thomas Henning, Chris Ormel and Andras Zsom for useful discussions. We also wish to thank the anonymous referee for his fast and insightful report which helped to improve the paper.

References

Andrews, M. R., & Williams, J. P. 2007, ApJ, 671, 1800 [NASA ADS] [CrossRef] (In the text)
Blum, J., & Muench, M. 1993, Icarus, 106, 151 [NASA ADS] [CrossRef] (In the text)
Blum, J., & Wurm, G. 2008, ARA&A, 46, 21 [NASA ADS] [CrossRef] (In the text)
Brauer, F., Dullemond, C. P., Johansen, A., et al. 2007, A&A, 469, 1169 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Brauer, F., Dullemond, C. P., & Henning, T. 2008a, A&A, 480, 859 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Brauer, F., Henning, T., & Dullemond, C. P. 2008b, A&A, 487, L1 [NASA ADS] [CrossRef] [EDP Sciences]
Dominik, C., & Dullemond, C. P. 2008, A&A, 491, 663 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Dullemond, C. P., & Dominik, C. 2005, A&A, 434, 971 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Furlan, E., Hartmann, L., Calvet, N., et al. 2006, ApJS, 165, 568 [NASA ADS] [CrossRef]
Hueso, R., & Guillot, T. 2005, A&A, 442, 703 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Johansen, A., Oishi, J. S., Low, M.-M. M., et al. 2007, Nature, 448, 1022 [NASA ADS] [CrossRef]
Kessler-Silacci, J., Augereau, J.-C., Dullemond, C. P., et al. 2006, ApJ, 639, 275 [NASA ADS] [CrossRef]
Kretke, K. A., & Lin, D. N. C. 2007, ApJ, 664, L55 [NASA ADS] [CrossRef]
Lynden-Bell, D., & Pringle, J. E. 1974, MNRAS, 168, 603 [NASA ADS] (In the text)
Lyra, W., Johansen, A., Klahr, H., & Piskunov, N. 2009, A&A, 493, 1125 [NASA ADS] [CrossRef] [EDP Sciences]
Nakagawa, Y., Nakazawa, K., & Hayashi, C. 1981, Icarus, 45, 517 [NASA ADS] [CrossRef]
Nakamoto, T., & Nakagawa, Y. 1994, ApJ, 421, 640 [NASA ADS] [CrossRef] (In the text)
Ormel, C. W., & Cuzzi, J. N. 2007, A&A, 466, 413 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Piétu, V., Dutrey, A., & Guilloteau, S. 2007, A&A, 467, 163 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] (In the text)
Stewart, S. T., & Leinhardt, Z. M. 2009, ApJ, 691, L133 [NASA ADS] [CrossRef] (In the text)
van Boekel, R., Min, M., Leinert, C., et al. 2004, Nature, 432, 479 [NASA ADS] [CrossRef] (In the text)
Weidenschilling, S. J. 1977, MNRAS, 180, 57 [NASA ADS] (In the text)
Weidenschilling, S. J. 1980, Icarus, 44, 172 [NASA ADS] [CrossRef]
Weidenschilling, S. J. 1984, Icarus, 60, 553 [NASA ADS] [CrossRef]
Weidenschilling, S. J. 1997, Icarus, 127, 290 [NASA ADS] [CrossRef]
Whipple, F. L. 1972, From Plasma to Planet, Proc. 21st Nobel Symp., ed. A. Evlius, 211 (In the text)
Youdin, A. N., & Lithwick, Y. 2007, Icarus, 192, 588 [NASA ADS] [CrossRef] (In the text)

Footnotes

... package: Available from www.engineers.auckland.ac.nz/~snor007

All Figures

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{12452fg1.eps} \end{figure}$	Figure 1: Snapshots of the vertically integrated dust density distribution of model 1 (only coagulation), model 2 (fragmentation at 10 m s^-1) and model 3 (fragmentation at 1 m s^-1). The grain size is given by $a=(3m/4\pi \rho _{\rm s})^{1/3}$ . The dashed line shows the evaporation radius within which no coagulation is calculated. The solid line shows the particle size corresponding to a Stokes number of unity.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{12452fg2.eps} \end{figure}$	Figure 2: Dust-to-gas ratio (solid line, linear scale) and 0.55 $\mu$ m optical depth at 10 AU (dashed line, logarithmic scale) as function of time for model 1 ( top), model 2 ( center) and model 3 ( bottom).
Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.