© ESO, 2013

1. Introduction

Since the first confirmed discovery of brown dwarf (BD) Teide 1 (Rebolo et al. 1995) and Gliese 229B (Nakajima et al. 1995) in 1995, several hundred BDs have been identified and many efforts have focused on understanding these objects, which are considered to be an intermediate step between planets and stars. Observations of BD show near-infrared excess emission (e.g. Muench et al. 2001; Liu et al. 2003) that reveals material around young BDs. Moreover, typical fluxes measured with millimeter observations (e.g Klein et al. 2003; Scholz et al. 2006; Joergens et al. 2012) are in most of the cases lower than few mJy at 1 mm, implying that these circumstellar disks have masses of few M_Jup or even lower. Determining whether these low-mass disks can be the scene of the formation of planetesimals or even planets is still a subject of debate.

Observationally, some aspects of disks around BD are different from those around T-Tauri and Herbig Ae/Be stars. They have a lower accretion rate (~10^-12  M_⊙  yr^-1, Herczeg et al. 2009), a flat tendency for the disk geometry inferred from Spitzer observations and spectral energy distribution (SED) modeling (e.g. Apai et al. 2005; Allers et al. 2006; Guieu et al. 2007; Scholz et al. 2007; Morrow et al. 2008; Pascucci et al. 2009; Szűcs et al. 2010), and a longer lifetime (Carpenter et al. 2006; Riaz et al. 2012; Harvey et al. 2012b).

To study dust growth in protoplanetary disks, different mechanisms should be taken into account, for instance turbulent mixing, settling, aerodynamical drag with the gas, collision rates and fragmentation. Essentially, when the dust particles are small, they are well coupled to the gas, move along with it, and grow as a consequence of surface forces. However, when the particles increase in size, they start to decouple from the gas and the relative velocities between particles increase, leading to fragmentation collisions (Weidenschilling 1977; Brauer et al. 2008). In addition, before any meter-sized object can be formed at Earth-Sun distances, dust drifts toward the central star because of the sub-Keplerian velocity of the gas. Millimeter (mm) grains experience the same rapid inward migration in the outer regions of the disks, even though mm-size particles have been observed in those regions of protoplanetary disks (e.g. Wilner et al. 2000; Natta et al. 2004; Rodmann et al. 2006; Ricci et al. 2010a,b, 2011; Ubach et al. 2012). Local pressure maxima in disks have been proposed as a solution of this rapid inward drift (Klahr & Henning 1997; Johansen et al. 2009; Pinilla et al. 2012a).

Birnstiel et al. (2010b) showed that under typical T-Tauri parameters, neglecting radial drift and considering different parameters for the disk such as turbulence, particles could reach millimeter sizes, meaning that the spectral index α_1−3  mm can have low values. Birnstiel et al. (2010b) models predict that for a disk mass lower than 5 M_Jup, such as the mass of a BD disk, α_1−3  mm would be close to 3. However, recent millimeter observations confirmed low values α_1−3  mm for two BD disks (Ricci et al. 2012, 2013), with α_1−3  mm ≈ 2.3. In addition, as we show below, the radial drift and fragmentation barriers are different for particles in disks around BD than for those around more massive and luminous stars. Radial drift indeed has a stronger influence on particles in BD disks, and as a consequence, any mechanism in a disk that may allow a reduced inward migration of grains has to be more extreme in BD disks to lead to an effective trapping of particles. Explaining how the first pebbles are formed from interstellar dust in BD disks is therefore a very intriguing topic.

The purpose of this paper is to investigate whether the dust growth models that were successful for T-Tauri disks are, when applied to BD disks, consistent with a short set of observations at millimeter wavelengths (two measurements from recent ALMA and CARMA observations). As mentioned above, infrared studies have on the other hand focused on larger surveys, showing interesting trends such as a flatter BD disk geometry. However, we do not aim to reproduce these trends, because infrared observations yield no strong constraints on the dust evolution models (other than those studied in e.g. Szűcs et al. 2010; or Mulders & Dominik 2012).

For this work, dust coagulation/fragmentation models were considered in two main cases: first, the extreme case of setting the radial drift to zero, and second, taking radial drift into account with strong inhomogeneities in the gas surface density that mimic long-lived pressure bumps. In both cases, different disk parameters were considered, to analyze which scenarios are the best incubators of the first pebbles found in BD disks with millimeter observations. In Sect. 2, we describe the drift and fragmentation barriers for the specific case of BD disks and the physical parameters of the dust coagulation/fragmentation model. Numerical results, observational perspectives, and comparisons with current mm-observations are presented in Sect. 3. A summary of the results with the corresponding discussion of this work is recapitulated in Sect. 4. Finally, the main conclusion of this paper is given in Sect. 5.

2. Dust evolution model, drift, and fragmentation barriers in BD disks

The interaction between the gas and the dust is fundamental for understanding how the particles evolve within the protoplanetary disk. When particles are well coupled to the gas, the dust relative velocities are mainly generated by Brownian motion and settling to the midplane. Considering these two sources for the velocities, particles stick by van der Waals forces and collisional growth is very efficient, producing mm-size grains in the outer regions (≳50 AU) on timescales of ~10⁵ years in typical T-Tauri disks (Birnstiel et al. 2010a). Unfortunately, when particles grow, the relative velocities substantially increase through turbulent motion and radial drift; the collision energies are therefore high enough to cause fragmentation, as was experimentally shown by Blum & Wurm (2008). From the theoretical point of view, stellar and disk properties strongly influence these first steps of planet formation, meaning that effects such as the usual inward migration of dust grains may vary for BD disks. In this section, we first explain the main characteristics of the numerical model for the dust evolution, followed by an explanation of the drift and fragmentation barriers in BD disks, and finally we describe the set-up for the numerical simulations that are considered in this work.

2.1. Dust evolution model

