© ESO, 2011

1. Introduction

The gas and dust surface densities of proto-planetary disks appear as one of the key parameters in the formation of planetary systems. For example, the formation mechanism of giant planets remains a debated problem. Competing models are the core-accretion mechanism (e.g. Hubickyj et al. 2005), which faces apparent timescale difficulties, and the gravitational instability (e.g. Boss 1997; Rice et al. 2005), which requires massive disks. Determining the dust and gas densities as a function of age of the proto-planetary disks would be a major step to decide the relative importance of the various processes that potentially lead to planet formation.

However, there is no ideal way to measure these densities. H₂ remains the more abundant molecule in proto-planetary disks but is difficult to observe because it only possesses quadrupolar rotation lines in the mid-IR. The gas column density is thus usually estimated from molecular tracers such as CO or less abundant molecules (Piétu et al. 2007). Uncertainties linked to a poor accuracy on the molecular abundance and its variation across the disk owing to the chemical behavior of the observed molecule usually affect the results (Dutrey et al. 2007). The dust surface density can, in theory, be directly derived from the dust brightness temperature. However, the dust emissivity is still poorly known and the accuracy on the surface density depends on the knowledge of the dust properties (composition, size, etc.) and its radial and vertical variations through the disk. Finally, the dust-to-gas ratio may also vary with radius.

In all cases, high angular resolution is required to derive the surface density profile because the typical size of disks range from 100 AU to 1000 AU. Attempts have also been made in the optical, using scattered light images (Burrows et al. 1996), but they are hampered by the need to extrapolate the density structure from the upper layers to the disk mid-plane. Other methods include silhouette disks against the bright background of HII regions: McCaughrean & O’dell (1996) showed that steep edges (power law exponent  ~ −4.5, or exponential taper) were needed to reproduce the “proplyds” in Orion, but this cannot be extrapolated inward because of the high opacities.

The mm domain is better suited to sample the bulk of the disk. However, the high angular resolution required, at least better than 1′′, implies the use of large mm/submm interferometers. For the dust emission, early attempts include the 3 mm study of Dutrey et al. (1996) with the IRAM array, the 2 mm survey of Kitamura et al. (2002) using NRO, and more recently the 1.3 − 0.8 mm study performed by Andrews & Williams (2007) with the SMA. These studies were interpreted in a simplified framework of truncated power laws for the surface densities.

High-resolution studies for the gas are even more difficult. Using the same simplified model, the CO outer radius is in general found to be much larger than the dust-derived outer radius (e.g. Dutrey et al. 1998; Simon et al. 2000; Isella et al. 2006). This is confirmed through CO isotopologue studies in several sources, such as AB Aur (Piétu et al. 2005), DM Tau, LkCa 15, and MWC 480 (Piétu et al. 2007). Although this may be interpreted as changing dust properties with radius, Hughes et al. (2008) suggested this could also be caused by a different surface density distribution, with an exponentially tapered fall-off of the density with radius. At the resolution of their observations,  ≃1′′, the truncated power law and the softened-edge version are indistinguishable.

A similar approach has been used by Isella et al. (2009) to interpret a  ≃ 0.7′′ resolution 1.3 mm survey with CARMA, and by Andrews et al. (2009) for SMA observations at 0.8 mm.

All these analysis were based on single frequency imaging, although the overall SED is often used to provide additional constraints on the disk parameters. For thermal emission, the only observable is the brightness distribution of the dust at frequency νwhere J_ν(ν,T) is the Planck function. At least, measurements at three different frequencies are required to independently constrain Σ(r), T_d(r) and κ(ν,r). In the mm domain, the dust is mostly optically thin and the Rayleigh-Jeans approximation valid in many cases, (3)To first order, this allows the separation of the evolution of κ(ν,r) from that of Σ(r)T(r) with measurements at two frequencies, only. Resolved images are needed at both wavelengths to remove the degeneracy between an optically thick core and possible radial variations of the spectral index β. Recently, Isella et al. (2010) reported a first resolved dual-frequency study of RY Tau and DG Tau, while Banzatti et al. (2011) published a resolved multi-frequency study of CQ Tau.

In this paper we report on a high angular resolution (0.4 to 1′′), dual-frequency survey of  ~20 of circumstellar disks located in the Taurus-Auriga complex, 8 of which have sub-arcsecond angular resolution at both 2.7 and 1.3 mm. Observations are described in Sect. 2. Section 3 presents the disk models that we used and the analysis we performed using our specifically developed method. In Sect. 4 describes the results of this analysis. The consequences and interpretations are presented in Sect. 5. We then conclude in Sect. 6.

2. Sample, observations, and data reduction

Table 1 indicates the properties of the sources in the sample. The sample contains classical T Tauri stars or late-type HAe stars, single or multiples (in italics), and a few embedded sources with optical jets and molecular outflows like DG Tau, DG Tau-b, HL Tau, and HH 30 (in boldface). Properties were obtained from the quoted literature. For homogeneity, all ages were derived using the Siess et al. (2000) evolutionary tracks, directly from the work of Bertout et al. (2007) when available, or re-derived using the cited estimates of luminosity and spectral types. These stellar ages tend to be somewhat higher (factor 1.5) than derived from other evolutionary tracks (D’Antona & Mazzitelli 1997; Palla & Stahler 1999), although even higher ages can be obtained using the Baraffe et al. (1998) tracks. Note that the evolutionary tracks remain ill constrained, and no available set reproduces the constraints derived from the kinematic masses, see Simon et al. (2000) and the small corrections brought by more accurate measurements of Piétu et al. (2007) and Dutrey et al. (2008). However, all evolutionary tracks produce a similar ordering of the ages, at least in the 0.5 − 1.5 M_⊙ range of masses, which dominate our sample. Because the DG Tau-b luminosity is unknown, its age is completely uncertain. Since it still displays an active molecular outflow, we have tentatively assumed it to be 1 Myr old, but with large uncertainties. For GM Aur, the mass derived by Bertout et al. (2007) is somewhat larger than the kinematically derived value 1.00    ±    0.05   M_⊙from Dutrey et al. (2008). Accordingly, its age may be overestimated by about 50%.

Part of the survey was made by simultaneously observing at 2.7 or 3.4 mm and 1.3 mm in the winter seasons between Nov. 1995 and Oct. 1998 using the dual frequency receivers on Plateau de Bure (see Simon et al. 2000, for a description of these observations). Sources were observed in track-sharing mode, typically six to eight at a time. In all cases, the intensity scale was calibrated by using MWC 349 as flux calibrator. This method ensures an homogeneous calibration across the sample, specially for the spectral index determination as MWC 349 has a precisely characterized spectral index of 0.6.

Additional high angular resolution with 750 m baselines data was collected from Piétu et al. (2006) for MWC 480 and LkCa 15, simultaneously at 110 and 220 GHz. For HH 30 we used the data from Guilloteau et al. (2008).

Higher angular observations (baselines up to 760 m) were also obtained in Feb. 2007 at 1.3 mm, and Feb. 2008 at 2.7 mm, again in track-sharing mode among six to eight sources, with the new dual-polarization, single frequency receivers. MWC 349 served as flux calibrator, but in addition MWC 480 was used as an internal flux-scale consistency check, because it is compact, bright enough and independently measured.

The main survey reaches angular resolution of 0.5 × 0.3′′ at 1.3 mm and a factor of 2 lower at 2.7 mm. Phase stability was good during the main survey observations: most observations are noise-limited, rather than dynamic-range limited. Dynamic range only limits the brightest sources HL Tau, T Tau (which were observed only during the first period) and, to a lesser extent DG Tau and MWC 480, for which the effective noise is twice the thermal noise.

Some sources also have 2.7 mm data from previous studies (Dutrey et al. 1996). In addition, more limited angular resolution data from Schaefer et al. (2009) for Haro 6-13 and Haro 6-33 (1.2 × 0.7′′ resolution) and Chapillon et al. (2008) for MWC 758 and CQ Tau (about 1.3′′ resolution) are also included for completeness.

Table 2 indicates the resulting beam sizes for each source at 230 GHz. The positions indicated in Table 2 are those determined from this study, and are the reference positions for Figs. 1, 2, 5, 7, 8 and G.1−G.22. Because the data span more than 12 years of time, correction for the star proper motions is important. The proper motions were taken from the Ducourant et al. (2005) catalog when available, or determined from our own measurement, as the astrometric accuracy of the Plateau de Bure is high enough to allow measurements to about 2 mas/yr in each direction over a 10 year span when sufficient signal-to-noise is available. The positions are in the J2000.0 system and referred to epoch 2000.0 after correction for proper motion. The positional accuracy is better than 0.05′′.

Figure 1 is a montage of the 1.3 mm images of the survey sources, presented in terms of fraction of the peak flux. Figure 2 is as Fig. 1, but for 2.7 or 3.4 mm, depending on the sources. Robust weighting was used to produce these images. Despite the fairly wide range of angular resolutions (from 0.5    ×    0.3′′ to about 1.5′′), clearly some objects are much more centrally condensed than others. In particular, the most compact sources are the two circumstellar disks in the Haro 6-10 binary.

3. Modeling

3.1. Simple analysis

The measured flux densities at 1.3 mm and around 3 mm are given in Table 3 (considering only baselines shorter than 100 m). They result from a simple elliptical Gaussian fit to the uv data. For the orientations and apparent sizes, all baselines were included. Short baseline data, although adequate to measure the overall flux densities and apparent spectral index α, are not suitable to derive characteristic sizes and even position angles. This is because, to first order, disks have power law distributions of the surface density and temperature and are optically thin at such wavelengths. Thus, when seen at low inclination, (<45° or so), the surface brightness is a power law of the radius and has no characteristic size. This can bias the apparent position angle, since the apparent half-power size only depends on the angular resolution and the exponent of the power law. For nearly edge-on disks (i > 75°), the disk thickness introduces a characteristic size, because the brightness falls off like a Gaussian in this direction, so the position angle is properly recovered. Thus, in general, reliable position angles can only be derived with sufficient angular resolution, i.e. from long baseline fits. These properties can explain the different position angles found by previous authors using lower resolution data (e.g. Dutrey et al. 1996; Kitamura et al. 2002). Note that these biases on the position angles can also affect analysis made with more elaborate disk models: only sufficiently high angular resolution can provide an unbiased determination of this parameter.

