C. Pinte^1,2 - L. Fouchet^3,4 - F. Ménard¹ - J.-F. Gonzalez³ - G. Duchêne¹

1 - Laboratoire d'Astrophysique de Grenoble, CNRS/UJF UMR 5571, 414 rue de la Piscine, BP 53, 38041 Grenoble Cedex 9, France
2 - School of Physics, University of Exeter, Stocker Road, Exeter EX4 4QL, UK
3 - Université de Lyon, 69003 Lyon, Université Lyon 1, 69622 Villeurbanne, CNRS, UMR 5574, Centre de Recherche Astrophysique de Lyon; École Normale Supérieure de Lyon, 46 allée d'Italie, 69364 Lyon Cedex 07, France
4 - Institute of Astronomy, ETH Zurich, Schafmattstrasse 16, HPT D19, 8093 Zurich, Switzerland

Abstract
Aims. Our objective is to study the vertical dust distribution in the circumbinary ring of the binary system GG Tau and to search for evidence of stratification, one of the first steps expected to occur during planet formation.
Methods. We present a simultaneous analysis of four scattered light images spanning a range of wavelength from 800 nm to 3800 nm and compare them with (i) a parametric prescription for the vertical dust stratification, and (ii) with the results of SPH bi-fluid hydrodynamic calculations.
Results. The parametric prescription and hydrodynamical calculations of stratification both reproduce the observed brightness profiles well. These models also provide a correct match for the observed star/ring integrated flux ratio. Another solution with a well-mixed, but "exotic'', dust size distribution also matches the brightness profile ratios but fails to match the star/ring flux ratio. These results give support to the presence of vertical stratification of the dust in the ring of GG Tau and also allow us to predict the presence of a radial stratification.

Key words: radiative transfer - stars: circumstellar matter - methods: numerical - polarisation - scattering

1 Introduction

The first step of planet formation is the coagulation of dust grains, during low speed collisions (Dominik et al. 2006, and references therein). Micrometer size grains grow to form aggregates that will finally give birth to kilometer size planetesimals and eventually planets. In parallel to grain growth, dust grains are expected to settle towards the disk midplane due to the conjugate actions of the stellar gravity and gas drag.

Dust settling arises because gas is pressure-supported (and thus has a subkeplerian velocity) whereas dust is not. Therefore, gas slows down the dust resulting in the dust settling towards the midplane and spiraling inwards to the central star because of angular momentum conservation (see Stepinski & Valageas 1996; Weidenschilling 1977). This settling is highly dependent on the grain size with the larger grains (>1 m) almost unaffected by gas drag and the smaller ones (<1 $\mu$ m) so strongly coupled to the gas that they follow its motion. For the intermediate sizes, the dust distribution decouples strongly from that of the gas (Barrière-Fouchet et al. 2005).

The result is the formation of a dust sub-disk close to the equatorial plane (e.g., Safronov & Zvjagina 1969; Dubrulle et al. 1995), that may fragment and produce planetesimals when the density of the dust layer exceeds the gas one (Goldreich & Ward 1973). Despite remaining questions and uncertainties - how dust grains overcome the "meter size barrier'' without being accreted onto the star (e.g. Weidenschilling 1977) or destroyed by high speed collisions (Blum et al. 2000; Jones et al. 1996), how Kelvin-Helmholtz instabilities prevent the dust sublayer to fragment - models predict that dust grains can grow up to very large sizes and settle towards the midplane, leaving a population of small grains close to the surface and larger grains deeper in the disk.

Such a stratified structure should have consequences on disk observables like their spectral energy distributions (SEDs) and scattered light images (Dullemond & Dominik 2004). So far however, most studies of dust settling have been based on SED analysis only (Furlan et al. 2005; Apai et al. 2004).

An alternative approach to probe the disk structure is through high resolution scattered light images. Such images of protoplanetary disks are becoming available over a wide wavelength range in a few cases, from the optical with HST, to mid-IR with ground-based adaptive optics. Analysis of these images with model fitting were usually restricted to a single wavelength but multi-wavelength studies are becoming more and more common (see Watson et al. 2006, and references therein). Because scattering properties are highly dependent on grain size and because different wavelengths probe different depths into the disk, scattered light images over a wide range of wavelengths are potentially a powerful tool to detect dust settling. Indeed, attempts to simultaneously fit all wavelengths with a single model under the assumption of well mixed gas and dust indicate that the grain size distribution is likely more complex than that of the ISM (Mathis et al. 1977). In particular, a dependence of the maximum grain size with vertical distance above the disk midplane (vertical settling) sometimes appears to be needed (Duchêne et al. 2004, hereafter D04).