For the dust evolution, we used the coagulation/fragmentation model explained in Birnstiel et al. (2010a). The dust evolution is described by the advection-diffusion differential equation of the dust surface density Σ_d, which for a single dust size can be written in cylindrical coordinates as (1)where D_d is the dust diffusivity and Σ_g is the gas surface density. Because the timescales for gas viscous evolution are longer than the dust growth timescales, we considered that for the dust evolution models, the gas surface density remains constant with time. This equation was solved for each size using the flux-conserving donor-cell scheme (see Birnstiel et al. 2010a, Appendix A). The radial velocity of the dust u_r,d has two contributions: the first one u_drag comes from the drag with the gas that depends on the radial gas velocity u_r,g and on the size of the particles; and the second one arises from radial drift u_drift, which is proportional to the radial pressure gradient ∂_rP, such that (2)where the coupling constant St – the Stokes number – is defined as the ratio between the largest eddy turn-over time (1/Ω, with ) and the stopping time of the particle within the gas. In the Epstein regime, where the ratio between the mean free path of the gas molecules λ_mfp and the sizes of the particles a is λ_mfp/a   ≥   4/9, St is defined at the disk midplane as (3)where ρ_s is the volume density of a dust grain, usually of about ~1 g cm³.

In addition, the turbulent gas viscosity is considered as (Shakura & Sunyaev 1973), where c_s is the isothermal sound speed, which is given by (4)where k_B is the Boltzman constant, m_p the proton mass and μ the mean molecular mass. Magnetorotational instability (MRI) is the most likely source of turbulence in disks (e.g. Johansen & Klahr 2005). MRI essentially depends on the disk temperature and penetration of cosmic, X- and UV-rays to the midplane. The dust diffusivity in Eq. (1) can be defined in terms of St when the gas diffusivity is considered to be the turbulent gas viscosity ν (Youdin & Lithwick 2007), hence (5)The Stokes number St describes how well the particles are coupled to the gas. When St ≪ 1, the first term dominates in Eq. (2), i.e., u_drag, and as a result these particles move along with the gas, meaning that they have sub-Keplerian velocities. When the particle size increases, the second term, i.e., u_drift starts to dominate and reaches the highest value when St = 1, and as a consequence these are the particles that react strongest to the sub-Keplerian velocity of the gas. The bodies with St ≫ 1 move with their own velocity, i.e., Keplerian speed.

Fig. 1 Vertically integrated dust density distribution (Eq. (6)) after different times of evolution and including radial drift. Case of the BD disk ρ-Oph102 parameters (MBD = 0.05  M⊙, LBD = 0.03  L⊙ and TBD = 2880 K) and Rout = 30 AU, Σ = Σ0  r-0.5, Mdisk = 2  MJup, vf = 10 m s-1 and αturb = 10-3. The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density. The cyan line corresponds to the largest size that particles can reach given a fragmentation velocity vf.

Open with DEXTER

In addition, dust particles grow, fragment, and crater depending on the relative velocities between them. For this the Smoluchowski coagulation equation (Smoluchowski 1916) was solved for the dust grain distribution n(a,r,z) (see Birnstiel et al. 2010a, Eqs. (35) and (36)), considering three different physical processes: coagulation, fragmentation, and erosion. For the relative velocities, we assumed Brownian motion, settling, turbulent motion (Ormel & Cuzzi 2007), and drift velocities in the azimuthal and radial direction. The fragmentation velocities v_f were estimated based on laboratory experiments and theoretical work of collisions, which are of the order of few m s^-1 for silicates (Blum & Wurm 2008) and several m s^-1 for ices (e.g., Wada et al. 2009, 2011). To describe the dust grain distribution, we refer to the vertically integrated dust surface density distribution per logarithmic bin, which is given by (6)where m is the mass of a single particle of size a. Therefore the total dust surface density can be written as (7)

2.2. High radial drift and fragmentation in BD disks

In a protoplanetary disk, large bodies feel no pressure and move with Keplerian speed, while the gas is slightly sub-Keplerian. As a result, particles that have grown to a size where they start to decouple from the gas, i.e., St ≈ 1, feel a strong headwind, lose angular momentum, and move inward. Solving the radial Navier-Stokes equation, the azimuthal velocity of the gas is given by (Nakagawa et al. 1986) (8)where u_k is the Keplerian speed and (9)with ρ_g(r,z) the gas density, such that the gas surface density is given by . For an ideal gas, the pressure P(r) is defined as (10)and therefore (11)where the term in parenthesis depends only on the exponents that characterize the power-law gas density and temperature radial profiles in a flared disk. For the temperature, T(r) can be approximated by (Kenyon & Hartmann 1987) (12)where r/R_⋆ ≳ 10². α_inc is the angle between the incident radiation and the local disk surface. For the last proportionality of Eq. (12) the weak dependence of α_inc with R_⋆ is neglected. Considering Eq. (12) for the disk temperature and assuming a simple power-law function for the gas density profile, the difference between the orbital gas velocity and the Keplerian speed scales as (independently of the radial location r) (13)This difference determines how much angular momentum the particles lose, and as a consequence, how fast they drift to the central star. Taking Eq. (2) for the dust radial velocity, and considering grains with St ≳   1, the drag term (first term of Eq. (2)) can be neglected and u_r,d can be written in terms of the difference u_φ,g − u_k as (14)For comparison, taking the parameters of the BD know as ρ −Oph 102 (Ricci et al. 2012) and comparing them with a young solar-type star, i.e., M_{ρ − Oph102} = 0.05   ×   M_sun and L_{ρ − Oph102} = 0.03   ×   L_sun, this leads to (15)The fact that the radial drift is higher for BD disks than for Sun-like disks also implies that destructive collisions due to radial drift are more likely because it contributes to increase the relative velocities between the particles. Therefore rapid inward drift and the potential fragmentation due to these high drift velocities are major problems for particles in BD disks.

In addition, the highest turbulent relative velocity of particles with a given Stokes number is given by (Ormel & Cuzzi 2007) (16)When the radial drift is set to zero, the largest size of particles a_max is calculated when the turbulent relative velocities (Eq. (16)) are equal to the fragmentation velocities v_f. For particles with St ≲ 1, a_max is approximately given by (17)This implies that considering only turbulent velocities, a_max will depend on disk parameters such as α_turb, Σ_g, and T(r), which strongly vary between BD and T-Tauri disks. Whether turbulence or radial drift is the cause for destructive collisions for the dust within the disk, fragmentation would occur differently for BD disks, as we discuss in the following sections.