On the other hand, for sources with an apparent core-halo structure, such as DM Tau or CI Tau, the long baseline fit tends to represent only the central part and misses substantial flux. The spectral index α₁₀₀ derived from long baseline data (Table 3, Col. 8) is systematically smaller than that from the short baseline fit only (Col. 7). This indicates either a contribution of an optically thick core and/or dust grain evolution.

3.2. Model description

Because the apparent size, orientation, and spectral index may depend on the uv coverage when using a simple Gaussian model, we must analyze the data with more realistic brightness distributions. Because a direct inversion of the brightness profile is impossible, due to the combination of insufficient resolution and the limited signal-to-noise, a global fitting technique using some a-priori model must be used. We therefore analyzed the continuum emission in terms of two “standard” disk models that differ only in the surface density distribution. Model 1 uses a simple truncated power law, Model 2 an exponentially tapered power law with an arbitrarily large outer radius. The surface density is characterized by three parameters plus an inner radius in each model. Our approach is to keep the model parametric and simple to avoid as much as possible biases towards a specific physical model for disks.

In Model 1, the surface density is a simple power law with a sharp inner and outer radius: (4)for R_int < r < R_out.

In Model 2, the density is tapered by an exponential edge: (5)Note that Model 1 derives from Model 2 by simply setting R_c → ∞ and p = γ in the above parametrization. Model 2 is a solution of the self-similar evolution of a viscous disk in which the viscosity is a power law of the radius (with constant exponent in time γ).

With the inner (R_int) and outer (R_out) radii, the disk mass is given by (6)which for small R_int and large R_out yields (7)The simple power law case is recovered for R_c → ∞, by developing to first order in (r/R_c)^2 − γ, (8)One can also used M_d as a free parameter instead of Σ₀, like in Andrews et al. (2009). Equation (6) can also be used to show that R_c is the radius which contains 63% of the disk mass, because M(r < R_c) = M_d(1 − 1/e) = 0.63M_d provided R_out is large enough.

An equivalent parametrization is that described by Isella et al. (2009)(9)The parameterizations using R_t or M_d become ill defined for γ = 2, which makes them less suited for use in a minimization scheme than the simple parametric expression of Eq. (5) (for which only R_c is unconstrained, as the surface density becomes a power law). R_t is related to R_c by (10)R_t/R_c is close to 0.5 for all values of γ below 1, reaches 1 for γ = 1.5, then diverges for γ → 2. In the framework of self-similar viscous evolution (Lynden-Bell & Pringle 1974; Hartmann et al. 1998), it can be shown that R_t is the radius at which the net mass flux changes sign.

In both models, the temperature is assumed to be a simple power law of the radius (11)The disks are thus vertically isothermal. To allow a homogeneous comparison, we used T₁₀₀ = T(100   AU) = 15 K and q = 0.4, except when those parameters can be constrained by the observations. The validity and impact of this assumption will be discussed in Sect. 4.1.

Similar analyses have been used by Kitamura et al. (2002) and Andrews & Williams (2007) for their 2 mm and 0.8 mm data respectively. Most previous studies (Kitamura et al. 2002; Andrews & Williams 2007; Isella et al. 2009) used the thin disk approximation to compute visibilities. Here, because our sample includes highly inclined objects, we assume that the disks are flared, with a scale height varying as a power law of the radius h(r) = H₁₀₀(r/100   AU) ^− h. For all but the two highly inclined objects (HH 30 and DG Tau-b), we used H₁₀₀ = 16 AU and h = −1.25. These values agree with those derived using the gas temperature determined from CO observations whenever available, and the stellar mass, either from kinematic determination (Simon et al. 2000) or standard evolutionary tracks. The results are, however, completely independent of the assumed scale height, which justifies a posteriori the thin disk approximation used by previous authors. However, for the two highly inclined objects, H₁₀₀ and the exponent h had to be used as adjustable parameters.

The inner radius R_int is also not significant in general, except for a few special sources that display inner cavities, such as GM Aur, HH 30 and LkCa 15 (see Sect. 4.4.2). We fixed it to 1 AU, but in general, any value lower than about 3 − 4 AU would not change the results. For Model 2, we used for R_out the outer radius derived from CO observations when available. If not, we set it to 500 AU. These outer radii are large enough to have negligible influence on the results.

Each model has thus a priori five free intrinsic parameters: two for temperature T₀ and q, three for the surface density Σ₀, p or γ, and R_out or R_c, plus the inclination, orientation and position.

The dust opacity as a function of wavelength and radius completes the description. In a first step, we assume it to be independent of radius and described by the following prescription (12)with κ₂₃₀ = 2   cm²   g^-1 (per gram of dust). This introduces one additional parameter, the mean dust emissivity index β_m. This is similar to the Beckwith et al. (1990) results, but using a different pivot frequency to avoid further dependence of the derived disk mass on β_m. The dust model used by Andrews & Williams (2007) and Andrews et al. (2009) also results in β_m = 1, but with a slightly different absorption coefficient κ₂₃₀ = 2.2   cm²   g^-1. Finally, we also assume that the gas-to-dust ratio is constant and equal to 100. In a second step, we shall relax the assumption of constant κ(ν) as a function of radius r, see Sect. 4.5.

Appendix A (available on-line only) illustrates some of the possible degeneracy between the various models, in particular between constant dust properties with an optically thick inner region, and variable dust properties.

3.3. Fitting method

For the inclination and orientation, we used the accurate determination from the CO kinematics when possible. Values derived from optical observations (scattered light images, or optical jets) or molecular jets were used for some sources for which the disk kinematics is not known. Independent fits of these parameters from the dust emission were also performed to check the consistency of the results: see Table 4 and references therein. We stress, however, that the uncertainties on the disk inclination and orientation do not significantly affect the derived radial structure.

At each observed frequency, the radiative transfer equation is solved by a simple ray-tracing algorithm, and model images are generated. Great care has been taken to avoid numerical precision problems caused by finite grid effects. The numerical integration is typically performed on a 128 × 128 grid, with 512 points along the line of sight. Two oversampling techniques are used to enhance the accuracy while keeping computational costs reasonable. First, the overall image is interpolated (by bilinear interpolation) by a factor 2 before computing the model visibilities. Second, the inner 64 × 64 pixels are re-computed on this finer grid with a smaller step along the line of sight (64 × 64 × 1024). Larger numbers were used for the largest disk. This results in effective pixel sizes of 2 to 7 AU in (x,y), depending on the outer disk radius used in the model, and steps 4 to 8 times smaller along the line of sight.

A modified Levenberg-Marquardt method was used to derive the disk parameters by a non-linear least squares fit of the modeled visibilities directly to the observed uv data, as detailed by Piétu et al. (2007). Like all methods, L-M minimization can be trapped in local minima when the starting point is too far away from the solution. We alleviate this problem by using multiple re-starts when needed, and also by adapting the step size used to compute the gradient. We found empirically that using steps equal to half a sigma on each parameter provided stable results. Error bars were derived from the covariance matrix, except when the parameter coupling was too strong (e.g. between R_c and γ in Model 2 for γ larger than about 1.5). In that case, the multi-parameter fit was reduced to a one parameter problem by finding the best fit for several values of this parameter, and determining the error bars from the resulting χ² distribution. Data at several available wavelengths are fitted simultaneously by the same model, which allows us to constrain β_m. However, whenever data at very nearby frequencies (220 and 230 GHz, for example) exist, only one was considered in this process, because even small absolute calibration error could result in a strong bias on the value of β_m. In the dual frequency fit, the long wavelength (2.7 or 3.4 mm) data do not in general influence the derived surface density law, because of their lower angular resolution, but only serve to determine β_m. Because the geometric parameters are largely decoupled from the disk intrinsic parameters, the simultaneous fit of dual-frequency data sets used (in general) four parameters: Σ₀, p (γ in Model 2), R_out (R_c in Model 2), and β_m. Additional parameters (T₀,q or H₀,h) were also fitted simultaneously when needed. Separate fits were also made at 1.3 and 2.7 mm for the few sources were the angular resolution at 2.7 mm is sufficient, or when data sets at 1.4 mm also existed: in these cases, β_m was set at the value found from the dual frequency analysis, and only the three remaining parameters were fitted together.

The choice of the pivot radius R₀ in Eqs. (4, 5) is important. There is always an optimal value that minimizes the error on Σ₀, which depends on the angular resolution and source surface density profile (see discussion in Piétu et al. 2007). Using a non-optimal value results in a coupling of Σ₀ with p for Model 1, and γ,R_c in Model 2. Another different pivot radius is also required for T₀ when the source is sufficiently optically thick and resolved to constrain T₀,q.

Two stars required a specific treatment: the binary Haro 6-10 and the quadruple UZ Tau, which have two disks in the field of view. For Haro 6-10, a simple Gaussian model of the emission from the other disk was subtracted before the analysis of each disk. For UZ Tau, the procedure was more elaborate. First a Gaussian model of the emission from the companion (UZ Tau W) was subtracted, and the remaining emission from UZ Tau E was analyzed. Then, the best-fit model of UZ Tau E was subtracted from the original data, and the emission from UZ Tau W analyzed separately.

All results are presented in Table 5. A comparison of the results obtained from independent data sets at similar wavelengths is shown in Table 6, which shows the excellent agreement of the constrained parameters (see also Fig. 4). In addition, the good agreement of geometric parameters with determinations from other studies is a further proof of the data quality (see Table 4, and Fig. 3).

Simple power law.

Results for the surface density parameters, Σ₁₀₀,p and R_out, are presented in Cols. 7 − 10 of Table 5. For most sources, the emission is largely optically thin, so the derived surface density will scales as roughly 1/T₀, but the outer radius remains essentially unaffected by the choice of the temperature. The only exceptions are the T Tau and Haro 6-10, which are essentially optically thick disks. FT Tau and Haro 6-13 may also be attributed to thick disks.