GG Tau is of special interest to study dust settling. Its geometry, well known from millimeter emission maps (Guilloteau et al. 1999), is simple: the dust is distributed in a circumbinary ring between 180 and 260 AU and the inclination is known, around 40 $^\circ$ , allowing a wide sampling of the scattering phase function. Multiple wavelength studies of the disk (D04) have shown that I-band scattered light images probe sub-micron grains located about 50 AU above the disk midplane, whereas at longer wavelengths, in the L' band, observations probe deeper regions of the disk (25 AU above the disk midplane) and reveal the presence of grains larger than 1 $\mu$ m, strongly suggesting the presence of a stratified structure.

In this paper, we investigate (i) whether multiple wavelength scattered light images of the GG Tau ring can be interpreted in terms of a disk model with a stratified structure and (ii) whether such a structure is required to explain the observations. In Sect. 2, we present model fitting of synthetic scattered light images. Two families of models are explored: (i) models where dust and gas are perfectly mixed in the disk and where optical properties of the grains are varied and (ii) models that include a stratified structure described in a parametric way. In Sect. 3, we explore if the model with a stratified structure, in better agreement with observations, is compatible with a more physical (hydrodynamical) description of dust spatial distribution. The effect of realistic dust settling on scattered light images of the GG Tau ring is studied by coupling hydrodynamical simulations of the dust distribution with a radiative transfer code. We discuss our results in Sect. 4 and conclusions are presented in Sect. 5.

2 Parametric models of the ring

2.1 Model description

Model fitting is performed on the multiple wavelength scattered light images obtained in the Hubble Space Telescope filters F814W, F110W, F160W by McCabe et al. (2002) and Krist et al. (2002) and on the Keck-AO image in L' band by D04.

Scattered light images are calculated by means of the Monte Carlo radiative transfer code MCFOST (Pinte et al. 2006). The general concept is to follow individual photon packets that randomly scatter in the disk, with probabilities given by the radiative transfer equations and dust optical properties, track them, and produce images when they exit the model boundaries.

In this paper, we focus our analysis on the stratified structure of the ring and we adopt the same ring geometry and dust grain composition as D04. The disk is assumed to be in hydrostatic equilibrium with a Gaussian vertical profile. The model parameters are the inner and outer radii (190 and 260 AU respectively), the total dust mass (0.0013 $M_\odot$ , Guilloteau et al. 1999), the vertical scale height (21 AU at 180 AU), the flaring index $H(r) \propto r^{1.05}$ and the surface density exponent $\Sigma (r) \propto r^{-1.7}$ . The density distribution extends inside the inner radius with a radial Gaussian decrease with a 1/2e width of 2 AU:

We assume dust grains to be homogeneous spheres (Mie theory) and adopt the porous dust grain optical properties from Mathis & Whiffen (1989, model A), with an 0.5 g cm^-3 average grain density. The grain size distribution $n(a,\vec{r})$ is defined locally, in each cell of the computation grid, and may then depend on the position within the disk, allowing to model disks with a stratified structure. Integrated over the whole disk, the grain size distribution is described by a power-law, ${\rm d} n(a) \propto a^{p} {\rm d} a$ . The minimum grain size is taken to be $a_{\rm min}=0.03~\mu$ m and $a_{\rm max}$ and p are considered as free parameters. With $a_{\rm max}=0.9~\mu$ m and p=-3.7, this dust model reproduces the interstellar extinction law up to 8 $\mu$ m (Mathis & Whiffen 1989). D04 showed that neither a grain distribution typical of interstellar medium nor the distribution used by Wood et al. (2002) to reproduce the SED of HH 30 are able to account simultaneously for both the scattered light images and I-band polarisation rates. In this section, we try to construct a grain size distribution that would reproduce both observed quantities simultaneously