As a general illustration of this problem in BD disks, Fig. 1 shows the vertically integrated dust density distribution (Eq. (6)) after different times of evolution, including radial drift and using the dust evolution model described in Sect. 2.1. For this case, taking the parameters ρ-Oph 102 (M_BD = 0.05  M_⊙, L_BD = 0.03  L_⊙ and T_BD = 2880 K), a truncated power law is taken for the gas density profile with Σ = Σ₀  r^-0.5 and r   ∈   (1,30) AU such that the mass of the disk is M_disk = 2  M_Jup. The largest size that particles can reach is represented by the solid cyan line, which was computed considering the fragmentation velocity as v_f = 10 m s^-1 and the turbulence parameter as α_turb = 10^-3. The solid white line corresponds to St = 1, which reflects the shape of the gas surface density based on Eq. (3) and particles that feel the strongest radial drift. The dust particles initially grow in the disk, allowing to have mm size grains in the outer regions (r   ≳   10 AU) after 0.1 Myr of evolution (left panel of Fig. 1). However, when particles grow to sizes close to St ~   1, turbulence and radial drift lead to fragmentation and inward migration of particles, and in only 0.5 Myr of evolution (middle panel of Fig. 1), the outer region is empty of mm-grains. This scenario does not change significantly with time (right panel of Fig. 1), but the inner region (r ≲ 10 AU) is even more depleted of mm grains, since they continue fragmenting and drifting toward the star, leading to a dust-poor disk after 1.5 Myr. As a result, explaining how micron-size dust grows to pebbles and how those are retained in the outer regions of BD disks is very challenging.

In addition, if perfect sticking and micron-sized compact particles are considered to be formed, the mean growth time of monomers with mass m at a fixed distance from the star is (Brauer et al. 2008; Okuzumi et al. 2011) (18)where dm/dt = ρ_dσ_collΔv, with ρ_d is the dust density, σ_coll the collision cross-section and Δv the collision velocity. For a narrow-size distribution, perfect mixing and settling assumed in the Epstein regime, a given dust-to-gas ratio, the growth timescale is proportional to (19)

with Σ_d as the dust surface density. The growth timescale in disks of BD and T-Tauri disks can be compared considering, for example, regions of similar temperature. Following Mulders & Dominik (2012), we re-scaled the distances with the stellar luminosity as well as the Keplerian frequency, to obtain (20)where and correspond to the growth timescales in BD and T-Tauri disks respectively. This implies that for micron-size particles and considering for simplicity only settling motion, the mean collision time τ_grow is generally shorter for BD than for T-Tauri disks. As an illustration, if we again take the BD ρ-Oph 102, the growth timescale due to settling is approximately one order of magnitude shorter for dust particles around this BD than for particles in a Sun-like disk. As a consequence of the result that settling and radial drift occur faster in BD disks, the very early stages of dust growth, when particles just stick and grow through molecular interactions, is more efficient at the location of a given temperature in BD disks than in T-Tauri disks.

Based on the low spectral indices (α_{1 − 3  mm} < 3) measure with recent ALMA-Cycle 0-observations at 3.2 mm and 0.89 mm of the young BD ρ-Oph 102 (Ricci et al. 2012) and CARMA observations of the BD disk 2M0444+2512 at 3 mm (Ricci et al. 2013) and 0.850 mm (Bouy et al. 2008), we aim to explain how mm-grains can form and be retained in disks around BD. We divide this work in two parts: In the first part, we do not allow particles to drift to study whether we can create mm-grains considering only grain growth and fragmentation. In the second part, we aim to explain how these grains are trapped in the outer regions of BD disks when radial drift is included.

2.3. Set-up

For the simulations of this work, we considered two different scenarios. In the first case, there is no dust radial evolution, i.e., the drag and drift terms in Eq. (2) are neglected. In the second case, we considered dust radial evolution with pressure bumps. For the gas density, we considered truncated power-law functions Σ_g = Σ₀  r^− p, with p =  { 0.0,0.5,1.0 } and Σ₀ was computed such that the mass of the disk was always M_disk = 2  M_Jup. For the outer radius of the disk, four possibilities were considered: 15, 30, 60 and 100 AU, since the exact typical disk spatial extent is unknown (Luhman et al. 2007; Ricci et al. 2012).

The turbulence parameter α_turb was taken to have values between 10^-5−10^-3. This parameter directly influences the maximum size of particles (Eq. (17)), if turbulence is the cause of fragmentation. Mulders & Dominik (2012) inferred α_turb ~ 10^-4 from SED modeling, assuming a fixed grain size distribution and gas-to-dust ratio. The fragmentation velocities were assumed to be v_f =  { 10,30 } m s^-1, which are the values expected for ices (Wada et al. 2009). For all simulations, the dust-to-gas ratio was initially considered to be 1%, the initial size of the particles was taken to be 1 μm and the maximum size that particles can reach in the simulations was fixed to be 1 km, because for these sizes gravitational effects start to play a role, which are not included in the present dust evolution model.

Fig. 2 Vertically integrated dust density distribution after 1 Myr of evolution without radial drift. Case of Rout = 30 AU and vf = 10 m s-1 and different values of the gas density slope p = 0 (top panels), p = 0.5 (middle panels) and p = 1.0 (bottom panels); with the turbulence parameter αturb = 10-3 (left panels), αturb = 10-4 (center panels) and αturb = 10-5 (right panels). The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density, while the cyan line represents the maximum size that particles can reach given a fragmentation velocity vf.

Open with DEXTER

Fig. 3 Vertically integrated dust density distribution after 1 Myr of evolution without radial drift. Case of Rout = 30 AU and vf = 30 m s-1 and different values of the gas density slope p = 0 (top panels), p = 0.5 (middle panels) and p = 1.0 (bottom panels); with the turbulence parameter αturb = 10-3 (left panels), αturb = 10-4 (center panels) and αturb = 10-5 (right panels). The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density.

Open with DEXTER