Exponential edge.

We generally used Eq. (5) to first locate the minimum. However, because of the direct dependency between the parameters, the errors on R_t and M_d were obtained by re-fitting the data using these parameters as primary parameters rather than R_c and Σ₀. Note that while the error on Σ_t may become very large, Σ₀ is generally constrained with a very similar accuracy as in Model 1.

Results are presented in Cols. 2−5 of Table 5. It was difficult to adjust this model to a few sources, among which were the apparently optically thick sources T Tau and Haro 6-10, and the single stars DQ Tau, DG Tau-b, and GM Aur. For the three latter stars, the best-fit power law has an index of p = 2. In this case, the expression in Eq. (5) attempts to fit γ = 2 and diverges. Finding the best fit requires the determination of the best transition radius R_c and its errorbar for all values of γ ranging from 0.6 to 1.9 (by steps of 0.1). The relative errors on R_t are generally larger than for R_c, because R_t/R_c diverges for γ → 2. No constraint on R_t is possible for DQ Tau. For GM Aur, only a lower limit is obtained, while for DG Tau-b, R_t is very marginally constrained: at the 2σ level, any value is acceptable. The error bars should be taken with care in those cases. A similar procedure was used for Haro 6-13 and Haro 6-33, for which R_t remains unconstrained at the 2σ level.

Column 6 of Table 5 indicates the difference in χ² between Model 1 and Model 2. A positive value indicates Model 2 (the viscous disk) provides an apparently better fit than the truncated power law. The significance of this result will be discussed in Sect. 5.2.

Deprojected, circularly averaged visibility profiles are displayed in the middle column of Fig. 5 for DM Tau and Figs. G.1 − G.22 for the others sources (in Appendix G). These deprojected visibilities only serve as an illustration of the fit results, but not to determine the parameter values and their errors.

4. Results

4.1. The dust temperature

For the assumed temperature law, our treatment differs quite significantly from those of Kitamura et al. (2002) and Andrews & Williams (2007), who assumed that the temperature derived from IR-emitting dust by fitting the SED also applies to the mm emitting dust. However, strong vertical temperature gradients are expected in disks (e.g. D’Alessio et al. 1999).

Because the mm emission comes from cold dust around the disk mid-plane, using a power law for the dust temperature distribution is an oversimplification. The dust temperature is expected to follow three different regimes, depending on whether the disk is optically thick or thin for absorption of the incident radiation and re-emission of its own radiation. The two extreme regimes predict temperature dependence, and are connected by a nearly constant temperature (or even slightly rising) region (“plateau”), whose extent depends on the source radial opacity profile (D’Alessio et al. 1999; Chiang & Goldreich 1997). A more self-consistent approach was taken by Isella et al. (2009), who derived dust opacities from the Mie theory assuming a specific dust composition and grain size distribution, and solve for the dust temperature in the two-layer approximation of Chiang & Goldreich (1997). However, in this case the derived dust temperature depends (by an unknown amount) on the assumed dust composition. Furthermore, using a single temperature for all grain sizes is an oversimplification. The dust thermal balance is largely dominated by the IR radiation (see Chiang & Goldreich 1997). Because the opacities are not gray, the temperature of dust grains is expected to depend on their size. The details will depend on the exact behavior of the dust emissivity as a function of wavelength, but generally larger grains are expected to be colder (Wolf 2003; Chapillon et al. 2008). Yet, these grains dominate the mm emission that we are observing.

Our approach of keeping the dust temperature as a parametric law allows us to directly measure the effective temperature of the emitting grains whenever the angular resolution is sufficient to resolve the optically thick core of the disk. Furthermore, we can estimate the impact of the temperature uncertainty on the derived surface density parameters. Such a step-by-step approach allows us to understand and quantify the existing couplings between the dust parameters, the disk temperature and the disk surface density.

Because the flux scales as T × Σ, the assumed values for the temperature may affect the derived shapes of Σ(r). In Model 1, the exponent p will be directly affected, because p + q is preserved for pure optically thin emission. This is confirmed by our analysis for both models (see Appendix C). However, the effects are small because our adopted value for q = 0.4 is a good first order approximation of most (reasonable) temperature profiles. In Model 2, an inappropriate temperature profile may affect R_c, because this parameter is constrained by the steepening of the emission as function of increasing radius. Again, Appendix C shows the effect is limited, R_c being affected by at most 20%.

In a few sources, Isella et al. (2009) derived dust temperature as a function of radius from a joint modeling of the SED and 1.3 mm images. We used the temperatures displayed in their Fig. 7 as an input in our modeling to check the magnitude of the effects in all sources we have in common. The results are presented in Table 9. The temperature law has no visible influence on the pivot radius, R_t, and affects γ by at most 0.1 − 0.2. Our used temperature laws are displayed on top of those of Isella et al. (2009) in Fig. 6. From Table 9 and Appendix C we conclude that the uncertainties in our assumed dust temperature distribution do not significantly affect the shape of the derived surface density distribution.

However, the disk masses are sensitive to the assumed dust temperature, since they scale to first order as 1/T. Furthermore, the dust emissivity index β_m can also be affected, because the contribution of the optically thick core depends on the dust temperature. The differences in the analysis of the MWC 480 performed by Hamidouche et al. (2006) and Piétu et al. (2007) illustrate the importance of the effect. A similar effect can be seen for DG Tau in Table 9: β_m changes by 0.5 between the two hypotheses on the temperature.

From the best-fit values, a few sources in our sample display partially resolved cores that may be interpreted as optically thick cores, and thus allow a direct determination of the temperature. As detailed in Appendix A, these “thick cores” satisfy two conditions: i) they have the same brightness at both wavelengths, and ii) their brightness distribution is relatively flat, because the temperature is expected to decrease as r ^{− 0.4 − 0.7} at most. The fitting process indicates that this happens for DG Tau, DG Tau-b, HL Tau, T Tau, and MWC 480. The derived values are presented in Table 8. Because the Model 1 and 2 have different opacity distributions (see Figs. G.1 − G.22), they predict different optically thick zones, and thus the temperature slightly depends on the assumed density model. For T Tau, the apparent difference is largely an artifact, because the source is basically a completely optically thick disk, for which the “viscous” disk model is poorly constrained. The measured values and slopes justify a posteriori our simple hypothesis for the temperature law. The dependence is small for DG Tau, DG Tau-b, and HL Tau, though. In the power law model, the extrapolated temperature at 100 AU for DG Tau is 17 K, close to our adopted value of 15 K for all other sources. HL Tau is slightly warmer, 19 K. For DG Tau-b, the temperature at 100 AU is 15 K, but the exponent is slightly lower than 0.4.

Formally, FT Tau has both a flat enough brightness distribution and a low apparent β_m to be consistent with optically thick dust, but would require a very low dust temperature to match the observed flux densities. A dust temperature of 10 K at 40 AU would just provide adequate flux (the brightness can be obtained from the (apparent) opacities displayed in Fig. G.9). Such a low value seems inconsistent with the relatively luminous and massive central star, so the warmer, optically thin solution with β_m ≃ 0 is to be preferred.

Among the observed sources, MWC 480 deserves specific comments concerning the temperature. In this bright source, the “thick core” is quite large, 50 − 80 AU. However, its brightness is moderate, which means that when this is interpreted as being optically thick, the derived temperatures are very low (see Piétu et al. 2006). The large size of the “thick core” results in substantial opacity corrections for β_m, which in turn leads to unrealistic values for Model 2.

An alternate explanation for the relatively flat brightness distribution in the inner part is a warmer, optically thin region with β ≃ 0. This is not consistent with the value of β_m derived from the integrated flux, and can only happen if β varies with radius (see Appendix A). This is studied in Sect. 4.5 and MWC 480 will be rediscussed in more detail in Sect. 4.6.

4.2. Surface densities

Isella et al. (2009) have published a high-resolution (0.7′′) survey at 1.3 mm of the Taurus region, with several sources in common to our study. It has been analyzed in terms of the viscous disk model, and Table 10 shows a comparison of the results. Note that in this analysis, we assumed no inner hole for Lk Ca 15 to provide a consistent comparison, and its apparent deficit of emission in the center is purely explained by a negative value for γ. In general, our data have a higher resolution and are slightly more sensitive, which results in error bars that are lower than in Isella et al. (2009), the only exception is GM Aur, for which our resolution is moderate.

The agreement between both studies is reasonable, typically within 2σ. The most notable exception is DG Tau. DG Tau was further studied at higher resolution by Isella et al. (2010); the agreement on R_t is reasonable, but they find γ = 0.28 ± 0.05 instead of γ = 1.6 ± 0.1 in our study. The difference between the two results may be due to the widely different uv coverage, linked to a non symmetric source. Our data are dominated by fairly moderate baseline lengths (up to 300kλ), while Isella et al. (2010) find a substantial contribution to the imaginary part of the visibilities at 1.3 mm up to 200kλ (see their Fig. 2 and image in Fig. 10).

We also note that the agreement is better on R_c (or R_t) than on γ. This is to be expected, as R_c is a first order parameter (the radius which encloses 63% of the disk mass), while γ is a second-order parameter (the slope of the surface density distribution).

4.3. Emissivity index

β_m values have been reported for a number of sources in our sample by Rodmann et al. (2006) and Ricci et al. (2010). Their analysis is different from ours, because β_m is derived from spatially unresolved multi-wavelength data, from a fit of the SED. Rodmann et al. (2006) use a simple power law to derive the spectral index α of the  mm SED between 7 and 1 mm. Overall, the agreement with our results is poor, most likely as a result of several effects. First, Rodmann et al. (2006) apply a uniform correction for opacity, while we have shown that the existence of optically thick cores affect β_m very inhomogeneously, with corrections ranging from 0 to 0.5. Second, the different frequency span must also affect β_m, because using a power law for the dust emissivity is only an approximation; in particular, the emissivity is expected to steepen at long wavelengths (e.g. Draine 2006). The agreement with the results of Ricci et al. (2010) is much better, most likely because they use a more elaborate procedure for the SED fit, in which some estimate of the disk size and surface density slope is used to account for the optical depths effects.