Models are monochromatic and are calculated at the central wavelength of the different filters (0.81, 1.0, 1.6 and 3.8 $\mu$ m for the F814W, F110W, F160W and L' filters). Each model is computed with a total number of 128 million photon packets. The pixel size is 0.04 $^{\prime \prime }$ and images are convolved with a Gaussian kernel with a resolution similar to those of observations, i.e., ${\it FWHM}= 0.065$ $^{\prime \prime }$ , 0.080 $^{\prime \prime }$ , 0.128 $^{\prime \prime }$ and 0.091 $^{\prime \prime }$ for I, J, Hand L' bands, respectively.

The quality of the models was determined by a $\chi ^2$ minimization method. The disk geometry was validated in previous studies (e.g., Guilloteau et al. 1999; Duchêne et al. 2004; Silber et al. 2000) and in this paper we focus instead on the scattering properties, in particular the scattering phase function. In a first step, we only consider the azimuthal brightness profiles in our fitting procedure. These profiles were constructed with the same method as described in D04 (see their Sect. 4.2) using $10^\circ$ sectors. For each wavelength, we define a reduced $\chi ^2$ , based on the uncertainties in the images. Model uncertainties are negligible with respect to the observations due to the large number of photon packets generated for each model. We finally define a total $\chi ^2$ as the sum of the $\chi ^2$ at each wavelength. From this total $\chi ^2$ , a Bayesian analysis is performed to estimate the probability of occurrence of each parameter value. It also gives an estimate of the validity range of each "best-fit'' value (see for instance Press et al. 1992 and also Lay et al. 1997 for an application to millimeter observations of HL Tau). Because our model is simplistic (the observations reveal systematic departures from the axisymmetry assumed in the model), this $\chi ^2$ should not be interpreted in statistical terms but as defining a metric to compare the general behaviour of different models.

2.2 Models without stratification

Modification of the grain size distribution or use of complex grain size distributions results in a dramatic increase of the parameter space. Here, our aim is rather to verify how one can, by increasing in a reasonable manner the parameter space, modify the grain distribution used by D04 so that it reproduces simultaneously all observations, without invoking any stratification.

Scattered light images at multiple wavelengths probe different grains inside a single distribution. Thus, for a grain size distribution typical of the ISM with $a_{\rm max}=0.9~\mu$ m, the median scatterer has a size of 0.5 $\mu$ m in the I band and 0.7 $\mu$ m in the L' band (D04).

A simple way to play on the characteristic size observed at each wavelength is to modify the relative fraction of each grain size via the slope of the size distribution. The goal is to introduce large grains, susceptible to reproduce the L' band image, and therefore larger than 1 $\mu$ m, but in a small enough number so that their contribution at small wavelengths is negligible.

To do so, we explore a 2 dimensional parameter space, by varying the maximum grain size between 1 $\mu$ m (size slightly smaller than the size required to reproduce the L' band image in the work of D04) and 100 $\mu$ m (size large enough that the contribution from these grains is negligible on scattered light, even at 3.8 $\mu$ m), and the index of the grain size distribution between -3.5 (as in the ISM) and -6.0, by steps of 0.5.

$\begin{figure} \par\psfrag{Angle}[c][c]{Position angle ($^\circ$ )} \psfrag{Int... ...Normalized intensity} \includegraphics[width=8.8cm]{7137fig1.eps} \end{figure}$	Figure 1: Azimuthal brightness profiles of the best model without stratification (solid line) superimposed to observed profiles (triangles). The model parameters are p = -5.5 and $a_\textrm{max} = 100~\mu$ m. The position angle $0^\circ$ corresponds to the disk semi-minor axis, on the front face of the ring.
Open with DEXTER

Figure 1 presents azimuthal brightness profiles for the best model. Results are in agreement with observations at all wavelengths with a global $\chi ^2$ of 57.28 (Table 1). The best model reproduces very well the images in the I, J and H bands but slightly overestimates the flux on the back side of the ring in the L' band.

Table 1: $\chi ^2$ of the best models without stratification. Each line gives the $\chi ^2$ in the I, J, H and L' bands, as well as the total $\chi ^2$ defined as the sum of the individual ones for the models minimizing one of these $\chi ^2$ , respectively.

The Bayesian analysis shows that a single index for the grain size distribution, equal to -5.5, can reproduce the observations (Fig. 2). While our crude sampling prevents us from accurately defining the allowed range for the slope index, we can readily exclude values smaller than -6.0 and larger than -5.0, leading to a tight constraint. Regarding the maximum grain size, we obtain that a minimum value of 3 $\mu$ m is required but all values beyond that are equiprobable.

2.3 Models with stratification

We model the vertical stratification of the disk with a simple parametric law where the vertical scale height of the grains depends on the grain size:

The particular geometry of the GG Tau ring in which the back side of the inner rim of the ring is seen directly suggests that such a stratified model is too simplistic to describe scattered light observations. In Fig. 3 a synthetic map obtained in the I band with the stratification described above is shown. The inner edge of the ring, observed on the back face, presents a vertical brightness gradient (the surface is about twice as bright as the midplane) which is not seen in observations (Krist et al. 2005,2002).

$\begin{figure} \par\includegraphics[angle=190,width=4.5cm,clip]{7137fig3.eps} %\end{figure}$	Figure 3: Synthetic image of the GG Tau ring in the I band with only vertical stratification. The color map is in logarithmic scale and is inverted with bright regions in dark. The image size is 3.7 arcsec.
Open with DEXTER

This gradient is caused by the fact that, as a function of altitude at which photons hit the ring, they do not scatter on grains of the same size. Photons that scatter in the upper layers of the ring interact with small grains, which are isotropic scatterers and send a significant fraction of energy back to the observer. On the contrary, photons that scatter close to the midplane interact with larger grains, which are more forward-scatterers and send only a smaller fraction of energy back to the observer.

Thus, the different scattered light images cannot be explained only with vertical stratification. A complementary stratification, in the radial direction appears needed to explain the globally uniform brightness observed on the back side of the ring. A configuration where small grains would form a layer enclosing larger grains, as presented in Fig. 4, might produce images in agreement with observations. At optical wavelengths, photons would scatter in the surface layer and meet small grains, whatever their height relative to the disk midplane. In the infrared, photons would penetrate deeper in the ring and scatter on larger grains.

$\begin{figure} \par\includegraphics[width=6.8cm,clip]{7137fig4.eps} \end{figure}$	Figure 4: Sketch of the vertical and radial stratification proposed to account for the scattered light images of GG Tau.
Open with DEXTER

To represent such a stratification, we introduce an analytical description, similar to the one used for the vertical stratification, so that the density remains a continuous function. The characteristic length of the radial decrease of the surface density, $\Delta r$ , becomes a function of the grain size, that we choose to describe with a power-law, whose index is the same as the vertical stratification one:

We then explore a 2D parameter space, playing with maximum grain size between 1 and 100 $\mu$ m as previously and the stratification index, $\xi$ , that we sample between 0 (no stratification) and 0.5. In this section, we set the index p of the integrated grain size distribution to a constant value of -3.5, following models from D04.

$\begin{figure} \par\psfrag{Angle}[c][c]{Position angle ($^\circ$ )} \psfrag{Int... ...Normalized intensity} \includegraphics[width=8.8cm]{7137fig5.eps} \end{figure}$	Figure 5: Azimuthal brightness profiles of the best model with stratification (solid line) superimposed to observed profiles (triangles). The model parameters are $\xi = 0.15$ and $a_\textrm{max} = 60~\mu$ m.
Open with DEXTER

Figure 5 presents azimuthal brightness profiles obtained for the best model, compared to observed profiles. Results are in good agreement with observations at each wavelength, with an azimuthal profile similar at all wavelengths. The intensity on the back side of the ring, in the L' band is slightly larger, by about 1 $\sigma$ , with respect to the observations. In order to quantify further the quality of the fit, we compare $\chi ^2$ values obtained for the best model to those that we get by fitting one single wavelength (Table 2). The global $\chi ^2$ , 57.70, is comparable to the one obtained from models without dust stratification (57.28, see Table 1).

Table 2: $\chi ^2$ of the best models with stratification. The different lines have the same meaning as in Table 1.

We note that, in the I, J and H bands, the best global model is only marginally poorer that the best models at each of these wavelengths but that the difference is more important for the L' band similarly to the models without stratification.

Results of the Bayesian analysis are presented in Fig. 6. The stratification index is well constrained, with a probability curve peaking between 0.125 and 0.15 and a dispersion of about 0.025. Regarding the maximum grain size, it appears that values smaller than 5 $\mu$ m are excluded but that beyond this limit, all values are roughly equiprobable, grains larger than 5 $\mu$ m contributing very little to the scattered flux at wavelengths smaller than 3.8 $\mu$ m.

2.4 Polarisation and flux ratio

Both models, with radial and vertical stratification or with a slope for the size distribution equal to -5.5, present equally good fits to the azimuthal brightness profiles. In order to discriminate between these two solutions, we compare their predictions of the polarisation and of the disk-to-star flux ratio with observations.

$\begin{figure} \psfrag{pola}[c][c]{Linear polarisation level} \psfrag{angle}[c... ...on angle ($^\circ$ )} \includegraphics[width=8.8cm]{7137fig7.eps} \end{figure}$

Figure 7: Observed (from Silber et al. 2000) and synthetic polarisation at 1 $\mu$ m as a function of the azimuthal angle. Polarisation maps have been binned by a factor 4, resulting in a pixel size of 0.172 $^{\prime \prime }$ and then divided in 20 $^\circ$ sectors where the polarisation is measured, for each sector, in the pixel of maximum brightness. Black circles represent the polarisation extracted from observations. The solid red line corresponds to the best model with stratification and a slope p=-3.5 and the dashed blue line to the model without stratification and a slope p=-5.5. The position angle $0^\circ$ is defined as in the previous figures and corresponds to the disk semi-minor axis, on the front side of the disk.

Open with DEXTER

Figure 7 presents the observed polarisation in the J band (Silber et al. 2000) as a function of the azimuthal angle compared with the polarisation predicted by the two best models and shows the differences between them. The model with stratification presents a polarisation rate around 15% on the front side and peaking around 65-70% on the back side of the disk whereas for the model without stratification but with a grain size distribution slope of -5.5, the polarisation rate is 10% on the front side and a maximum polarisation level of 80-90% on the back side. None of the two models provides a good fit to the observations of Silber et al. (2000). The model with stratification performs slightly better ( $\chi^2 = 16$ against a $\chi ^2$ of 33 for the model without stratification) and reproduces reasonably well the polarisation in the regions where it is the best defined, the front face of the disk. The model without stratification results in somewhat lower polarization rates (on the front side). However, both models overestimate the polarisation on the back side of the disk. This behaviour is due to the high polarisability of the porous grains used here. An exploration of the dust properties, beyond the scope of the paper, is needed to obtain a better agreement with observations and, at this point, polarisation cannot be used as a strong constraint to discriminate between the two models.

Table 3 gives all disk-to-star flux ratios for observations and models. These measurements depend on the ring geometry, in particular to its projected surface on the plane of the sky. Because we did not explore this parameter, we focus here on the variations of these ratios with wavelength. The model with dust stratification is in better agreement with observations: the flux ratio is roughly constant between the I and H bands and decreases by about 15% in the L' band. Even though it does not exactly reproduce the observations where the ratio slightly increases from the I to H bands, and then decreases in the L' band ( $\approx$ $35\%$ smaller than in the I band), the global behaviour is respected, which is not the case for the model without stratification where a steady flux ratio decrease with wavelength is witnessed, with a factor around 6 between the I band and the L' band.

Results presented in this section were based, for both cases, on the best model according to $\chi ^2$ minimizations. The study of brightness profiles only gives a lower limit for the maximum grain size. Polarisation and colours show the same insensitivity to the precise value of $a_{\rm max}$ as the phase function. Results presented are therefore valid for all values $a_{\rm max} \ga 5~\mu$ m and, in particular, cannot help to better constrain this parameter.

3 Hydrodynamical modeling

The parametric modeling of GG Tau's circumbinary ring presented in the previous section points towards a stratified disk where both vertical and radial segregation of grains are needed.

In this section, we aim at testing the validity of this hypothesis, using a more physical modeling of the processes of vertical settling and radial migration of dust grains. Results of hydrodynamical simulations, where dust settling and migration are computed more realistically, are used as an input to the radiative transfer code to produce synthetic scattered light images.

3.1 Hydrodynamical code

To compute the hydrodynamical evolution of gas and dust, we use Smoothed Particle Hydrodynamics (SPH), a Lagrangian technique described by Monaghan (1992). The equations and approximations are rigorously established by Bicknell (1991).

Our code was derived from Murray's code (Murray 1996) and is described in Barrière-Fouchet et al. (2005). It is fully 3D, bi-fluid, non self-gravitating and locally isothermal with the temperature law given by $T(r)\propto r^{-3/4}$ , where T is the temperature and r the distance to the star projected onto the midplane. The initial surface density is given by $\Sigma \propto r^{-3/2}$ , following estimations of the Minimum Mass Solar Nebula (MMSN) and of Guilloteau et al. (1999).

The dust phase, representing one single grain size at a time, is treated as a pressureless fluid and the coupling between gas and dust is done through aerodynamic drag (see Monaghan 1997; Monaghan & Kocharyan 1995). The drag exerted by dust on gas is neglected (see details in Barrière-Fouchet et al. 2005).

In the initial state, dust and gas are considered as well mixed, which means that the dust disk is flared. This is correct for the smaller grain sizes when dust is strongly coupled to the gas, but for larger grain sizes, we can expect grains to start close to the midplane because they are produced by the growth of smaller grains that have already settled. Still, up to 1 mm, grains settle quickly enough to the midplane that this initial state does notreally matter. The contribution of larger grains to the synthetic images is negligible (see Sect. 2.2).

3.2 Parameters

We run hydrodynamical simulations for a circumbinary disk around two stars of 0.65 and 0.5 $M_\odot$ respectively, separated by 35 AU and evolving on orbits with 0.25 eccentricity (Dutrey et al. 1994). We start from a disk with a 0.5 AU central hole and let the binary dig the hole during the 10 000 first time-steps ( $\simeq$ 8 orbits at 100 AU) while the gas relaxes to a quasi equilibrium state before dust is injected. We then let the whole evolve for 40 000 more time-steps ( $\simeq$ 32 orbits at 100 AU).

The disk mass is 0.13 $M_\odot$ with 1% dust in mass, its outer boundary is initially fixed to 350 AU and doesn't significantly spread beyond 400 AU. Each simulation includes about 25 000 particles (gas+dust). The grain sizes range from 1 to 10 $\mu$ m, we omit larger grains because their contribution to the scattered light images is negligible. Table 4 summarizes the simulation parameters.

The inner edge radius obtained in the simulations, about a hundred AU, is appreciably smaller than GG Tau's (180 AU). This is due to the uncertainty on the binary parameters for GG Tau and the depletion of the disk's central regions is therefore different. This simply introduces a scaling factor and does not alter in any way the results presented here. More recent binary parameters that could better account for the location of ring's inner edge have been obtained but still need to be confirmed by further astrometric monitoring, as noted by Beust & Dutrey (2005,2006).

3.3 Simulation results

$\begin{figure} \resizebox{8.8cm}{!}{ \includegraphics{7137fig8.eps} } \end{figure}$	Figure 8: Density contours for dust ( left) and gas ( right) at the end of the SPH simulations of the circumbinary disk for runs including 10, 5 and 1 $\mu$ m grains from top to bottom.
Open with DEXTER

Figure 8 displays the dust and gas densities in the (r, z) plane at the end of the SPH simulations. The settling of dust appears clearly compared to the gas disk. The outer dust density contours show it to be more efficient for larger grains (10 $\mu$ m) than for smaller ones (1 $\mu$ m) that remain strongly coupled to the gas. The radial distribution also differs from one grain size to the next. The 1 $\mu$ m grains move closer to the star than the 5 and 10 $\mu$ m grains (see also Fig. 10), which provides physical support to our parametric model of radial and vertical stratification (Sect. 2) in which the small grains form a layer surrounding the central parts of the ring where the large grains are present.

The reason for this radial segregation is that larger grains respond more strongly to the pressure gradient and are more efficiently concentrated in the gas pressure maxima. Therefore, the radial extent of the density distribution decreases as the grain size increases (this is more easily seen in Fig. 10). Any asymmetry between the inner and outer pressure gradients will produce a shift in the peak of the density distribution according to the grain size. In the case of GG Tau, the pressure gradient is stronger in the inner side than in the outer side, which results in the density peak being shifted outwards for larger grains.

3.4 Interfacing the two codes: analytical fits to the data

The SPH data are unfortunately too noisy to be used directly in MCFOST, even after smoothing on a regular grid. We therefore derive analytical fits that must consistently reproduce this differential settling. To do this, we first perform a vertical fit of the dust density at the end of the SPH simulations at a given radius, the parameters of this fit are then fitted in the radial direction.

3.4.1 Vertical fit

We first decide of a shape for the fit in the vertical direction. According to Garaud & Lin (2004), the time evolution of the dust density with height is self-similar, with initial conditions satisfying the hydrostatic equilibrium. This leads to the classical expression of density for a geometrically thin disk (in the vertically isothermal case):

We produce a vertical fit at a given radius r₀ with two parameters $\rho_{\rm max}$ and $\sigma_0$ so that

3.4.2 Radial fit

$\begin{figure} \par\includegraphics[width=8.8cm,clip]{7137fig9.eps} \end{figure}$	Figure 9: Fits for the 5 $\mu$ m dust particles in the circumbinary disk. a) Scaleheight. b) Surface density. In both cases, the dotted line corresponds to the vertical fit only and the solid line to the subsequent radial fit. c) SPH density after smoothing on a regular grid. d) Vertical fit of the density. e) Complete fit of the density in both vertical and radial directions.
Open with DEXTER