When we considered pressure inhomogeneities, we assumed the same prescription as in Pinilla et al. (2012a) to simulate long-lived pressure bumps, with the unperturbed density as a simple power law (21)where the wavelength of the perturbation L(r) is taken to be a factor f of the vertical disk scale-height h(r) i.e. L(r) = f  h(r), with (22)Pinilla et al. (2012a) showed that under these assumptions, bumps with an amplitude of ~30% and a wavelength of one scale-height are necessary for T-Tauri disks, to reduce the radial drift and keep the particles in the outer regions of the disk after several Myr. However, in Sect. 2.2, we showed that radial drift effects are stronger for BD disks, independently of the distance of the particles from the star. In addition, since the wavelength of the perturbation was considered proportional to the disk scale-height h(r), at a given distance r, L(r) would by definition be larger for BD than for T-Tauri disks (Eq. (22)). This implies that for a given amplitude of the perturbed density, the pressure gradient would be lower in BD disks, making the trapping of particles even more difficult. We used a fixed bump width equal to one scale height, i.e., f = 1.0, since the scale of turbulent structures from MRI are of the order of the scale height of the disk (Flock et al. 2011), and they may be the origin of pressure inhomogeneities in the disk. The amplitudes considered are therefore higher than for T-Tauri disks and these are taken A = { 0.4,0.6 }. All stellar and disk parameters are summarized in Table 1.

3. Results

In this section, we present the results of the numerical simulations setting the radial drift to zero, followed by the case of non-zero radial drift and pressure bumps in the disk. For each case, the observational perspectives are also presented.

3.1. No radial drift

3.1.1. Dust density distribution

Figures 2 and 3 show the vertically integrated dust density distribution (Eq. (6)) after 1 Myr of dust evolution when the radial drift is set to zero, for the case of R_out = 30 AU and two different values of the fragmentation velocity: 10 m s^-1 and 30 m s^-1. In each case, the three different values of the gas density slope p and the turbulence parameter α_turb are plotted. For this case, the dust is considered to be in a steady state, therefore u_drag and u_drift are neglected in Eq. (2).

Effect of turbulence:

when particles grow and reach sizes such that they are more weakly coupled to the gas, the main sources of relative velocities are turbulent and azimuthal velocities, since radial drift velocities are set to zero. Because the highest turbulent relative velocities (Eq. (16)) depend on α_turb and the disk temperature is assumed to be very low in BD disks (~10 K at 10 AU), fragmentation is not caused by turbulent motions. Instead, it mainly happens because the azimuthal dust velocities u_d,φ are as high as the limit beyond which particles fragment (fragmentation velocity, v_f). The azimuthal drift velocity u_d,φ is given by (Birnstiel et al. 2010a) (23)Hence, as we noticed in Fig. 2, the fragmentation limit and consequently a_max are independent of the α_turb and for a given gas density slope, the largest size that particles can reach is the same for each α_turb and only depends on the gas density profile. On the other hand, we notice that the concentration of mm- and cm-size grains is more evenly distributed in the whole disk for α_turb =  { 10^-4,10^-3 } than for α_turb = 10^-5, for which they are mainly in the inner part of the disk.

The first remarkable difference of grain growth between T-Tauri and BD disks is that the turbulent mixing strength α_turb does not play as important a role for dust fragmentation in BD disks as it does for T-Tauri disks. Without radial drift, fragmentation in T-Tauri disks mainly occurs because of turbulent motion and the effect of different α_turb values is significant (see e.g. Birnstiel et al. 2010b, for the effect of turbulence in the spectral index), while for BD disks the azimuthal relative velocity dominates, which is independent of α_turb.

	Fig. 4 Predicted fluxes at 1 mm (F1  mm) and the spectral index between 1 and 3 mm (α1 − 3  mm) after 1 Myr of dust evolution, without radial drift, vf = 10 m s-1 and for all other parameters discussed in Sect. 2.3. Black dots with error bars are millimeter observations of the young BD ρ-Oph 102 (Ricci et al. 2012) and 2M0444+2512 (Ricci et al. 2013).
Open with DEXTER

Effect of the gas density slope:

for the same turbulence parameter and different gas density slope p, the dust density distribution changes (see Figs. 2 and 3). When p = 1 (steep surface density), large particles (a > 1 cm) are well coupled to the gas in the dense inner region and grow to even larger sizes (a ~ 50 cm) before they fragment. Therefore, with p = 1 the mm-size particles are less concentrated in the region where the gas density is lower, i.e., in the outer region of the disk (r ≳ 10 AU), than for p =  { 0.5,1.0 }.

Effect of the fragmentation velocity:

when the velocity at which particles fragment increases to 30 m s^-1 (Fig. 3), grains can grow to very large sizes such that St > 1. Fragmentation does not occur and this implies that bodies that are no longer affected by gas (a   ≳   1 m) grow and even in some cases reach the maximum size of particles considered in these numerical simulations (~1 km). Small particles are not replenished and as a result the disk can be almost empty of small grains (≲1 cm). For example, when p = 0 and α_turb = 10^-5, most of the disk consists only of big grains (a > 10 cm). However, when p increases and the gas density is higher in the inner regions, larger particles are coupled to the gas and can still be affected by gas turbulence. For this reason, there is some mm-dust in the inner part for the same α_turb = 10^-5 and p = 1.

As we discussed above, because of the low disk temperatures, destructive collisions due to turbulence are less likely in BD disks, and even more with these high fragmentation velocities considered for ices (such as v_f = 30 m s^-1). With ices, the material properties vary, leading to fragmentation velocities of around ~10–50 m s^-1 (Wada et al. 2009, 2011), but it is still a matter of debate if ices can have these high fragmentation velocities. Some scenarios produce disks with only large grains (Fig. 3), which may not allow a good agreement between these theoretical models and the recent submillimeter observations of disk around two BD. For this reason, we focus the following results on the case of v_f = 10 m s^-1, where the fragmentation allows a constant replenishment of small particles.

Effect of the outer radius:

when the outer radius of the disk decreases or increases, the same amount of dust is distributed in a smaller or larger region with the same M_disk and Σ_g profile, and the possibility of having cm- mm-size particles in the outer regions changes, we discuss this effect below.

3.1.2. Comparison with observations

Fig. 5 Vertically integrated dust density distribution after 5 Myr of evolution including radial drift and a bumpy gas surface density (Eq. (21)) with A = 0.4 (left panel) and A = 0.6 (right panel). Case of Rout = 30 AU, vf = 10 m s-1, p = 0.5 and αturb = 10-3. The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density, while the cyan line corresponds to the fragmentation limit.

Open with DEXTER