4.4. Individual objects

4.4.1. Multiple stars

Haro 6-10 stands out as exceptional. Although the formal fit gives marginally optically thin disks and β_m ≃ 0, this is likely to be an artifact caused by seeing limitation. Indeed, any small “seeing” effect spreads out a little emission and makes the source slightly more extended than in reality. This mimics an (optically thin) halo. Thus, Haro 6-10 is best represented by (two) optically thick disks of radii around 15 AU (scaling as 1/ since only the total flux is constrained). This result indicates that the amplitude and phase calibrations are sufficiently accurate to determine sizes as small as 30 AU (total), or about 1/5th of the synthesized beam in this case. The inclination cannot be derived for Haro 6-10. The minimum mass of each disk is 10^-3   M_⊙(see Appendix F).

T Tau was already studied by Hogerheijde et al. (1997) and Akeson et al. (1998) in the mm domain. As in these studies, only the northern member of the multiple system is detected. Like Haro 6-10, the emission can be explained by a nearly optically thick disk. Because of the larger size, the seeing effect is negligible and only the optically thick solution is found to be viable. Our best-fit inclination of 40    ±    4° is somewhat larger than the  ~ 20° derived by Ratzka et al. (2009) from IR studies. However, this only influences the apparent opacities by the ratio of the cos(i), i.e. about 20%. The minimum mass of the disk is 0.007   M_⊙, assuming the disk is optically thick.

The quadruple system UZ Tau shows emission from two regions: one around the spectroscopic binary UZ Tau East, the other near the optically resolved wider binary UZ Tau West (separation 0.34′′ at PA ≃ 0, Simon et al. 1992). Given the disk inclination of UZ Tau East (Simon et al. 2000) which is confirmed by our new measurements, and assuming disks and orbits are coplanar, the true deprojected separation would be  ~ 100 AU. Interpreting the emission around UZ Tau W as a single disk yields a similar orientation (consistent with coplanar disks) and an outer radius of 120 ± 45 AU. This is fairly large compared to the binary separation, and may be difficult to reconcile with tidal truncation. This result, however, could be an artefact of improper subtraction of the UZ Tau East emission because any small (positive) residual emission left around UZ Tau East could bias the derivation of the position angle and size. A solution with two circumstellar disks is not totally excluded by our data. From the images, we find that the emission centroid is in between UZ Tau West A and B (see Fig. 7). The displacement observed between 1.3 mm and 2.7 mm suggests that the disk around West B disk is more optically thick than that around West A. Under the interpretation of circumstellar disks, their minimum mass is 6 × 10^-4   M_⊙.

The small size of circumstellar disks in known binaries suggests that tidal effects are responsible for their truncation, although a firm conclusion cannot be drawn because the inclination of Haro 6-10 is unknown.

Mathieu et al. (1997) found DQ Tau to be a non-eclipsing, double-lined spectroscopic binary, comprised of two relatively equal-mass stars M ≈ 0.65   M_⊙with spectral types in the range of K7 to M1 and an orbital period of 15.804 days. The orbit is eccentric, but with an apastron around 0.28 AU, the tidal cavity should be much smaller than 1 AU. DQ Tau has been recognized as variable in the mm domain by Salter et al. (2008). The variability is caused by interactions between the magnetospheres when the two stars are near periastron, so that flares happen periodically. The observation dates and derived total flux for each date are given in Table 11. No evidence for variability is found in our data as expected, since none of our observations happened close to periastron. The measured emission is thus coming purely from the dusty (circumbinary) disk.

Another star is possibly affected by binarity: FT Tau, which displays a weak, but significant (6σ) emission 1.3′′ west of the main star (see Fig. 8), and a very small ( ≃ 60 AU radius) disk with β_m ≃ 0 (see Table 5 and Fig. G.9). The position of the secondary peak of mm emission is, however, different from that of the near IR source found by Itoh et al. (2008).

The case of HH 30 is somehow unusual. Anglada et al. (2007) suggested that HH 30 is a binary based on the precession of its optical jet, but could not decide between a close binary and a  ≃ 15 AU separation. Guilloteau et al. (2008) showed that the deficit of mm emission could be interpreted as a central hole consistent with the tidal truncation in the wide binary model. Here, in Model 2, the inner radius becomes unsignificant: any value below about 45 AU is acceptable for R_int, because of the very steep decrease of the surface density profile for this high negative value of γ ≃  − 2. In essence, this means γ is constrained by the apparent sharp decrease of the emission near 120 AU, and not by the central deficit.

4.4.2. Sources with holes

For DM Tau, modeling the near and mid-IR SED (Calvet et al. 2005) indicates an inner hole of about 3 AU. Although this small hole is below the detectability limit of our observations, we used it in our analysis.

A deficit of emission at the center of the disk of LkCa 15 was discovered by Piétu et al. (2006), who interpreted it as a 45 AU radius hole. This central dip was also observed at lower resolution by Isella et al. (2009), but they suggested that it could be due to a negative value of γ. Our higher angular resolution data allow us to test which hypothesis best represents the observed brightness distribution. Results are reported in Table 12. The no-hole hypothesis is rejected at the 3σ level, and the best fit is obtained with an inner hole of 38 ± 4 AU. The near-IR imaging of Thalmann et al. (2010) confirms the sharp nature of the rim of the inner hole and indicates a radius of 46 AU. The transition radius R_t remains relatively unaffected by the presence or absence of a hole, but the value found for γ strongly depends on the hole size: the best-fit solution is compatible with γ = 0.

For GM Aur, the lack of 10 μm emission suggested a central hole of R_int = 25 AU (Calvet et al. 2005). The hole has also been detected in the gas traced by CO, through spectro-imaging of the J = 2−1 transitions of the ¹²CO, ¹³CO, and C¹⁸O isotopologues indicating very low gas surface densities in these regions: Dutrey et al. (2008) indicate a size of R_int = 19 ± 4 AU. This size has been confirmed by direct imaging of the dust emission at mm wavelengths (Hughes et al. 2009). Like for DM Tau, we thus assumed R_int = 20 AU. The strong dependence of γ upon the possible existence of a central hole also exists for GM Aur. Indeed, assuming no hole, we recover a very similar solution to that found by Isella et al. (2009, see Table 10), although it is somewhat worse (near the 2σ level) than our nominal solution obtained for R_int = 20 AU.

In conclusion, with the exception of HH 30 which was discussed in Sect. 4.4.1, allowing for central holes offer better solutions, and brings the surface density exponent γ back to “standard” values between 0 and 1.5.

4.4.3. Young sources

HL Tau is a Class II object, for which the central star is not directly visible. Our measured position angle is consistent with that of the jet and of previous high-resolution studies of the mm and centimeter emission from this region (Looney et al. 2000; Anglada et al. 2007; Carrasco-González et al. 2009). The inclination of the source is more debated: early work from Cohen (1983) assumed a nearly edge-on disk, while i = 56 ± 10° can be derived from the 7 mm deconvolved size from Wilner et al. (1996). Our result better agrees with the submm data obtained by Lay et al. (1997), 42    ±    5°, and is also consistent with the obscuration of the redshifted jet (Pyo et al. 2006). At the observed scale, the envelope that dominates the submm flux is filtered out (Looney et al. 2000). Our major finding is the substantial opacity at mm wavelengths in the inner 40 AU, which allows us to constrain the temperature, but this significant opacity does not prevent structures from becoming visible at longer wavelengths, 1.3 cm or 7 mm. Our angular resolution is insufficient to separate the possible enhancement reported near 65 AU at 1.4 mm by Welch et al. (2004) and 1.3 cm by Greaves et al. (2008), but not confirmed at 7 mm by Carrasco-González et al. (2009).

DG Tau is a bright embedded star driving an optical jet at PA 226° (Eislöffel & Mundt 1998). It is surrounded by a large ¹³CO disk orthogonal to the jet (Sargent & Beckwith 1989; Kitamura et al. 1996), whose kinematics indicate a stellar mass around 0.7 M_⊙ (Testi et al. 2002). Inclinations of 45° and 38° are found by Pyo et al. (2003) and Eislöffel & Mundt (1998), respectively. Our measured inclination of 32 ± 3° is in favor of lower values. The results quoted in Table 5 only slightly depend on the assumed orientation and inclination: γ can be decreased by 0.1 and R_t increased to 19 AU for the best fit orientation. The higher resolution data of Isella et al. (2010) also give lower inclinations and small (22 AU) R_t, but with a very different value for γ (see Sect. 4.2).

DG Tau-b is a young, totally obscured, star at the apex of a wide angle cavity seen in scattered light (Padgett et al. 1999). It drives an optical jet and a molecular outflow (Mitchell et al. 1997). Although the position angles derived from the jet and disk agree, we find the disk inclination to be only 64    ±    1°, while the jet inclination is estimated to be higher than  >75° from proper motion measurements (Eislöffel & Mundt 1998). We also note that the DG Tau-b disk is best fitted with a higher flaring index h than assumed for the other objects of our sample. We used h = 1.35, for which a scale height H₁₀₀ = 27 ± 8 is required to reproduce the observed continuum emission. The high flaring index is consistent with the fairly flat temperature distribution (q = 0.3 ± 0.1) also found in this source.

4.5. Radial dependency of the dust properties

Most previous studies assumed that the dust properties were uniform across the disk. The dual-frequency resolved images allow us to test the validity of this hypothesis, and eventually constrain the variations of dust properties as a function of radius.

4.5.1. Emissivity Index β

In Table 7 smaller transition radii R_t are found from 2.7 mm data than from 1.3 mm data for four out of eight sources: CI Tau, CY Tau, DL Tau and DM Tau. For the other sources, the combination of sensitivity and resolution at 2.7 mm data is insufficient to distinguish. Equivalently, in the truncated power law analysis (Table 7) the slope of the surface density p is systematically steeper at 2.7 mm than at 1.3 mm. A similar result was recently obtained for CQ Tau by Banzatti et al. (2011).