In the case of a circumbinary disk, the surface density, plotted against the logarithm of the radius ( $\ln\Sigma = f(r)$ ), looks like a Gaussian because the disk is actually a wide ring and there is a density drop towards the center as well as towards the outside of the disk. We choose a fitting formula where we look for the value of the three parameters $\Sigma_{\rm max}$ , $r_{\rm0max}$ , $\sigma'_0$ , verifying

Finally, with the fits of H, of $\Sigma$ and the theoretical expression of the density (Eq. (4)), we get an analytical expression for it:

The fitted scaleheights and surface densities for different grain sizes are plotted in Fig. 10, which clearly shows the differential settling, on the left panel, with the larger grains having the more important settling (because we are in the small particle regime, see Barrière-Fouchet et al. 2005). The right panel shows that the dust density maximum shifts towards larger radii as the grain size increases, and therefore as grains settle more efficiently in this size range (see Sect. 3.3).

3.4.3 Interpolation between grain sizes

Each hydrodynamical simulation only includes one grain size at a time. To synthesize scattered light images, we need to reconstruct a grain size distribution at each point of the disk. Because we do not consider interactions between dust grains and because all simulations, whatever the grain size used, predict consistent density for the gas phase, we can do this by superimposing the calculated densities for each grain size. The relative fraction of grains is fixed by assuming a grain size distribution, integrated over the disk, described as a power law with an index of -3.5. Finally, an interpolation, logarithmic with respect to grain sizes, is performed between the spatial distributions of grains resulting from the SPH calculations, to get a well sampled grain size distribution. The smallest grains ( $a_{\rm min}=0.03~\mu$ m) are assumed to follow perfectly the gas distribution.