To compare model predictions with disk fluxes of low-mass disks, we considered the dust distribution from the simulations after they almost reached a quasi-static state ~1 Myr. We calculated the opacities κ_ν for each grain size and at a given frequency ν. We used for simplicity optical constants for magnesium-iron silicates (Jaeger et al. 1994; Dorschner et al. 1995) from the Jena database¹ and followed the Mie theory. Once the opacities were calculated, the optical depth τ_ν was computed as (24)where i is the disk inclination, which we took to be zero. The flux of the disk at a given frequency F_ν is therefore (25)with d being the distance to the source, which was taken to be 140 pc as in the young disks in Taurus and Ophiuchus star-forming regions. B_ν(T(r)) is the Planck function for a given temperature profile T(r). At mm wavelengths, the flux is proportional to the dust mass in the outer region of the disks. Hence, the flux F_ν could be approximated as a power law F_ν ∝ ν^α_mm, where the spectral index α_mm gives information about the size of the grains and is expected to be lower than 3 when the dust reaches mm sizes (Natta et al. 2007). Figure 4 shows the predicted fluxes at 1 mm (F_1  mm) and the spectral index between 1 and 3 mm (α_{1 − 3  mm}) for v_f = 10 m s^-1 and all other parameters considered in Sect. 2.3.

First of all, comparing F_1mm and α_1 − 3mm for different R_out, we notice that the spectral slope increases when the disk is more extended. This is a natural result of distributing initially the same amount of dust in an extended region, decreasing the possibility of having mm-size grains. The obtained fluxes are between ~1–20 mJy for all cases. In addition, a single case (the only diamond-point of Fig. 4) where fragmentation does not happen is plotted, with v_f = 30 m s^-1, R_out = 30 AU and p = 0.0. It is possible to see that most of the dust in the disk has grown to large sizes (a > 10 cm). When the particles have such a large size, they have very low opacities, which results in low millimeter fluxes. The predicted fluxes are very low (close to 1 mJy) compared with the other fluxes (see top-left panel of Fig. 3 and the corresponding diamond-point in the second-left panel from the top to the bottom of Fig. 4).

In general, when fragmentation still occurs, i.e., for v_f = 10 m s^-1, the integrated flux weakly depends on the gas density slope p, therefore we focus in the subsequent results on the intermediate value, i. e., p = 0.5. In most of the cases, the spectral slope and the flux are slightly sensitive to the turbulence parameter α_turb, and α_{1 − 3  mm} increases for low turbulence. This is because there are fewer mm-size grains that contribute to decrease the spectral index when α_turb is low (Figs. 2 and 3), as was discussed in Sect. 3.1.1.

Comparing the millimeter observations (dot points with error bars in Fig. 4) of two BD disks (2M0444+2512 and ρ Oph 102), it is important to notice that when the radial drift is set to zero, the predicted and the observed values of F_1  mm and α_{1 − 3  mm} match well for the brightest BD disk and in particular for R_out = 15 AU. The error bars come from the the optical depth uncertainties for different grain composition (Beckwith et al. 1990). When the extension of the disk increases to R_out = 30 AU, the spectral slope has values that agree with observations (especially for higher values of α_turb), but the fluxes are slightly higher than those detected by mm-observations. This leads to the conclusion that in this case, a combination of fragmentation with a minor drift is necessary to reduce the number of mm-grains and have a better agreement between theoretical predictions and recent millimeter observations of 2M0444+2512. For the faintest disk ρ Oph 102, radial drift is indeed needed in any case. For R_out = 60 AU, only few cases allow to have α_{1 − 3  mm} ≲ 3, nevertheless, these cases match the millimeter observations less well. Finally, for most of the cases of R_out = 100 AU, the model predictions for the millimeter fluxes and spectral indices are inadequate to explain the observations.

As a conclusion, theoretical models of dust evolution in which the radial drift is set to zero and considering BD disk conditions such as low-mass disks (M_disk = 2  M_Jup), low radial extension (R_out = 15 AU), presence of ices that allow fragmentation velocities of about v_f = 10 m s^-1 and average turbulence strength α_turb =  { 10^-4,10^-3 }, are the models with the best agreement to mm-observations.

3.2. Radial drift and pressure bumps

We focus our attention on the most favorable cases: R_out =  { 15,30,60 } AU (particularly 15 and 30 AU) and α_turb = 10^-3. The gas density slope was taken to be p = 0.5. The amplitudes considered for the long-lived pressure bumps were A =  { 0.4,0.6 } (Eq. (21)). Pinilla et al. (2012a) showed that these inhomogeneities are comparable with global simulations of zonal flows with an amplitude of ~25% by Uribe et al. (2011).

For the simulations with radial drift, u_drag and u_drift for Eq. (2) were taken into account. This implies that radial drift also contributes to the total relative velocities of dust particles.

3.2.1. Dust density distribution

In a region where the pressure gradient is positive, gas moves with super-Keplerian velocity and particles with sizes corresponding to St ~ 1 would move outward. As a consequence, inside the pressure bumps the amount of dust naturally increases, allowing an increase in the frequency of sticking collisions. Figure 5 compares the dust density distribution after 1 Myr of evolution for the case of R_out = 30 AU, v_f = 10 m s^-1, p = 0.5, α_turb = 10^-3 and two different values of the amplitude of the perturbation A = 0.4 and A = 0.6. In both cases, with the pressure bumps considered, it is possible to reduce the strong radial drift and retain mm-size particles even after million-year timescales. The efficiency of the pressure bumps is evident, with many more mm grains than in the test case of Fig. 1 (empty disk after 0.5 Myr of dust evolution, Fig. 1). In addition, the effectiveness of the amplitude strength, i.e., the pressure gradient strength, for the trapping of particles is almost the same for both amplitudes after 1 My of evolution.

	Fig. 6 Dust-to-gas ratio at different times of dust evolution including radial drift and a bumpy gas surface density (Eq. (21)) with A = 0.4. Case of αturb = 10-3, vf = 10 m s-1, p = 0.5, Rout = 30 AU (top panel) and Rout = 60 AU (bottom panel).
Open with DEXTER

Inside the bumps, the radial drift is reduced and fragmentation is also less efficient than in the case of Fig. 1. With a fragmentation velocity of v_f = 10 m s^-1, particles grow to sizes corresponding to St ~ 1 for the cases in Fig. 5, and because of the positive pressure gradient, grains are trapped inside the bumps. However, since fragmentation due to turbulent and azimuthal relative velocities is still happening, when particles collide and become smaller such that St < 1, they are more difficult to trap in the pressure bump (Pinilla et al. 2012a) and would finally drift to the star.