A possible cause for this effect is contamination by free-free emission, which adds a point-like source at lower frequencies. However, none of these sources have sufficient free-free emission to significantly contaminate the 2.7 mm flux (see Rodmann et al. 2006, for the measurements). From Rodmann et al. (2006), the contamination does not exceed 3% near 2.7 mm. Removing a point source of this intensity from our 2.7 mm data does not affect our results.

Thus, the different solutions found at the two wavelengths indicate a change of dust properties, at least in the spectral index of the emissivity β, with radius. The larger p and smaller R_t at 2.7 mm than at 1.3 mm imply that the ratio of T_b(1.3 mm)/T_b(2.7 mm) increases with radius, hence β increases with radius. The apparent β(r) as a function of radius can be derived from (13)where β₀ is the constant value used to derive the apparent surface densities Σ_a(r) and Σ_b(r) at both wavelengths, i.e. β₀ = β_m (see also Isella et al. 2010).

The increase of β(r) with radius is most easily understood in the framework of the truncated power law analysis, because it simply turns into a logarithmic dependence of β(r) as a function of radius Δp = p(2.7   mm) − p(1.3   mm) is systematically positive in our sample (see Table 7). However, the apparent significance level is low for each source, as Δp apparently exceeds its 2σ uncertainty in only two sources (CY Tau and DM Tau). Better constraints can be obtained by fitting the logarithmic dependence of β(r) directly to the data (16)The values of β_r are reported in Col. 2 of Table 13 (for Model 1, but similar values are obtained for Model 2). It becomes now clear that the radial dependence is highly significant, because the weighted mean value is β_r = 0.34 ± 0.04 (ignoring FT Tau, which has a negative β everywhere). A Student’s T-test applied to the distribution of values of β_r reported in Table 13 (including FT Tau) indicates less than 2% chances of being compatible with β_r = 0.

For the softened-edge model, the β(r) function implied by Eq. (13) is more complex, and an illustration of the shape of this function is given in Fig. 10, which displays this apparent β(r) for two of the sources, CI Tau and DL Tau. The hatched areas indicate the approximate range of allowed values, obtained by adding and subtracting 1σ to each of the parameters defining the opacity function at the two wavelengths (R_c, γ and Σ₀ from Eq. (5)). Because some of these parameters are actually correlated, this is only an estimate of the error on the profile. The apparent index β is large (>1.7) in the outer disk parts (r > 150 − 250 AU), while it is smaller than about 0.6 near 50 AU.

The shape of the radial dependence of β in Fig. 10, and the logarithmic dependence in Eq. (16), are simple results of the choice of shape of the surface emissivity distribution, and have no physical constraints attached. In particular, apparent values of β below 0 or above 1.7 can result from such an analysis.

Because of the limits in angular resolution and sensitivity, some prescription of the evolution of the dust properties as a function of radius, assuming realistic conditions, must be specified to obtain better insights on the dust properties versus radius. A poor choice could make the radial dependence apparently non significant. To illustrate the problem, we used (17)which varies between 0 (large grains) and 1.7 (small ISM-like grains). With this functional, we obtain significantly better fits, at least by 3σ, but up to 8σ in DL Tau (see Table 13). Furthermore, the improvement does not depend on the assumed shape of the surface density: power laws or tapered edges yield identical results for the pivot R_b and width R_w, although the errorbars on these parameters are typically 30% lower in the power law hypothesis. Note that there is a fairly strong correlation between R_b and R_w, and their errorbars are in general not symmetric. To better illustrate the variations of β(r), the resulting range of allowed values for β(r) for each source is given in Fig. 11. The logarithmic dependence found from Eq. (16) is also indicated. Both functionals give approximately the same values in the regions where β(r) is actually constrained, that is from 30 AU to the R_out of the power law. However, the log dependence fails to characterize the sharpness of the transition from low to high values of β.

Finally, although our analysis of β excludes the flux calibration uncertainty, it is worth emphasizing that this does not affect the radial variations of β(r), but only the mean value β_m. It also does not affect the relative differences in β_m between sources, because all observations were made in an homogeneous way, with all spectral index measurements based on an assumed index of 0.6 for MWC 349.

4.5.2. Absorption coefficient κ_ν1.3   mm

If grains vary in size with radius, the absorption coefficient κ(ν,r) will also vary. The surface density laws that were derived so far were derived assuming κ(ν,r) = κ(ν₀)(ν/ν₀)^β(r) with ν₀ = 230 GHz. In practice, there is no physical justification for any value for ν₀, because for essentially all models of grain growth the absorption coefficient and the apparent emissivity index vary simultaneously in a more complex way. One can attempt to use a more physical approach to the grain properties, using dust absorption coefficients derived from a physical model (e.g. Draine 2006, and references therein). For example, Isella et al. (2009; see also Natta et al. 2004) derived the absorption coefficient from a fixed grain composition, with a size distribution controlled by a single variable parameter. The size distribution is a power law with a fixed minimum and maximum radius and an exponent g. The absorption coefficient κ and apparent emissivity index β at 1.3 mm are plotted as a function of g in Fig. 9. From this dust model, we can derive a function κ(β), which can be used in our model with the same assumptions about the radial dependency of β(r) as previously done.

With the prescription of the opacity law described by Fig. 9, and β(r) as in Eq. (17), the pivot radii R_t/R_c are not changed very significantly. The largest changes are for CY Tau and UZ Tau, where R_c decreases by 50%, DM Tau, where it increases by 50%, and MWC 480. Effects on γ are negligible except for MWC 480 and UZ Tau (see Table 13). The relatively small effect on R_c, and γ is explained because β(R_c) is close to 1 in most of the sources studied, and for this value κ_ν(230   GHz) has an extremum. Thus the variations of κ_ν(r) around R_c are relatively moderate, and accordingly the shape of the derived surface density is mildly affected by the radial variations of κ_ν(r).

However, these small apparent changes may be misleading, because they implicitly depend on the assumed shape of the surface density law. As β(r) is getting close to 0 in the disk center, and thus the absorption coefficient κ(1.3   mm) could be much smaller at small radii, it is also completely possible to have a much steeper surface density gradient in the inner 40 AU. This remains hidden from our study because of the angular resolution, but also because the inner 20 − 40 AU become optically thick in some sources. If steep gradients like this exist, longer wavelength images should be able to reveal them. The strong changes observed in R_c,γ for MWC 480 and UZ Tau are also manifestations of this effect, although at larger scales.

In our sample, only HL Tau was studied with sufficiently high angular resolution at 7 mm and 1.3 cm to confront images with the above prediction. Although surrounded by a diffuse halo, the 7 mm image of Carrasco-González et al. (2009) is indeed very centrally peaked, but a quantitative comparison with our results is not directly possible because it displays complex structures.

4.6. MWC 480 revisited

In the simple β(r) = β_m hypothesis, the disk of MWC 480 appears sufficiently optically thick at 230 GHz to allow the derivation of the dust temperature (Piétu et al. 2006). The optically thick region is even large enough to constrain the exponent q to some extent. Leaving both T₀ and q as free parameters, Model 1 and Model 2 give different best fits for the temperature, because in the best fit for Model 2, the radius at τ(1.3   mm) = 1 is much larger (80 AU) than for Model 1 (35 AU, see Fig. G.19). In addition, β_m is larger by about 0.3 in Model 2 than in Model 1, because the optically thick core is much larger. Furthermore, since the extrapolated temperatures in the best Model 2 are very low (7 K at 100 AU, and 2.7 K at 400 AU), the emission is no longer in the Rayleigh-Jeans domain, and β_m increases because the corrections are larger at 1.3 mm than at 2.7 mm. In practice, Model 2 finds a low temperature with a steep exponent ( ≃ 0.6) because of two effects: i) the brightness is identical at 2.7 mm and 1.3 mm in the inner 40 AU, and ii) the imposed shape of the surface density is too flat in the inner 80 − 100 AU (in order to provide sufficient opacities beyond 100 AU). To account for these two constraints, an optically thick core of 80 AU is fitted, with a steeply decreasing temperature. High temperatures can only be found by allowing the surface density to fall faster than Model 2 allows between 40 and 80 AU.

Clearly, in this case, although low temperatures are needed in the inner regions, extrapolating the same power law introduces non-physical biases on the disk mass and on β_m. Leveling the temperature to a minimum value of 12 K beyond 45 AU provides a better fit to the observations, and allows us to bring back β below 2. This may be an indication of the temperature rise with radius that is expected to begin when the opacity for re-emission drops below 1. However, the very low apparent temperatures in this object are surprising because of the luminous central star. This may be linked to the geometry of that source. From its IR SED, MWC 480 is a Group II Ae disk, which is interpreted as a self-shadowed disk with small flaring (Meeus et al. 2001). Indeed, it has never been detected in scattered light, despite a fairly favorable inclination. Yet, the temperature derived from ¹³CO is  ≃ 23 K at 100 AU, with an exponent q = 0.4 ± 0.1 (Piétu et al. 2007), and if the disk remains optically thick even at 3 mm, we would expect dust and gas to be thermalized at the same temperature.

Allowing β(r) to change with radius also offers a much more attractive solution to the continuum emission of MWC 480. The flattening of the emission in the inner 50 − 80 AU is no longer ascribed to an optically thick core at low temperatures, but to a flattening of the surface density distribution, while the ratio of 2.7 to 1.3 mm emission is matched by allowing β(r) to become small in the inner 30 AU. Although it is equivalent in χ² to the constant β(r) = β_m solution, this new model agrees with less extreme dust temperatures. In fact, the dust emission is largely optically thin in this case, and there is a substantial degeneracy between the dust temperature and the derived disk mass/surface density. A lower limit to the dust temperature at 100 AU is 23 K (assuming q = 0.4), which is consistent with the temperature derived from ¹³CO line emission by Piétu et al. (2007). This lower limit was used to derive the surface density.

If we use κ(β) as implied in Fig. 9, the fit quality is slightly degraded, but most importantly, the derived shape for the surface density and the temperature profile are significantly affected (see Table 13). We find γ ≃ 1.5, and a large transition radius R_t > 70 AU, much like for GM Aur. The best-fit temperature profile is flat, q = 0.0 ± 0.1, with T > 25 K. The higher χ² value derived under these assumptions may be related to an oversimplified temperature profile, as in the simpler analysis q ≃ 0.5 was found in the inner regions.