3.5 Synthetic images

Figure 11 compares the azimuthal brightness profiles in the I, J, H and L' bands resulting from our simulations of the circumbinary disk, with and without grain stratification. In the latter case, the aerodynamic drag of gas on dust is switched off and the dust distribution follows that of the gas for all grain sizes. The differences between models with and without stratification are especially apparent on the back side of the disk. As for the analytic stratification we used in Sect. 2.3, we observe the smaller grains at the inner rim of the disk. They are isotropic scatterers and send back a larger fraction of light than the larger grains, of which we detect the scattered light when all grains are well mixed. We want to stress here that the models are not fits to the data, which are plotted as a reference in order to study the qualitative behaviour of our models.

$\begin{figure} \par\psfrag{intensite normalisee}[c][c]{Normalized intensity} \p... ...e L'}[c][c]{L' band} \includegraphics[width=8.8cm]{7137fig12.eps} \end{figure}$	Figure 11: Azimuthal brightness profiles of the circumbinary disk with (solid line) and without (dashed line) stratification superimposed to GG Tau's observed profiles (triangles).
Open with DEXTER

The inclusion of dust segregation is a step in the right direction to explain GG Tau's observations. In particular, we note that the back side of the disk is particularly sensitive to the radial stratification. The light coming from the front surface of the ring, which corresponds to angles between 0 and 90 $^\circ$ on one side and 270 and 360 $^\circ$ on the other side, is in fact little affected by the segregation and both brightness profiles (with and without stratification) are almost identical. However, the back side is noticeably brighter in the stratified case, for the reasons mentioned above.