Fig. 7 Time evolution of predicted fluxes at 1 mm (F1  mm) and the spectral index between 1 and 3 mm (α1 − 3  mm) including radial drift and a bumpy gas surface density (Eq. (21)) with Rout =  { 15,30 }  AU, αturb = 10-3, vf = 10 m s-1, A = 0.4 (left panels), A = 0.6 (right panels), and p = 0.5. Black dots with error bars are millimeter observations of the young BD ρ-Oph 102 (Ricci et al. 2012) and 2M0444+2512 (Ricci et al. 2013).

Open with DEXTER

3.2.2. Dust-to-gas ratio

Figure 6 shows the radial dependence of the dust-to-gas ratio at different times of evolution for the cases of α_turb = 10^-3, v_f = 10 m s^-1, and p = 0.5, in a bumpy gas surface density with A = 0.4 and R_out =  { 30,60 } AU. For this case, the dust-to-gas ratio remains almost constant with time until 1 Myr of evolution. It increases inside the pressure bumps and varies between ~10^-2 − 10^-1, but there is a remarkable reduction after 5 Myr of evolution, where it can reach values of 10^-3 in the inner part and 10^-4 − 10^-5 in the outer parts of the disk.

When the outer radius increases, the dust-to-gas ratio does not remain constant with time in the outer regions of the disk, and it is possible to distinguish that it decreases with time, reaching even values of ~10^-6 in those regions (bottom panel of Fig. 6).

It is important to note that the gas density profile remains constant for our dust evolution models, since viscous evolution timescales are longer than dust growth timescales. In addition, we did not consider any mechanism that might selectively disperse the gas, such as photoevaporation (see e.g. Gorti et al. 2009; Owen et al. 2011).

3.2.3. Comparison with millimeter observations

Focusing the attention on the cases where fragmentation is effective (i.e., when v_f = 10 m s^-1 and α_turb = 10^-3), Fig. 7 shows the time evolution of F_1  mm and α_{1 − 3  mm} including radial drift and a bumpy gas surface density with R_out =  { 15,30 }  AU, p = 0.5, and A =  { 0.4,0.6 }. For A = 0.4, we notice a good match between the theoretical predictions and the observations, even after 1 Myr of evolution, there are still enough mm-grains in the outer regions that contribute to have high values of the flux at 1 mm. The fluxes become lower at 5 Myr, because a significant amount of grains continues to fragment and drift toward the star. However, F_1  mm and α_{1 − 3  mm} at intermediate time steps like 1–2 Myr coincide well with the observed values, and are in the range of the estimated ages of the observed BD disks. The effect of increasing the pressure gradient, i.e., higher A, is notable, since the fluxes increase because there are more trapped grains, in particular at long timescales, for instance 5 Myr.

To achieve a better agreement of the predicted fluxes F_1  mm and the spectral index α_{1 − 3  mm} with the current mm-observations of 2M0444+2512 and ρ Oph 102, it would be necessary to increase the mass of the disk. Figure 8 shows the comparison of F_1  mm and the α_{1 − 3  mm} for two different disk masses (M_disk =  { 2,5 }   M_Jup) and the specific case of R_out = 15 AU, p = 0.5, α_turb = 10^-3, v_f = 10 m s^-1 and A = 0.4. Because there is more dust for M_disk = 5  M_Jup than for M_disk = 2  M_Jup, there are more mm-grains that contribute to the integrated flux, and therefore the spectral slope can be lower.

4. Summary and discussion

One of the most important problems in the core accretion theory for planet formation is “the meter-size barrier”. This phenomenon results from the combination of high radial drift and fragmentation. A meter-sized object at 1 AU drifts toward the star on timescales much shorter than the growth timescales. In the colder regions of a disk, the same problem occurs when translated to mm-sized grains. Both radial drift and fragmentation appear to be different for particles in BD disks.

We have studied how dust evolves and explored different scenarios where mm-grains survive in BD disks. For the results presented here, we considered typical BD parameters. For instance, the disk outer regions reach temperatures of ~10 K. Assuming that MRI is the source of turbulence, which depends basically on the disk temperature and considering α_turb-prescription, low and moderate values of α_turb were taken (10^-5 − 10^-3). In addition, the presence of ices in the disk was assumed and fragmentation velocities were taken to be v_f =  { 10,30 } m s^-1, based on numerical simulations of collisions with ices (see e.g. Wada et al. 2009). Furthermore, we considered low-mass disks (2  M_Jup) and typical values for masses and luminosities of BD (see Table 1). For the radial extension of BD disks, we considered R_out =  { 15,30,60,100 } AU. Luhman et al. (2007) presented observations combining Spitzer spectroscopy and high-resolution imaging from Hubble Space Telescope (HST) of a circumstellar disk that is inclined close to edge-on around a young BD in Taurus, estimating a disk radius of ~20–40 AU. ALMA-Cycle 0 observations of ρ-Oph 102 (Ricci et al. 2012) do not resolve the disk, but they estimate an outer radius for the dust of ~15–40 AU. From 3-mm CARMA observations of 2M0444+2512 (Ricci et al. 2013), the outer radius for the dust is estimated to be ~20–75 AU for p = 1 and ~15–30 AU for p = 0.5.

	Fig. 8 Comparison of the predicted fluxes at 1 mm (F1  mm) and the spectral index between 1 and 3 mm (α1 − 3  mm) for two different disk masses (Mdisk = 2  MJup-star points and Mdisk = 5  MJup-dot points). Time-evolving data from right to left (from 0.1–5.0 Myr as Fig. 7). Case of Rout = 15 AU, p = 0.5, αturb = 10-3, vf = 10 m s-1, and A = 0.4.
Open with DEXTER

4.1. Fragmentation and drift barriers

We have shown how several dust growth aspects in BD disks differ from T-Tauri and Herbig Ae/Be disks. In the first stages of dust evolution, when particles move along with the gas and coagulate, relative velocities are mainly due to Brownian motion and settling to the midplane. It was discussed in Sect. 2.2 that growth timescales due to settling are shorter at the location of a given temperature in BD than in T-Tauri disks.