5. Discussion

5.1. Dust properties

From the above results, we find large grains (β < 0.5) in the inner 60 to 100 AU, and small grains beyond for seven sources (CI Tau, CY Tau, DL Tau, DM Tau, DG Tau-b, MWC,480 and UZ Tau-E). Two other sources in our sample have very low β: the FT Tau disk is truncated at 60 AU, while DQ Tau has not been observed with sufficient resolution at 2.7 mm, so that β(r) is not constrained in the outer regions. A third source may have low β up to 60 AU: T Tau N, although we interpreted it as being optically thick.

On the other hand, β is not constrained in the inner region for the two other sources observed with sufficiently high resolution in our sample, because in DG Tau, the inner 50 AU may be optically thick, while for LkCa 15, the inner 50 AU are (largely) devoid of dust. For HH 30 we find a low β below 120 AU, while it is known from the scattered light images that small grains exist at least up to 250 AU or even 430 AU (Burrows et al. 1996), the outer radius of the gas distribution (Pety et al. 2006). Finally, in HL Tau, large grains exist in the inner 20 AU, as shown by the 1.3 cm and 7 mm images (Carrasco-González et al. 2009). Thus, in essence, all sources in our high-resolution sample show large grains (low β) below 100 AU and small grains beyond, although the detailed shape of the radial dependence cannot be characterized by our data.

The apparent variations of p with wavelength observed by (Banzatti et al. 2011) for CQ Tau also points out towards an increase of β(r) with radius in that source. Moreover, although they considered it to be insignificant, the same trend is found in RY Tau and DG Tau by Isella et al. (2010). Thus, the radial dependency of β(r) appears to be a general property of disks.

Our findings that β is low in the inner 60 to 100 AU of all disks in which we can constrain the radial dependency also sheds new light on the results quoted by Ricci et al. (2010). Ricci et al. (2010) found a lower average value for the spectral index α for disks with low 1.3 mm flux than in disks that show strong emission. A possible interpretation is that these weaker disks are optically thick and very small, like those surrounding the binary Haro 6-10. These weak disks may just miss the extended, low brightness parts with high values of β that we found in bright sources. In our sample, a clear example for this behaviour is FT Tau. Given our measured R_b, a testable prediction is that these faint disks should be smaller than about 100 AU in radius. Note that this does not address the origin of the small size for these disks, although tidal truncation is an obvious candidate. On the other hand, in AB Aur, which has an inner hole around 100 AU (Piétu et al. 2005), the mm emission is coming from the small grain regions, which results in a mean β_m = 1.4 ± 0.2, which is different from all other sources. Such a high β_m is not an indication of different grain growth in this source, but just a side effect of the radial dust distribution. We further stress that the β_m values derived in all previous analyses represent an ill-defined average over the disk structure.

5.2. The shape of the surface density distribution

Given the high resolution and sensitivity, can we decide which model fits the data better? The lowest χ² is the usual indicator, but care must be taken that the χ² is not affected by different biases between the two models owing to numerical effects in the model computation. The precision required for this is always higher than the precision required to obtained converged parameters and errors within a given model, because the discretization effects impact models differently (see Appendix B). For the models considered, the problem is somewhat relaxed because they both derive from a generic one (see Eq. (5)). We nevertheless checked by using oversampled grids that the χ² results were converged.

From Table 5, the softened-edge model does not appear superior to the power law model to represent the observations. In this process, the compact optically thick sources should be ignored. For these sources, the data are insensitive to the true shape, but can be significantly affected by small instrumental effects. For example, the seeing that results in flux spreading because of atmospheric phase variations tends to produce a small halo around the compact core. In our sample of 23 individual objects, this may affect five sources. Of the remaining objects, four sources are best represented by a power law: DG Tau, DQ Tau, HL Tau and (marginally) DM Tau. On the other hand, six sources are better fitted by the exponential-edge model: CI Tau, CY Tau, DL Tau, UZ Tau E, and marginally LkCa 15. Both models fit equally well the last seven sources, which were observed with lower angular resolutions except for DG Tau b.

Despite the high resolution (projected baselines above 500kλ) and sensitivity, the shape of the surface density remains difficult to constrain. The truncated power law was initially used because it provides the simplest parametric model. It is furthermore not linked to any specific physical disk model, a property which can be seen either as an advantage (by providing no specific bias) or handicap (as having no physical ground). Its principal failure was its inability to represent continuum and spectral line emission with the same outer radii (Piétu et al. 2007). The softened-edge model has recently gained favor because, as suggested by Hughes et al. (2008), it may provide a framework that can explain both the continuum and optically thick CO emission. The exponential taper is often referred to as having a physical background, because viscosity is expected to spread out initially small disks. However, the exponential taper is only a specific solution of self-similar evolution of a viscous disk with a power law distribution of the viscosity (with a constant in time exponent). In practice, self-similarity and time independence are unlikely to strictly apply to real disks, so the resulting specific shape is also an approximation. Any core+halo structure would essentially yield the same result, provided the halo is just dense enough to explain the molecular emission, but tenuous enough to have little continuum emission from dust. This core+halo structure was invoked by Dutrey et al. (1994) and Guilloteau et al. (1999) to interpret the circumbinary environment of GG Tau.

In the strict framework of a viscous disk model, we find values of γ that are somewhat larger (and with higher dispersion) than those derived in the previous studies by Isella et al. (2009) and Andrews et al. (2009). The discovery of a radial dependence of the dust properties brings additional complexity to the problem. Clearly, the surface density of the gas is not well traced by the continuum emission at a specific frequency in this picture. It all depends on how the dust emissivity κ(ν,r) changes as a function of radius, so that the derived γ is expected to depend on the assumed dust properties.

It is also important to realize that the derived dust masses of the disk may be significantly affected by the variations of the dust properties with radius. Table 14 indicates the disk masses obtained for Model 2 using Eq. (17) with the two different hypotheses on the dust absorption coefficient (κ_ν(1.3 mm) constant or tied to β(r) as from Fig. 9). In our sample, although the effect is small for the other sources (about 20%), the masses of the CY Tau and DM Tau disks are strongly modified when using the dust properties from Isella et al. (2009). In particular, the mass of the DM Tau disk becomes quite significant (0.2 M_⊙) compared to the stellar mass (0.5 M_⊙Dartois et al. 2003). Such a large mass would have significant effect on the rotation curve of the gas, while it is known to be Keplerian with high accuracy (velocity exponent 0.50 ± 0.01, Piétu et al. 2007).

This result is under the assumption of a “normal” gas-to-dust ratio of 100. However, the gas-to-dust ratio itself is expected to change as a function of time and radius in the disk. For DM Tau, the potentially large dust mass suggests that the gas-to-dust ratio must be decreased. Molecular tracers may help to constrain the gas surface density more directly, but then a good understanding of the chemistry is required to recover the hydrogen content from the (very few) trace molecules that display strong enough lines to be observable: CO and its isotopologue ¹³CO, HCO⁺, CN, HCN, CS and H₂CO (Dutrey et al. 1997).

A simpler alternative is that the adopted dust properties are inappropriate. For example, with similar grain size distributions, but using a different dust composition (in particular porous grains), Ricci et al. (2010) derive dust opacities on the order of 3 − 20 cm² g^-1 instead of 0.4−2 cm² g^-1 from Isella et al. (2009). The overall dependencies of κ and β upon the grain size distribution display the same characteristic behavior. Note, however, that it is possible to obtain β values above 2, provided the grain size distribution as a relatively steep cutoff near a₊ = λ/2π, i.e. about 0.5 mm for λ = 3 mm, because the emissivity of a grain size a has a pronounced maximum at wavelengths  ~2πa (see e.g. Natta et al. 2004), before dropping as 1/a at longer wavelengths. In the following, we scale down the surface densities of sources analyzed with radial dust opacity gradients by a factor 3, to avoid a different bias in the comparison with sources for which this analysis was not possible.

Figure 12 displays the resulting surface densities (of gas+dust) for the sources in the sample. Uncertainties were omitted for clarity in this figure, but can be recovered from Figs. G.1−G.22. We note that the younger sources have higher surface densities in the inner 50 AU than other objects. They are also more centrally peaked, on average. This picture is qualitatively similar to the predictions from viscous disk evolution. Figure 12 also displays the profiles derived for the Solar Nebula, the MMSN (Hayashi 1981, gray line) and the solution proposed by Desch (2007), which accounts for the early planet migration as proposed by the Nice model (Tsiganis et al. 2005; Gomes et al. 2005) (yellow range). The solution proposed by Desch (2007) (see Appendix E) for the Solar Nebula is a steady-state solution, which allows for sufficient time (few Myr) for the giant planets to reach isolation mass. In comparing with our results, it is important to realize that our observations constrain essentially the slope and surface density between 50 to 150 AU, while the other regions are obtained by extrapolation of the analytically prescribed shape. In our sample, only the youngest objects have sufficiently high surface densities to be compatible with the MMSN.

5.3. Towards an evolutionary model?

5.3.1. Viscous evolution of the gas disk

Figure 13 displays the characteristic radius R_c as a function of estimated stellar ages. The figure apparently suggests an increase of R_c with age. Performing a Spearman rank-order correlation test indicates a correlation coefficient of 0.60, and a small probability of random distribution (0.7% only). This correlation study does not include the error bars on age and R_c, however. Furthermore, the correlation coefficient is heavily influenced by the two youngest objects, DG Tau and HL Tau and the two oldest ones, GM Aur and DM Tau, all sources for which the power law model gives a better fit than the softened-edge model. Nevertheless, taken at face value, our data seem to confirm the trend suggested by Isella et al. (2009), which they have interpreted as evidence for the viscous evolution of disks.

In the framework of self-similar viscous evolution, the surface density depends on five intrinsic parameters: a normalization constant C, the initial disk radius R₁, a normalized age T, the viscosity ν₁ at radius R₁ and its radial exponent γ. We have three measurements from our study (M_d   or   Σ₀,R_c,γ), the stellar age t_∗ from evolutionary tracks as quoted in Table 1, and the mass accretion Ṁ, usually derived from the accretion luminosity (Gullbring et al. 1998). Appendix D details the relationship between these observable quantities and the primary parameters of the radial surface density evolution.