It appears that the stratification resulting from the present hydrodynamical model is not sufficient to explain the observed brightness profiles. This may be due to the relatively short time on which the dynamical evolution is followed, or to physical processes which also lead to a spatial segregation of grains in the disk but are not considered here, like a preferred grain growth in the denser parts of the disk. Nonetheless, our hydrodynamical simulations demonstrate that it is physically plausible that the ring possesses an "inner wall'' of small grains, as suggested by our parametric approach.

4 Discussion

The inclusion of stratification in a parametric model of GG Tau's disk (Sect. 2) allows to reproduce the scattered light images from the I to the L' band, in particular the azimuthal brightness distribution which is remarkably similar at all wavelengths. This model is also in good agreement with the 1 $\mu$ m polarisation data and globally reproduces the ring colours with respect to the central star. However, a simple parametric vertical stratification alone cannot explain all observations, especially the brightness of the inner rim. The presence of a radial segregation at the inner edge, with smaller grains closer to the stars than larger grains, also appears required to obtain synthetic images which reproduce the observations. This kind of distribution can result from a preferred grain growth in the central, denser zones of disks where the collision probability is higher. It is also predicted by simulations of the dust dynamics around a binary system (see Sect. 3.3), which presents synthetic scattered light images showing a good qualitative match with the observations, as shown in Sect. 3.5.