When particles reach sizes such that they start to decouple from the gas (i.e. St ~ 1), radial drift and turbulent motion become the main sources for relative velocities. On one hand, we demonstrated in Sect. 2.2 that radial drift can be even twice as fast for particles in BD disks than in T-Tauri disks, and as a consequence they would be depleted on shorter timescales. On the other hand, considering moderate α_turb and low disk temperatures, turbulent relative velocities are quite low. If turbulence is the main cause for fragmentation, destructive collisions are less likely in BD disks than in T-Tauri disks. When radial drift and fragmentation are considered, BD disks would be dust-poor after only few Myr because of the fast inward migration of dust grains (Fig. 1).

4.2. No radial drift

This drastic case, where radial drift is not included, is important for studying whether mm-grain can be indeed formed considering dust coagulation/fragmentation models. When radial drift was set to zero, we found that two different scenarios can happen under BD disk conditions. In one scenario fragmentation occurs because of the relative azimuthal velocities, considering v_f = 10 m s^-1 (Fig. 2). In this case, it is possible to form mm-size particles in BD disk and achieve a good match for the predicted mm-fluxes and spectral indices, in particular for the brightest BD disk 2M0444+2512 (Fig. 4) for R_out = 15 AU. However, for larger disk sizes, i.e., R_out =  { 30,60 }  AU, a slight radial drift is necessary to decrease the predicted mm-fluxes. For the other BD disk ρ Oph 102, radial drift is always needed to reduce the mm-fluxes by at least a factor of 3.

The second scenario is when fragmentation is less likely, i.e., v_f = 30 m s^-1 (see e.g. Fig. 3), when most of the disk consists of big grains (a > 10 cm), leading to low mm-fluxes that underestimate millimeter observations. Therefore, values such as v_f = 30 m s^-1 do not yield a good agreement between models and observations, and average values of v_f = 10 m s^-1 are needed.

4.3. Radial drift and pressure bumps

To reduce the rapid inward drift that dust particles experience in protoplanetary disks and explain the presence of mm-grains in the outer regions of disks, pressure bumps have been suggested as a possible solution (see e.g. Pinilla et al. 2012a). In this work, when radial drift was taken into account, we considered strong gas density inhomogeneities (40% and 60% of amplitude), that led to regions in the disk with a positive pressure gradient. These pressure bumps could be the result of MRI effects (Johansen et al. 2009; Uribe et al. 2011). However, there are no known mechanisms that could produce this type of pressure bumps with these high amplitudes globally present in the disk. Nevertheless, if the scale-height of a disk is higher, as in the case of BD disks, the scale of turbulent structures increases, leading to higher pressure scales (e.g. Flock et al. 2011).

In addition to global pressure inhomogeneities from MRI in the disk, a single and high-pressure bump can exist in disks. For instance, Rossby wave instability (RWI) generated by the presence of a dead zone (Regály et al. 2012) would create a high single-pressure bump. The possibility of a dead zone in a BD disk is still debated. On the other hand, the presence of a massive planet could be also the reason for trapping of particles in a single huge pressure bump (Pinilla et al. 2012b). However, from low-mass disks, it is unlikely to have planets massive enough to open a gap. The criterion to open a gap in a disk depends on the local scale height h(r), meaning that a clear gap would be less likely for BD than T-tauri disks for the same disk and planet parameters (Crida et al. 2006). Hence, these two possibilities were ruled out for the conditions that we considered for BD disks, and we only focused on global-synthetic pressure bumps that allow a broad radial distribution of millimeter dust grains in the whole disk.

For simulations with pressure bumps, we considered a sinusoidal perturbation for the gas surface density with two different amplitudes A =  { 0.4,0.6 } (Eq. (21)) and a wavelength equal to one disk scale-height. The parameters considered for these simulations were taken from the most optimistic no-drift cases. The pressure bumps cause an accumulation of dust in the location of pressure maxima, where radial drift is reduced and fragmentation due to radial drift decreases. Inside those regions, fragmentation would happen only because turbulent or azimuthal motions and a good match between theoretical perspectives and spectral slopes found with millimeter observations for the two BD disks is possible, especially for R_out =  { 15,30 }  AU. Increasing the mass of the disk from 2 to 5 M_Jup, theoretical predictions and mm-observations match even better. The disk mass predictions for ρ-Oph 102 (Ricci et al. 2012) are 0.2–0.6 M_Jup and for 2M0444+2512 (Ricci et al. 2013) are 2.0–5.0 M_Jup, assuming a gas-to-dust mass ratio of 100.

4.4. Comparison between T-Tauri and BD disks

Although T-Tauri and BD disks have very different properties, a first-order comparison between the models which agree well with millimeter observations in each case can be made. If we compare the best model described in Pinilla et al. (2012a) for T-Tauri disks (pressure bumps with an amplitude A of 30% and one scale-height of the disk as the width f of the bumps) with the best model for BD disks (A = 0.4, f = 1.0, R_out = 15 AU and M_disk = 2  M_Jup), considering the same disk and dust properties (fragmentation velocity (v_f = 10 m s^-1), viscosity (α_turb = 10^-3), initial dust-to-gas ratio (1/100), etc.), we find that 85% of the dust mass is in grains larger than 1 mm in the BD case, versus 70% for T-Tauri disks. As a consequence, we expect our models to produce slightly flatter SEDs for BD disks than for TTS disks. We stress that this conclusion, which is inferred from our 1D models that are only relevant in the mid-plane, should be taken with care. A proper vertical structure calculation is needed to assess this question, and see if the effect is strong enough to account for 15–20% variation reported by e.g. Mulders & Dominik (2012) on Chameleon I targets.

4.5. Additional improvements

The aim of this paper was to determine the best conditions for grain growth in the mid-plane of BD disks. Our 1D model is well-suited to this goal, but requires to azimuthally average the surface density profiles, which implies that we do not take into account possible local disk features that may locally affect dust growth. 2D and 3D calculations would allow the vertical disk structure and vertical grain distribution to be computed and provide us predictions at shorter wavelengths (e.g. shape of the SED and of the silicate feature) that can be directly compared with large surveys.