A perfect correlation between disk sizes and age is not expected. The initial characteristic sizes of disks will add significant scatter. In this respect, the most significant fact is perhaps the envelope of allowed R_c versus ages, which places an upper limit on these initial sizes. This limit is related to the initial specific angular momentum. Larger disks would fragment and lead to binary and/or multiple systems. In this respect, it may be relevant that the young object with the largest R_c is UZ Tau E, a spectroscopic binary member of a hierarchical quadruple system. Another source of scatter resides in the exponent of the viscosity γ. Figure 13 displays evolutionary curves of R_c = R₁T^{1/(2 − γ}) for three values of γ, starting with a common initial radius R₁ = 10 AU (see Appendix Eq. (D.4) for a derivation of the evolution of R_c versus time). Although the extreme values of γ = 0.5 and 1.5 appear, at first glance, to provide a good fit to the envelope of the distribution of R_c versus ages, the actual picture is more complex. In particular, a number of stars close to the γ = 0.5 curve have in fact γ = 1.5 from our data set, while the reverse is also true.

The viscous timescale is given by (see Eq. (D.11))

Unfortunately, because the stellar ages t_∗ are very uncertain and we expect in general t_s ≪ t_∗, t_s remains largely unconstrained by the observations. Rather, Eq. (D.11) provides a loose constraint on the allowed range of disk masses and ages. An alternate way to constrain the viscous timescale is to look at the younger objects, for which the viscous evolution may not have had time to erase the initial conditions. In our sample, younger objects are better represented by power laws. This suggests that young disks are still influenced by the history of infall from the original proto-stellar cloud, and that the viscous timescale is on the order of a few 10⁵ years, the age of these youngest objects.

With this rough estimate for t_s, we can in principle derive T and recover the distribution of R₁ in our sample from Eq. (D.4), but the propagation of errors leads to large uncertainties (see Fig. 14). This is to be expected, because the viscous evolution has largely erased the memory of the initial conditions.

Self-similarity would also imply that the exponent γ remains constant over age. The distribution of γ vs. age is given in Fig. 15. Our distributions of γ are somewhat different from those derived by Andrews et al. (2009) and Isella et al. (2009). The former is centrally peaked around 0.9. The latter exhibits values lower than 0.8; however, we have argued in Sect. 4.4.2 that some of the derived values are affected by the interpretation of the central deficit of emission. In our sample, although there is no obvious correlation, stars of ages 1 − 3 Myr have on average lower γ ( ≃ 0.3) than either younger or older objects. Note that from Eq. (D.7), for γ = 1.5, we expect Ṁ(t) ∝ t^-2 in good overall agreement with the empirical relation found by Hartmann et al. (1998). On the other hand, γ = 0.5, which corresponds to the so-called β prescription of the turbulence (see Appendix E for details), yields Ṁ(t) ∝ t ^− 4/3 only, somewhat too small to explain Hartmann et al. (1998) correlation. Note however that large values of γ are unlikely to apply to the whole lifetime of the disks: if we assume γ has been constant with time, the two old disks with large γ (GM Aur and MWC 480) would have started with exceptionally small radii (<2 AU, see Fig. 14). Thus, invoking some evolution of the viscosity exponent with age seems required.

An alternate vision on the viscosity is to look at its value at some arbitrary fixed radius, for example at R₁₀₀ = 100 AU. Using the α prescription of the viscosity, the α parameter at 100 AU is given by (Eq. (D.17)) (18)where c_sH scales as to first order (see Appendix D for the derivation). Using from Table 1, and our adopted values of T_g = 15 K and H = 16 AU at 100 AU for the median value of L/M² = 3 L_⊙/M_⊙², the resulting α are displayed in Fig. 16. There is substantial scatter, but the measurements suggest an overall decrease of α versus time, roughly as 1/t_∗.

5.3.2. Evolution of the dust

The radial dependence of β(r) and the behavior of R_c as a function age may be understood in a more complex scheme where viscous spreading plays a significant role. Indeed, only small dust grains are efficiently coupled to the gas, while the larger ones should drift quickly inward (e.g. Weidenschilling 1977). Hence, one naturally expects that large grains will remain confined to the inner regions, which leads to an apparent increase of β(r) with radius. Simulations of the grain-growth, dust-gas coupling and fragmentation processes have been performed by Brauer et al. (2008), and further expanded by Birnstiel et al. (2010a) to include the disk accretion phase and viscous evolution. There is no specific prediction for the evolution of the shape of the grain size distribution with radius which could be compared to our data. However, from Fig. 10 of Birnstiel et al. (2010a), the smaller grains have outward net velocities beyond about 80 AU. A similar result was found by Garaud (2007), although her approach neglects the fragmentation processes. This radius is similar to the transition radius between low and high values of the emissivity exponent β found in our study.

Despite being a rather ill-defined quantity, the average β_m has been used to characterize disks in most previous studies. Figure 17 displays β_m as a function of estimated stellar ages. Very young sources have high values of β_m, comparable to the value found for ISM grains, which could indicate that the dust grains have not yet significantly evolved in these objects, at least at the characteristic distances that we sample in these sources (100 − 300 AU). Typical T Tauri disks have β_m < 0.7, which indicates significant grain growth. However, we also find that the older disks display high values of β_m too, well above the characteristic value for the “middle-aged” T Tauri stars. The radial dependency of β provides an explanation for this distribution of average β_m with stellar ages. As disks get larger with time, the apparent average β_m increases, which leads to the secondary increase of β_m for older objects as shown in Fig. 17.

Birnstiel et al. (2010b) evaluate the effect of the grain growth and fragmentation on the apparent spectral index α_{1−3   mm} for disk masses ranging from 0.005 to 0.1 M_⊙ and compared them to the observed distribution obtained by Ricci et al. (2010). They use a fixed disk model with R_c = 60 AU and γ = 1. The growth and fragmentation model predicts an increase of β(r) at radii ranging from 40 to 100 AU (their Fig. 3), which broadly agrees with our finding. However, in their analysis the distribution of average α_{1 − 3   mm} vs. observed 1.3 mm flux density only roughly matches the strongest sources. This effect is related to the dependence of the “fragmentation barrier” grain size on the surface density: less massive disks are expected to have smaller grains, and thus larger α_{1 − 3   mm}. As mentioned in Sect. 5.1, a simple way to provide a better agreement with the data is to assume that the faint disks are actually not less dense than the strong ones, but truncated to smaller sizes. In this case, α_{1 − 3   mm} is expected to decrease with source flux because of the radial dependence of β(r).

Finally, the discussion of the viscous properties of the disk presented in Sect. 5.3.1 used γ and R_c derived from the dust content, i.e. they implicitly assume a constant dust-to-gas ratio. In reality, the dust-to-gas ratio is expected to change with radius because of the concurrent effects of accretion, viscous spreading, grain drift, growth and fragmentation. In general, it is expected to decrease with radius, because the coupling between dust and gas increases with density, see for example Fig. 6 of Birnstiel et al. (2010a). The average dust-to-gas ratio is also expected to decrease with time, as the largest particles are drifting inward, being eventually advected onto larger bodies, either embryos or the central star. Accordingly, the discussion on possible changes of viscosity (α parameter and/or γ) presented in Sect. 5.3.1 should be taken with some additional care.

6. Conclusions

We report here the results of the first dual-frequency and high-resolution study of dust disks in the mm domain where the dust is mostly optically thin.

• Independent data sets allowed us to verify the robustness of thederived parameters and of their error bars. The geometricparameters (inclination and orientation) agree well withdeterminations from other constraints, such as scattered lightimages, optical jets, and the Keplerian rotation of the disks.

• We derived proper motions for 10 sources in our sample.

• Tidal truncation is found to affect the disk sizes in binary systems.

• Despite the combination of high angular resolution and sensitivity, we found that the viscous disk model does not generally provide a significantly better fit of the continuum data only than the simple truncated power law description. Baselines well above 300 kλ are required to distinguish between these two descriptions. In very young sources, the simple power law model appears to work somewhat better, while the exponential edge is marginally better for evolved objects.

• Inner holes also appear to provide a better explanation than negative values of γ for sources showing a deficit of emission at the center like GM Aur and LkCa 15.

• We have strong evidence for radial dependence of the dust emissivity exponent β with radius. In all cases, β is found to increase with radius, i.e. we find grain size which decreases with distance from the star. High β values (1.7−2, typical for ISM grains, or even possibly higher) are found beyond 100 AU, while the inner regions may display values down to nearly 0. This result is obtained whatever disk model has been adopted (surface density shape and temperature profile).

• We have possible evidence for optically thick cores in a few sources, which provide a direct estimate of the temperature of large grains. However, in some cases, inner regions with β = 0 may be misinterpreted as thick cores at low temperatures.

• Despite the ambiguities introduced by the variable dust properties, the characteristic size of the disk appears to increase with stellar ages, which broadly agrees with the viscous evolution. A more detailed comparison with the models suggests a decrease of the α viscosity parameter with time, as well as changes in the exponent of the viscosity.

These observations provide the first evidence for the expected effect of the dust grain evolution in circumstellar disks resulting from grain growth, fragmentation and, viscous transport. The comparison with model predictions is limited by the angular resolution obtained at the longest wavelengths, about 100 AU, which requires some parametric approach to constrain the radial dependence of β(r) (and by inference, κ(r) using some specific dust model). With the advent of ALMA and e-VLA, a direct inversion of the β(r) profile at linear resolutions of order 10−20 AU will become possible, enabling us to derive much more accurate constraints on the dust properties as a function of radius. This will be possible not only using two wavelengths as here, but over more than a decade in frequency.

Acknowledgments

We acknowledge the IRAM staff at Plateau de Bure and Grenoble for carrying out the observations. We also thank Rowan Smith who made a part of the numerical tests during a training period at LAB in July 2005. We thank the national programmes PCMI “Physico-Chimie du Milieu Interstellaire” and PNPS “Physique Stellaire” from INSU/CNRS for providing funding. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France, and of the NASA ADS Abstract Services.