Adding dust stratification and settling into models that already contain several free parameters renders the problem of finding unique solutions even more acute. We explored a large parameter space and tried to quantify, through Bayesian analyses, the likelihood of occurrence of the models we presented. By doing so, we also found a solution with well-mixed dust that also matched the observed azimuthal profiles, at all wavelengths. However, this solution appears physically less realistic than the expected one with settling because it requires a very steep slope for the dust size distribution. Yet, the solution exists and we added observational constraints to confirm or reject its validity. When compared to polarisation data of the ring, the "well-mixed'' solution proves only marginally compatible with the data, and certainly worse at matching observations than the model with stratification. More interesting, the well-mixed models we explored are not able to reproduce the correct dependency of disk colours with respect to the star as a function of wavelength.

The best "well-mixed'' solution we identified is made of a distribution of spherical and homogeneous particles, with a slope of -5.5, noticeably steeper than that of interstellar medium grains. Models of dust coagulation (Ormel et al. 2007; Dullemond & Dominik 2005) show that, at least during the first steps of the process, the slope of the grain size distribution does not vary so much, in the size regime we are interested here. We succeeded until now in fitting a growing number of observable quantities by modifying the grain properties in a well-mixed model but this fit is still not perfect, and noticeably poorer than in the stratified case. These modifications are increasingly arbitrary and it becomes difficult to find physical causes to explain them. Therefore, we consider this solution with a steep dust size distribution more "exotic'' and less likely to occur. The model with stratification provides a better match to more data sets. Furthermore, stratification appears a more natural outcome of the disk evolution processes and for now it is our preferred model for the ring of GG Tau.