Our model can also be improved in the prescription of the dust-evolution process. Recent laboratory experiments using silicates have revealed that particles probably also bounce (Güttler et al. 2010), which is not included in our model so far. It is an open question, however, whether this happens in the case of collision between ices (Wada et al. 2011), and whether it has a strong impact on dust growth. In fact, the bouncing barrier (Zsom et al. 2010) can be overcome if Maxwellian velocity distributions of relative velocities among dust grains are considered (Windmark et al. 2012a). Windmark et al. (2012b) also showed that after including bouncing effects and considering mass-transfer collisions, km-sized objects can still form at few AU from the star, but only if radial drift is neglected and cm-sized seeds are inserted initially.

Finally, the main caveat of our study is that we assumed that the gas density profile remains static for million-year timescales, neglecting the impact of different physical processes, for example the potential time evolution of the pressure bumps, like the case of zonal flows (see e.g.  Uribe et al. 2011). Simultaneous modeling of gas and dust evolution is a step forward of our models.

5. Conclusion

The first steps of planet formation differ between disks around T-Tauri stars and disks around BD. In BD disks, settling is more efficient and they are expected to be flatter. In addition, when dust grows and particles start to decouple from the gas, radial drift is a major problem and fragmentation due to turbulence is less likely than in T-Tauri disks. Our models show that when BD disks are small (e.g., R_out ≲ 15AU), radial drift needs to be completely suppressed to account the mm-observations of one of the two BD disks observed so far. When disks are larger, a small degree of radial drift is necessary for both cases. Including radial drift, the most favorable conditions after million-year timescales of dust evolution are with strong pressure inhomogeneities of an amplitude of ~ 40%–60%. Global dust evolution models that include the vertical structure are the next step to reproduce all observational signatures of BD disks, for example the flatter geometry of BD compared with T-Tauri disks, or the large average grain sizes in BD disk atmospheres.

Acknowledgments

We thank the anonymous referee for his/her constructive report, which helped us to improve and clarify the main results of the paper. T. Birnstiel acknowledges support from NASA Origins of Solar Systems grant NNX12AJ04G. P. Pinilla acknowledges the CPU time for running simulations in bwGRiD, member of the German D-Grid initiative, funded by the Ministry for Education and Research (Bundesministerium für Bildung und Forschung) and the Ministry for Science, Research and Arts Baden-Wuerttemberg (Ministerium für Wissenschaft, Forschung und Kunst Baden-Württemberg).

References

All Tables

All Figures

	Fig. 1 Vertically integrated dust density distribution (Eq. (6)) after different times of evolution and including radial drift. Case of the BD disk ρ-Oph102 parameters (MBD = 0.05  M⊙, LBD = 0.03  L⊙ and TBD = 2880 K) and Rout = 30 AU, Σ = Σ0  r-0.5, Mdisk = 2  MJup, vf = 10 m s-1 and αturb = 10-3. The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density. The cyan line corresponds to the largest size that particles can reach given a fragmentation velocity vf.
Open with DEXTER
In the text

	Fig. 2 Vertically integrated dust density distribution after 1 Myr of evolution without radial drift. Case of Rout = 30 AU and vf = 10 m s-1 and different values of the gas density slope p = 0 (top panels), p = 0.5 (middle panels) and p = 1.0 (bottom panels); with the turbulence parameter αturb = 10-3 (left panels), αturb = 10-4 (center panels) and αturb = 10-5 (right panels). The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density, while the cyan line represents the maximum size that particles can reach given a fragmentation velocity vf.
Open with DEXTER
In the text

	Fig. 3 Vertically integrated dust density distribution after 1 Myr of evolution without radial drift. Case of Rout = 30 AU and vf = 30 m s-1 and different values of the gas density slope p = 0 (top panels), p = 0.5 (middle panels) and p = 1.0 (bottom panels); with the turbulence parameter αturb = 10-3 (left panels), αturb = 10-4 (center panels) and αturb = 10-5 (right panels). The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density.
Open with DEXTER
In the text

	Fig. 4 Predicted fluxes at 1 mm (F1  mm) and the spectral index between 1 and 3 mm (α1 − 3  mm) after 1 Myr of dust evolution, without radial drift, vf = 10 m s-1 and for all other parameters discussed in Sect. 2.3. Black dots with error bars are millimeter observations of the young BD ρ-Oph 102 (Ricci et al. 2012) and 2M0444+2512 (Ricci et al. 2013).
Open with DEXTER
In the text

	Fig. 5 Vertically integrated dust density distribution after 5 Myr of evolution including radial drift and a bumpy gas surface density (Eq. (21)) with A = 0.4 (left panel) and A = 0.6 (right panel). Case of Rout = 30 AU, vf = 10 m s-1, p = 0.5 and αturb = 10-3. The solid white line represents the particle size corresponding to St = 1 (Eq. (3)) and reflects the shape of the gas density, while the cyan line corresponds to the fragmentation limit.
Open with DEXTER
In the text

	Fig. 6 Dust-to-gas ratio at different times of dust evolution including radial drift and a bumpy gas surface density (Eq. (21)) with A = 0.4. Case of αturb = 10-3, vf = 10 m s-1, p = 0.5, Rout = 30 AU (top panel) and Rout = 60 AU (bottom panel).
Open with DEXTER
In the text

	Fig. 7 Time evolution of predicted fluxes at 1 mm (F1  mm) and the spectral index between 1 and 3 mm (α1 − 3  mm) including radial drift and a bumpy gas surface density (Eq. (21)) with Rout =  { 15,30 }  AU, αturb = 10-3, vf = 10 m s-1, A = 0.4 (left panels), A = 0.6 (right panels), and p = 0.5. Black dots with error bars are millimeter observations of the young BD ρ-Oph 102 (Ricci et al. 2012) and 2M0444+2512 (Ricci et al. 2013).
Open with DEXTER
In the text

	Fig. 8 Comparison of the predicted fluxes at 1 mm (F1  mm) and the spectral index between 1 and 3 mm (α1 − 3  mm) for two different disk masses (Mdisk = 2  MJup-star points and Mdisk = 5  MJup-dot points). Time-evolving data from right to left (from 0.1–5.0 Myr as Fig. 7). Case of Rout = 15 AU, p = 0.5, αturb = 10-3, vf = 10 m s-1, and A = 0.4.
Open with DEXTER
In the text