References

Online material

Appendix A: Optical depth vs. variable β

Because a direct inversion of the brightness temperature profile is impossible, the determination of the parameters is fully implicit. Figure A.1 illustrates two possible types of brightness temperature profiles that can occur in our analysis. The continuous and dashed lines represent brightness at 1.3 mm and 2.7 mm for a typical power law distribution, with constant dust properties β(r) = β_m. The outer region is optically thin, and the slope constrains p + q. β_m is derived from the brightness ratio O. The inner region is optically thick, and constrains the exponent q as well as the temperature at 20 AU, T₂₀. The small brightness difference between the two frequencies is caused by the Rayleigh-Jeans correction. The dotted and dot-dashed lines represent an optically thinner disk at 1.3 and 2.7 mm respectively, with a viscous type profile with R_c = 150 AU. In addition, β(r) is assumed to vary with radius following Eq. (17) with R_b = 60 AU and R_w = 20 AU, κ(1.3 mm) being constant. Here, the inner region is optically thin, and its slope is q + γ. Note that if the temperature of that disk would be 4 times higher, it would mimic reasonably well the previous power law, optically thick case, provided γ is not too large. Accordingly, sources displaying a wavelength-independent flattened (apparent exponent  ≃ 0.4 − 0.7) inner brightness distribution can be interpreted either as optically thick sources, or as variable β(r) with β(r) ≃ 0 in the inner region. Steeper apparent exponents are not realistic for the temperature dependence. Note that the typical noise level is around 0.05 − 0.1 K in our observations at both wavelengths.

Appendix B: Sampling effects and best model

Because of the fully implicit derivation of the model parameters, an objective determination of the “best” model is difficult. The same source may be (nearly) equally well represented by either Model 1 or Model 2. We use a χ² criterium to determine the best matching model. However, it is important to realize that our data consists in a large (several 10⁴) number of statistically independent visibilities, each with very little (essentially zero) signal-to-noise. The χ² is given by (B.1)where O_i are the (complex) observed visibilities (O² actually being used to note O × O ^∗ , O ^∗  being the complex conjugate of O), M_i the modeled visibilities. The weights are derived from the theoretical noise using the system temperature, antenna gain, observing bandwidth and integration time. In general, , where M^b is the best-fit model, so even the null model M_i = 0 yields a χ² on the order of N, the number of visibilities, as W_i is the inverse of the variance of . Thus, the reduced χ², is a poor evaluation of the fit quality, which is close to 1 even for a very poor (null) model. Only the relative differences Δχ² between models of equivalent number of parameters can reveal whether one is better than the other.

Another subtle effect in comparing absolute values of χ² is the impact of discretization. A numerical model M is an approximation of the theoretical model T, M = T + E, where E is a numerical error term. So Because the model fit the observations and the numerical errors are not correlated with the observations, the last term is negligible, consequently the final χ² is a sum of the true (no numerical errors) term plus an offset cause by numerical effects. To make numerical errors negligible requires to be much less than 1. This is especially important when comparing different theoretical models. However, within a given model, the best-fit parameters may be determined with sufficient precision even if the numerical error term is not small.

Appendix C: Impact of the assumed temperature law

In this appendix, we investigate the impact of the dust temperature profile on the derived disk parameters. We consider two different profiles. Profile (i) is a power law T(r) = T₁₀₀(r/100AU) ^− q. Profile (ii) is a broken power law: it has a constant temperature between R_i = 40 AU and R_f, R_f being a variable parameter, while for r < R_i or r > R_f, the temperature is a power law with exponent q = 0.5, with to T(r) = T₁ at r = 1 AU. The temperature law is continuous as a function of r, and we used T₁ = 200 K by default. We analyzed the observations of a few sources (DL Tau, DM Tau and MWC 480) to explore the dependency of the derived surface density parameters on T₀,q and R_f. Figure C.1 illustrates the main impact of the temperature law on the surface density parameters, which is applicable to all optically thin sources. Figure C.1 is for Model 2 (so p is to be interpreted as γ), but similar results are obtained for Model 1.

For Profile (i):

• Σ₀ is nearly proportional to 1/T, with small corrections at low T owing to deviations from the Rayleigh-Jeans behavior.

• Similarly, p + q is nearly constant. This is equally valid for Model 1 (power law) and Model 2 (tapered edge).

• In Model 1, R_out is only weakly affected by the changes in p

• In Model 2, R_c increases by 20 to 30% when q increases from 0 to 0.5.

For Profile (ii):

• In Model 1, R_out slightly decreases with R_f (by about 10%), and p changes by about 0.1. Variations are not fully monotonic, however.

• In Model 2, R_c decreases by about 20 to 30%, when R_f goes from 50 to 200 AU. This is similar to the effect of q in Profile (i), as increasing R_t flattens the temperature distribution.

For more optically thick sources, like MWC 480, the effect on p is larger, because of the opacity corrections. However, in this case, q can be determined from the observations, because the χ² significantly depends on its value. Restricting the range of q to within its typical uncertainty limits the impact on p to about 0.2.

Except for the absolute scaling of the density as 1/T₁₀₀ (or 1/T₁ in Profile (ii)), the derived density distribution are thus not significantly affected by the assumed temperature law.

More importantly, R_c and p are affected in the same proportions at both wavelengths. Thus, the uncertainties on the temperature law have no significant effect on the derivation of the radial dependence of β(r) (see Fig. C.2). Incidentally, we note that in DL Tau, a better fit to the observations is obtained using Profile (ii) with R_f = 100 AU.

Appendix D: Disk parameters from observable quantities in the viscous model

The shape of the surface density profile used in Model 2 corresponds to the self-similar solution of the viscous evolution of a disk under the assumption that the viscosity is constant in time and a power law of radius (see Lynden-Bell & Pringle 1974; Pringle 1981). Under these assumptions, the surface density as a function of time and radius is given by (Eq. (17) of Hartmann et al. 1998) (D.1)where r = R/R₁, T = (1 + t_∗/t_s) is a dimensionless time, t_∗ the disk/star age and t_s is the viscous timescale at R₁, defined by (D.2)Our observations (at unknown time T) are characterized by the surface density law described by our Eq. (5) (D.3)So by identification, we obtain (D.4)and (D.5)which, eliminating R₁ usinq Eq. (D.4) (D.6)A time derivative of Eq. (D.1) (taken for r = 0) further indicates that the mass accretion rate is (D.7)We have in principle five unknowns (C,R₁,T,ν₁,γ), and five measurements: three from our study (M_d   or   Σ₀,R_c,γ), the stellar age t_∗ from evolutionary tracks and the mass accretion Ṁ, usually derived from the accretion luminosity (see Gullbring et al. 1998). Although this formally yields a solution, it is nearly degenerate when one considers the uncertainties on the measured quantities. This can be realized by noting that the mass accretion rate can be rewritten as (Eq. (14) from Isella et al. 2009) (D.8)while from Eqs. (6) and (D.5), the time dependency of the disk mass is simply (Eq. (A7) from Andrews et al. 2009) (D.9)so, by simple elimination (D.10)which simply gives (D.11)This is the only equation involving t_s. R_c, and thus R₁ does not appear in this expression because R₁ only reflects the initial condition of disk size, not its future evolution.

The (time independent) viscosity at any arbitrary radius is given by

(D.12)which, using the expression of R₁ in Eq. (D.4), can be expressed in terms of the observable quantities as

(D.13)It is customary to express it in terms of the α parameter, ν(r) = α(r)c_s(r)H(r), where c_s is the sound speed, and H(r) the scale height (D.14)In hydrostatic equilibrium, (D.15)T_g being the gas temperature in the disk mid-plane. Approximating T_g(r) by a power law of exponent  − q (q ≃ 0 − 0.6), we derive R_c and γ are directly constrained by our observations, while t_∗ is derived from evolutionary tracks. The last term c_s(R₀)H(R₀) depends on the stellar properties. We note from Eq. (D.15) that , and thus αC_s scales to first order as .

Appendix E: Alternate disk models

With the alpha prescription of the viscosity (radially uniform and constant in time α) and a (time independent) power law temperature T_k = T₀(r/r₀) ^− q, , so γ = 3/2 − q, Eq. (D.1) can also be written as (E.1)At long times, T ≫ 1, the density profile evolves as p = 3/2 − q, or p + q = 3/2.

A similar formula can also be recovered for the β prescription of the viscosity, (Richard & Zahn 1999). It is equivalent to setting q = 1 in Eq. (E.1), and thus corresponds to γ = 0.5.

The self-similar solutions of the evolution equation for the disk surface density were obtained under several simplifying assumptions. Desch (2007) pointed out that accounting for the early planet migration as predicted by the Nice model (Tsiganis et al. 2005; Gomes et al. 2005), the initial exponent of the surface density for the Solar Nebula would be very close to p = 2.2. To explain this slope, Desch (2007) recovered a different shape for the surface density in steady state configuration. The general form of the surface density in the Desch (2007) solution is (E.2)where r_u is the radius at which the disk has an apparent slope p and x_u = (2 − p − q)/(p + q − 3/2). For p + q > 2, x_u < 0 and the surface density vanishes at radius . Note that the classical steady-state result Σ(r) ∝ r ^{− (3/2 − q)} corresponds to the asymptotic limit x_u → ∞, and is obtained by imposing different boundary conditions on the evolution equation of angular momentum.

Appendix F: Unresolved, possibly thick, sources

For unresolved sources, the outer radius can only indirectly be constrained from the observed flux. Assuming uniform opacity τ, and a standard power law for the temperature T(r) = T₀(r/R₀) ^− q, the outer radius is given by (i being the inclination) (F.1)

A lower limit is recovered or i = 0 and τ → ∞ (F.2)The disk mass is given by (F.3)With q ~ 0 − 0.5, a lower limit on M_d is obtained for τ ≃ 0.5. Solutions with density/opacity decreasing with radius will lead to higher masses.

Appendix G: Figures for individual sources

We display here the figures for individual sources.

All Tables

All Figures