Krist et al. (2005) show evidence that the azimuthal brightness profiles change with time with front/back amplitude variations around 20%. These fluctuations remain moderate and taking them into account only slightly affects the major observational property of the ring in scattered light: the similar morphology of the disk over a large range of wavelengths. Therefore, our general conclusion that dust stratification is required to interpret scattered light images still holds but a more precise determination of the degree of settling will probably require comtemporaneous multiple wavelengths data sets.

In all cases, the presence of grains larger than a few microns is required. Both models require a minimum value of $a_{\rm max}$ of the order of 3 $\mu$ m, indicating that the grains have evolved compared to those of the interstellar medium. The value inferred from our simulations assuming spherical grains is indeed a lower limit, a particle of complex shape being able to mimic the scattering properties of a smaller spherical and homogeneous particle (Rouleau 1996). This implies that grain growth towards sizes larger than one micron is happening in this circumbinary ring and that the first steps of planet formation are very probably being initiated.

5 Conclusions

The definitive proof of a stratified structure in GG Tau's disk will probably rely on additional observations. However, stratified dust models are able to match several observable quantities of the GG Tau ring: in particular the azimuthal brightness profiles from 800 nm to 3800 nm, but also polarisation measurements at 1 micron and the variation of star/disk flux ratios with wavelength.

We have shown the attractiveness of polarimetric observations to better constrain the grain properties. The broadening of the wavelength coverage also seems very promising, both in intensity and polarimetry. The detection of the disk at 10 $\mu$ m, where our models show that scattered light still dominates over the disk thermal emission, would enlarge the wavelength leverage at our disposal and probe layers even closer to the disk midplane than the current 25 AU. The inclusion of the thermal emission allows additionally to derive an estimate of the ratio between the absorption and scattering cross-sections and to further restrict the families of solutions compatible with the observations of a given disk.

The modeling of the circumbinary disk in Sect. 3 allowed us to give a physical basis to the parametric modeling of GG Tau's ring in Sect. 2. Hydrodynamic models validate our hypothesis of a radial stratification of grains and are in qualitative agreement with the observations. It is indeed to the grains radial segregation, rather than to their vertical stratification, that we are most sensitive in GG Tau's disk. The configuration we used does not exactly correspond to GG Tau's, due to uncertainties in the binary's orbital parameters. The use of better parameters should allow us to make more quantitative comparisons with the data and to start to better constrain the processes causing the stratification in GG Tau's disk.

We would like to stress here that validating the parametric approach by hydrodynamic computations is important because it allows the use of a faster tool to quickly estimate whether settling has occurred in a disk or not, before embarking on dedicated hydrodynamic calculations, which are more realistic but necessarily longer.

On the long term, it will be necessary to include new physical processes acting in addition to dust dynamics. Turbulence, in the calculations presented here, was only taken into account in a simplified manner. It would be interesting to see how a better treatment of turbulence limits the settling and affects the observable quantities.

Grain growth is intimately linked to their settling. Because it enlarges the relative velocities, settling increases the number of collisions, which favours grain growth, which in turns increases the settling. The effects on the different observable quantities are potentially enhanced. The simultaneous modeling of both processes appears to be the next step towards a better understanding of the first stages of planet formation. It is currently being implemented into the SPH code (Laibe et al. 2007) and the resulting simulations will be coupled to the radiative transfer code to continue our study.

On the stratified dust distribution of the GG Tauri circumbinary ring