Gap opening beyond dead zones by photoevaporation

R. Morishima
University of California, Los Angels, California, USA
Jet Propulsion Laboratory, Pasadena, California, USA
E-mail: ryuji.morishima@jpl.nasa.gov
Accepted 2011 September 30. Received 2011 June 28.
Abstract

We propose a new hypothesis for the origin of protoplanetary discs with large inner holes (or gaps), so-called transition discs. Our gas disc model takes into account layered accretion, in which poorly-ionized low-viscosity dead zones are sandwiched by high-viscosity surface layers, and photoevaporative winds induced by X-rays from the central stars. We find that a gap opens at a radius outside a dead zone, if the mass loss rate due to photoevaporative winds exceeds the mass accretion rate in the dead zone region. Since the dead zone survives even after the gap opens, mass accretion onto the central star continues for a long time. This model can reproduce large gap sizes and high mass accretion rates seen in observed transition discs.

keywords:
methods: numerical – protoplanetary discs.
pagerange: Gap opening beyond dead zones by photoevaporationGap opening beyond dead zones by photoevaporationpubyear: 2011

1 Introduction

Spatial distribution of dust in protoplanetary discs is measured from disc continuum spectra in near-infrared to millimeter wavelengths whereas the gas density is measured directly by line emission from such as CO and indirectly from H emission which is used for estimations of accretion rates of gas onto central stars. Infrared excess due to dust emission disappears in 3 Myr for roughly a half of discs and for almost all discs in 10 Myr (Haisch, Lada & Lada, 2001). Gas in inner discs disappears in similar time-scales (Fedele et al., 2010), but the life time of gas in outer discs is not clearly known.

Recent observations, such by the Spitzer infrared space-telescope, the Submillimeter Array, and the Subaru telescope, revealed discs in the dispersal phase, so-called transition discs (Merín et al., 2010; Thalmann et al., 2010; Andrews et al., 2011). They have optically thick outer discs whereas inner discs with sizes up to 70 AU are optically thin. This indicates depletion of dust in the inner discs. Some of transition discs show near-infrared excess indicating formation of gaps instead of holes and some others have weaker mid to far infrared emission indicating optically thinner outer discs (Muzerolle et al., 2010). The fraction of transition discs in protoplanetary discs increases with time (Currie & Sicilia-Aguilar, 2011): 15-20 per cent at 1-2 Myr and more than 50 per cent at 5-8 Myr. A puzzling thing is that a large fraction (75 per cent; Merín et al. (2010)) of transition discs exhibit gas accretion onto their central stars and their accretion rates are close to those for classical T Tauri stars, yr (Hartmann et al., 1998).

There are three proposed explanations for the origin of transition discs but all of them do not seem to be fully satisfactory (Andrews et al., 2011; Williams & Cieza, 2011). (1) Gap opening by photoevaporative winds (Gorti, Dullemond & Hollenbach 2009; Owen, Ercolano & Clarke 2011). If the mass loss rate of photoevaporative winds exceeds the mass accretion rate of viscous evolution, a gap can open. However, before a gap opens, the mass accretion rate largely decreases and becomes much less than yr. This is the case even with X-ray photoevaporation with the mass loss rate as high as yr (Owen et al., 2011). The sizes of the inner holes due to photoevaporation cannot be larger than 20 AU when mass accretion onto the central star still remains. (2) Dust growth (Dullemond & Dominik 2005; Tanaka, Himeno & Ida 2005). This is expected particularly in inner discs as the models show that the growth time of dust has strong dependence on distance from the central star. However, dust growth alone cannot be responsible for transition discs, because the models for dust growth predict smooth opacity changes with radius whereas observed transition discs show clear gaps of dust. (3) Gap opening by giant planets. Zhu et al. (2011) show that multiple giant planets can open a gap as wide as the largest holes observed. However, depletion of the surface density in the gap results in a large decrease of the mass accretion rate.

The most probable source of the disc viscosity is the turbulent viscosity caused by the magneto-rotational instability (MRI; Balbus & Howley (1991)). Since the ionization rate in the thick inner disc is not high enough, a poorly ionized layer with a low viscosity called a dead zone forms near the disc midplane sandwiched by MRI active layers (Gammie 1996; Zhu, Hartmann & Gammie 2010). This type of accretion is named layered accretion. In this paper, we develop a gas disc model which takes into account layered accretion and photoevaporative winds. Layered accretion has not been taken into account in previous works of photoevaporation. As we will show, our disc model can reproduce some of important properties of observed transition discs, such as large gap sizes and high mass accretion rates.

In Sec. 2, we introduce our gas disc model. In Sec. 3, results are shown for both discs with and without layered accretion over wide ranges of parameters. In Sec. 4, implications for transitions discs are discussed. In Sec. 5, the summary of this study is given.

2 Methods

2.1 Basic equations

We assume that mass accretion from a molecular cloud core has already completed and that the mass of the central star is fixed during disc evolution. Time evolution of the surface density is derived from the mass and angular momentum conservation equations:

 ∂Σ∂t=12πr∂˙M∂r−˙Σw, (1)
 ∂(Σj)∂t=12πr∂∂r(˙Mj+˙Jvis)−˙Σwj, (2)

where is the mass accretion rate (assumed to be positive for inward accretion), is the mass loss flux of photoevaporative winds, is the specific angular momentum with being the Keplerian frequency, and is the angular momentum flux transferred by the disc viscosity, . From Eqs. (1) and (2), the mass accretion rate is obtained as

 ˙M=−2πrΣvr=6πr1/2∂∂r(Σνr1/2), (3)

where is the radial velocity. Substituting this form to Eq. (1), we obtain the well known diffusion equation (Lynden-Bell & Pringle, 1974):

 ∂Σ∂t=3r∂∂r[r1/2∂∂r(Σνr1/2)]−˙Σw. (4)

2.2 Disc viscosity

Using the viscosity parameter (Shakura & Sunyaev, 1973), the viscosity is given as

 ν=αch=αRTμΩ, (5)

where is the isothermal sound velocity at the disc midplane and is the vertical scale height of the disc. The isothermal sound velocity is given as , where is the gas constant, is the disc temperature at the midplane, and we take the mean molecular weight as = 2.33 g mol. For models without dead zones, is radially constant. For models with dead zones, the parameters for active surface layers and dead zones are defined to be and . Then, the effective disc viscosity parameter is given as

 α=Σaαa+(Σ−Σa)αdΣ, (6)

where is the sum of the surface densities of both sides of active layers. In simulations, we adopt . The surface density is modeled after Suzuki, Muto & Inutsuka (2010):

 Σa=min[ΣCR+ΣX(r1AU)−2,Σ], (7)

where is the surface density of layers ionized by cosmic rays and is the surface density of layers ionized by X-rays from the central star at 1 AU. Equations (6) and (7) indicate that if . On the other hand, if , the region is entirely MRI active () so that .

The surface densities and particularly depend on dust amount in the disc, as recombination of ions and electrons occurs on the dust surface (Sano et al., 2000; Bai, 2011). For example, Suzuki et al. (2010) uses = 12 g cm and = 25 g cm, whereas Gammie (1996) adopts nominal values as = 200 g cm and = 0. The value of adopted in Gammie (1996) comes from the attenuation length of comic rays ( 96 g cm; Umebayashi & Nakano (1981)). If the dust abundance is negligible, MRI is sustained in the layers as thick as the attenuation length (Sano et al., 2000). On the other hand, in the MHD simulations of Suzuki et al. (2010), they adopt the ionization degree calculated in Sano et al. (2000) and Inutsuka & Sano (2005) assuming that the dust-to-gas ratio be 0.01. In the present study, we vary as a parameter. On the other hand, we fix to be 25 g cm, as we find that is less important than for overall disk evolution, in particular, for . MHD simulations show that dead zones have non-zero viscosities (Fleming & Stone, 2003; Suzuki et al., 2010); the value of is on the order of -.

2.3 Disc temperature

As heat sources of discs, we take into account viscous heating, irradiation from the central star, and background irradiation. We simply superpose all contributions which are independently modeled:

 T4=T4vis+T4irr+T4amb, (8)

where , , and are the temperatures contributed from viscous heating, irradiation from the central star, and background irradiation. If only viscous heating is considered, the radiative cooling from both sides of the disc balances with the viscous heating as

 2σSBT4vis,eff=94ΣνΩ2, (9)

where is the effective disc surface temperature. The relation between the midplane temperature and the effective temperature is given by Hubeny (1990) as

 T4vis=τeffT4vis,eff, (10)

with

 τeff=[3τ8(1+ΣaΣ−αaΣaαΣ)+√34+13τ], (11)

where is the optical depth at the midplane. Eq. (11) is applicable to all cases. For pure MRI active regions without dead zones, we set and . For the dead zone region with , we set so that Eq. (10) becomes the same with Eq. (7) of Gammie (1996) for . A detailed derivation of the first term in the square bracket of the right hand side of Eq. (11) is given in Wünsch et al. (2006). The opacity is represented by , where is a constant. We use from Stepinski (1998) shown in Table 1. Results of the present study do not sensitively depend on opacity law.

The disc effective temperature irradiated by the central star is given by (Ruden & Pollack, 1991)

 T4irr,eff=T4∗[23π(r∗r)3+12(r∗r)2(hr)(dlnhdlnr−1)], (12)

where and are the photospheric temperature and radius of the central star. We adopt K and . The scale height used in Eq. (12) is the photospheric disc scale height and we simply assume it to be the same with the pressure scale height. After Hueso & Guillot (2005), we assume . This is validated as at large radii where irradiation from the central star is usually the dominant heat source. The midplane temperature, , irradiated by the central star is calculated using the vertical temperature structure derived by Malbet & Bertout (1991) (see also Malbet, Lachaume & Monin (2001)): {dmath} T_irr^4 = T_irr, eff^4[34μ_0(1-exp(-τμ0)) + 12 + 14μ0exp(-τμ0)]. Here is the averaged cosine of the incident angle of the central star given as

 μ0=−H0J0, (13)

with the zeroth and first order moments, and , of the incident intensity at the upper most layer:

 J0=σSBT4∗4π[1−(1−(r∗r)2)1/2], (14)
 H0=−σSBT4irr,eff4π. (15)

Eq. (2.3) represents that of an optically thick disc is lower than by a factor of , as the upper most layers emit half of the received energy toward the midplane. In the optically thin limit (), the third term in the bracket dominates and we obtain , which corresponds to the equilibrium temperature irradiated by the half hemisphere of the central star without any obstacle. The temperature at should be higher by a factor of , but this difference does not affect our results as the stage with a low surface density is very short with photoevaporative winds.

The ambient temperature is simply assumed as K.

2.4 Mass loss due to photoevaporation

Among various types of photoevaporations, the one induced by X-rays from the central star is most strong and most important for gas disc evolution in the planet formation region (Ercolano, Clarke & Drake, 2009; Owen et al., 2010, 2011), unless discs are close to very massive external stars (Adams et al., 2004; Mitchell & Stewart, 2010). The total mass loss rate from a disc due to X-ray photoevaporation is estimated by Owen et al. (2011) as

 ˙Mw=6.4×10−9A(LX1030ergs−1)1.14M⊙yr−1, (16)

where is a constant of the order of unity (we fix it to be unity) and is the X-ray luminosity from the central star. For sub-solar to solar mass protostars, is - erg s (Güdel et al., 2007). Numerical simulations of Owen et al. (2010) show that is roughly proportional to (their Fig. 13). Thus, we give the mass loss rate per unit area as

 ˙Σw=˙Mw4π(r1/2X,out−r1/2X,in)r−3/2(forrX,in≤r≤rX,out), (17)

where and are the inner and outer edges of the range where mass loss occurs. We take 1 AU and 70 AU from Owen et al. (2010). Outside , other mechanisms such as EUV or FUV photoevaporation (Gorti et al., 2009) may work, although these effects are not taken into account in the present work. Results of disc evolution do not sensitively depend on -dependence of , as similar results are obtained even if we vary the power-law index between -1 and -2.

In some test simulations, we also adopt the mass loss due to MRI disc winds (Suzuki et al., 2010). The mass loss rate of MRI disc winds at a given radius is proportional to the surface density. We find that gas depletes very rapidly even near the inner edge, and that the mass accretion rate onto the central star becomes less than yr within 1 Myr, as shown in Suzuki et al. (2010). This accretion rate looks too low as compared with those for typical discs around classical T Tauri stars (Hartmann et al., 1998). Probably, the mass loss rate of MRI disc winds is much smaller, but this mechanism may still play an important role in disc evolution.

2.5 Numerical procedures

We numerically solve Eq. (4) using the method described in Bath and Pringle (1981), in which the radial grid size is proportional to . We use 1000 grids between 0.03 AU and 3000 AU. These numbers give the inner most grid size of 0.02 AU. At the inner boundary is fixed to be zero, and at the outer boundary we adopt the outward mass flux given by . We find that the evolution of the total mass is nearly the same even with larger outer boundary radii whereas artificial mass loss from the outer boundary is not negligible with smaller outer boundary radii. The temperature is calculated by the Newton-Raphson method.

For all simulations, the mass of the central star is the solar mass, , the initial mass of the disc is 0.1 , and the initial surface density is proportional to with the outer edge of 20 AU. The adopted time step size is 1 yr. A single simulation for 10 Myr takes roughly one cpu day. We conduct 8 runs for discs with radially constant ’s and 18 runs for discs with dead zones. The input parameters are shown in Tables 2 and 3 for simulations of discs without and with dead zones.

3 Results

Results for all simulations are summarized in Table 4.

Figure 1 shows time evolution of the surface density , the temperature , and the mass accretion rate for a disc with a radially constant (run N2). As the disc spreads with time, and decrease. The kinks seen in the profile correspond to the opacity transition temperatures (150 K, 180 K, and 1380 K) and the corresponding kinks are also seen in . The direction of the radial gas motion is inward in the inner disc and outward in the outer disc, and the zero radial velocity radius expands with disc expansion. The accretion rate is independent of in the inner disc and is equivalent to the mass accretion rate onto the central star , as long as is sufficiently larger than the mass loss rate due to photoevaporation . When becomes less than , rapidly decreases and eventually a gap opens at slightly beyond 1 AU, as well as simulations in Owen et al. (2011). Once a gap opens, the inner disc (inside the gap) disperses very quickly and mass accretion onto the central star stops. The hole size is only 1.3 AU at this time. The temperature increases when the disc becomes optically thin, but this phase (with a non-zero gas density) is very short.

The gap opening slightly outside is explained as follows. Let us assume that the inward mass flux due to viscous accretion at is . As gas moves inward by , the mass flux is reduced by . Therefore, the inward mass flux decreases with decreasing and eventually becomes zero at if . Thus, if is independent of without photoevaporation, a gap inevitably opens near . Since it takes a time for gas to radially move from to (the viscous time scale at ), becomes much less than at the time of actual gap opening (see the black curve in Fig. 1).

In some runs, large gaps can open (runs N3 and N6; see Table 4). This happens because the disc sizes are not sufficiently larger than and the mass accretion rates are not radially constant near the outer edges of the discs. Even in these cases, the mass accretion rates onto the central stars at the time of gap opening are much smaller than those seen in transition discs. Overall, a large accretion rate and a large hole size are not simultaneously reproduced as long as is radially constant.

3.2.1 Cases with no dead zone viscosity (αd=0)

Pictures of evolution of discs with dead zones are very different from those without dead zones. Figure 2 shows evolution of a layered disc with (run D1). The mass accretion rate discontinuously drops at K, because the opacity law changes. As a result, two dead zone regions appear. Splitting into multiple dead zones was theoretically predicted in Gammie (1996). Most of the disc mass remains in the dead zones and remains nearly constant as long as the dead zones exist. Since in the dead zone region increases with , increases with time (Gammie, 1996; Zhu et al., 2010).

For run D1, a gap opens at 48 AU roughly when the outer dead zone disappears, because of the outer dead is larger than while in the inner dead zone is comparable to . Since the inner dead zone survives for a long time, remains high even after a gap opens. The outward mass flux from the dead zone remains to be after a gap opens, as long as the dead zone exists. Once the dead zone disappears, the inner disc quickly dissipates and mass accretion onto the central star is halted.

3.2.2 Cases with finite dead zone viscosity (αd>0)

It is possible to retain a steady state mass accretion with a finite residual viscosity in the dead zone. The condition to retain a radially constant is given from Eqs (3), (5), and (6) as

 r1/2∂∂r[[Σaαa+(Σ−Σa)αd]Tr2]=const. (18)

If the radial profile of is given such as irradiated discs, the above condition is fulfilled when so that the radial profile of is adjusted. On the other hand, if viscous heating is the dominant heat source, the midplane temperature is given from Eqs. (5), (6), (9), (10) as

 T=[27128Rκ0σSBμΩ((Σ2−Σ2a)αd+Σ2aαa)]1/(3−β). (19)

Thus, if , the radial profiles of and are mutually adjusted so that can be independent of .

Figure 3 shows an example of evolution of a disc with (run D6). Since this case satisfies the condition , becomes almost independent of (see the orange line at 2.0 Myr in Fig. 3). Split into multiple dead zones is also suppressed. As well as the cases with , a gap opens at a radius outside the dead zone.

Figure 4 shows the radial location of a gap for various values of , , and . We define as a radius where becomes zero in the earliest time while in outer radii. If is independent of in the dead zone region, the condition for gap opening is given by . For small values of and , can be smaller than even when a dead zone exists. Thus, a gap opens beyond a dead zone in such a case. On the other hand, for a case with large and , a gap opens at a small radius only after a dead zone disappears and becomes sufficiently small, as well as the cases with radially constant ’s shown in Section 3.1. Not surprisingly, gap opening beyond dead zones is possible even with large values of and if (or ) is large.

4 Discussion

Figure 5 shows comparison between observed transition discs and discs from our simulations. In the upper panel of Fig. 5, gap sizes and mass accretion rates are plotted. For simulations, we use and the mass accretion rate onto the central star at the time of gap opening (see Table 4). Subsequent evolution curves are also shown. As can be seen, large gap (or hole) sizes and large accretion rates seen in observed transition discs are well reproduced in models with dead zones. On the other hand, models without dead zones are able to reproduce only discs with small mass accretion rates and small gap sizes. In runs N3 and N6, a gap opens at 60 AU, but the mass accretion rate is small (in these runs, another gap opens at 1AU, and the disk outside the outer gap quickly dissipates).

In the lower panel of Fig. 5, we compare between disc masses from observations and our simulations. For our simulations, we take the disc mass outside the gap ( in Table 4), not the total disc mass (), for consistency with observations. It should be noted that submillimeter - millimeter observations measure emission from dust primarily around 100 AU and the surface density of gas is estimated using the interstellar gas-to-dust ratio. We find that our model cannot reproduce the masses of the most massive transition discs with high accretion rates, although the total disc masses including the inner discs in our simulations are comparable to the observed values (see Table 4). This means that most of the mass is in the dead zone for a disc with layered accretion. In our model of layered accretion discs, at 100 AU is 1 g cm before a gap opens but it decreases to less than 0.1 g cm after a gap opens (see Figs. 2 and 3). On the other hand, observed transition discs show that 1-10 g cm at 100 AU (Andrews et al., 2011). This value is similar to or even larger than those seen in classical T Tauri discs (Kitamura et al., 2002; Andrews et al., 2010). The apparent contradiction between modeled and observed disc masses might be resolved if dust is accumulated near the inner edge of the outer disc after a gap opens, resulting in a large dust-to-gas ratio (Alexander & Armitage, 2007). The dust-to-gas ratio may also increase if photoevaporative winds remove only gas but not dust.

The expected fraction of gapped discs is much larger for discs with dead zones than that for discs without dead zones (Fig. 6), if the disc life times are similar. The fraction of transition discs in protoplanetary discs is more than 50 per cent at 5-8 Myr (Currie & Sicilia-Aguilar, 2011), and a large fraction of them are accreting (Merín et al., 2010). The fraction of transition discs with large holes ( AU) is at least 20 per cent among millimeter-bright disc population (Andrews et al., 2011) although discs with large holes seem to be rarer than those with small holes (Merín et al., 2010). For discs with dead zones, it is possible to reproduce a high fraction of discs with large holes, if and are small and is large. A detailed comparison may be able to constrain and with a given distribution of .

Dead zones in our model are optically thick even after gaps open in discs. This is inconsistent with observations, although some of transition discs have inner optically thick regions indicated from near-infrared excess (pre-transition discs; Espaillat et al. (2010)). This indicates that dust removal due to their growth or migration is necessary whereas our current model adopts a fixed opacity. Since dust is usually supplied from the outer discs (Brauer, Dullemond & Henning 2008), gap opening may make the inner discs or dead zones optically thin. Alternatively, giant planets may clean up dust in the inner discs (Merín et al., 2010; Zhu et al., 2011). Giant planets may also be responsible for transition discs with small holes and large accretion rates, as our model does not reproduce such discs (Fig. 5).

5 Summary

In this study, we developed a gas disc model which takes into account layered accretion and photoevaporative winds induced by X-rays from the central stars. We found that a gap opens at a radius outside a poorly-ionized dead zone, if the mass loss rate due to photoevaporation exceeds the mass accretion rate in the dead zone region. Since the dead zone survives even after the gap opens, high mass accretion onto the central star remains for a long time.

We found good agreements between modeled and observed transition discs in regards to gap sizes and mass accretion rates. However, our model shows disc masses (we take masses outside gaps) an order of magnitude smaller than those for the most massive observed transition discs. This may indicate that the dust-to-gas ratios are large in the outer discs of observed transition discs while our model assumes a fixed dust opacity.

Acknowledgments

We are grateful to anonymous reviewers for fruitful comments on our manuscript.

References

• Adams et al. (2004) Adams F. C., Hollenbach D. Laughlin G., Gorti U., 2004, ApJ, 611, 360.
• Alexander & Armitage (2007) Alexander R. D., Armitage P. J., 2007, MNRAS, 375, 500.
• Andrews et al. (2011) Andrews S. M., Wilner D. J., Espaillat C., Hughes A. M., Dullemond C. P., McClure M. K., Qi C., Brown J. M., 2011, ApJ, 732, 42.
• Andrews et al. (2010) Andrews S. M., Wilner D. J., Hughes A. M., Qi C., Dullemond C. P., 2010, ApJ, 723, 1241.
• Bai (2011) Bai X.-N., 2011, ApJ, in press (arXiv:1107.2935).
• Balbus & Howley (1991) Balbus S. A., Hawley J. F., 1991, ApJ, 376, 214.
• Brauer et al. (2008) Brauer F., Dullemond C. P., Henning T., 2008, A&A, 480, 859.
• Currie & Sicilia-Aguilar (2011) Currie T., Sicilia-Aguilar A., 2011, ApJ. 732, 24.
• Dullemond & Dominik (2005) Dullemond C. P., Dominik C., 2005, A&A, 434, 971.
• Ercolano, Clarke & Drake (2009) Ercolano B., Clarke C. J., Drake J. J., 2009, ApJ, 699, 1639.
• Espaillat et al. (2008) Espaillat C. et al., 2008, ApJ, 689, L145.
• Espaillat et al. (2010) Espaillat C. et al., 2010, ApJ, 717, 441.
• Fedele et al. (2010) Fedele D., Ancker M. E., Henning T., Jayawardhana R., Oliveira J. M., 2010, A&A, 510, A72.
• Fleming & Stone (2003) Fleming T., Stone J. M., 2003, ApJ, 585, 908.
• Gammie (1996) Gammie C. F., 1996, A&A, 457, 355.
• Gorti et al. (2009) Gorti U., Dullemond C. P., Hollenbach D., 2009, ApJ, 705, 1237.
• Güdel et al. (2007) Güdel M. et al., 2007, A&A, 468, 353.
• Hartmann et al. (1998) Hartmann L., Calvet N., Gullbring E., D’Alessio P., 1998, ApJ, 495, 385.
• Hubeny (1990) Hubeny I., 1990, ApJ, 351, 632.
• Hueso & Guillot (2005) Hueso R., Guillot T., 2005, A&A, 442, 703.
• Inutsuka & Sano (2005) Inutsuka S., Sano T., 2005, ApJ, 628, L155.
• Kim et al. (2009) Kim K. H. et al., 2009, ApJ, 700, 1017.
• Kitamura et al. (2002) Kitamura Y., Momose M., Yokogawa S., Kawabe R., Tamura M., Ida S., 2002, ApJ, 581, 357.
• Lynden-Bell & Pringle (1974) Lynden-Bell D., Pringle J. E., 1974, MNRAS, 168, 603.
• Malbet & Bertout (1991) Malbet F., Bertout C., 1991, ApJ. 383, 814.
• Malbet et al. (2001) Malbet F., Lachaume E., Monin J. -L., 2001, A&A, 379, 515.
• Merín et al. (2010) Merín B. et al., 2010, ApJ, 718, 1200.
• Mitchell & Stewart (2010) Mitchell T., Stewart G., 2010, ApJ, 722, 1115.
• Muzerolle et al. (2010) Muzerolle J., Alllen L., Megeath T. M., Hernández J., Gutermuth R. A., 2010, ApJ, 708, 1107.
• Najita, Strom & Muzerolle (2007) Najita J.R., Strom S. E., Muzerolle J., 2007, MNRAS, 378, 369.
• Owen et al. (2010) Owen J. E., Ercolano B., Clarke C. J., Alexander R. D, 2010, MNRAS, 401, 1415.
• Owen et al. (2011) Owen J. E., Ercolano B., Clarke C.J., 2011, MNRAS, 412, 13.
• Pollack et al. (1994) Pollack J. B., Hollenbach D., Beckwith S., Simonelli D. P., Roush T., Fong W., 1994, ApJ, 421, 615.
• Ruden & Pollack (1991) Ruden S. P., Pollack J. B., 1991, ApJ, 375, 740.
• Sano et al. (2000) Sano T., Miyama S. M., Umebayashi T., Nakano T., 2000, ApJ, 543, 486.
• Shakura & Sunyaev (1973) Shakura N.I., Sunyaev R.A., 1973, A&A, 24, 337.
• Stepinski (1998) Stepinski T. F., 1998, Icarus, 132, 100.
• Suzuki et al. (2010) Suzuki T. K., Muto T., Inutsuka S., 2010, ApJ, 718, 1289.
• Tanaka et al. (2005) Tanaka H., Himeno Y., Ida S., 2005, ApJ, 625, 414.
• Thalmann et al. (2010) Thalmann C. et al., 2010, ApJL, 718, 87.
• Umebayashi & Nakano (1981) Umebayashi T., Nakano T., 1981, PASJ, 33, 617.
• Williams & Cieza (2011) Williams J. P., Cieza L. A., 2011, ARA&A, in press.
• Wünsch et al. (2006) Wünsch R., Gawryszczak A., Klahr K., Różyczka M., 2006, MNRAS, 367, 773.
• Zhu et al. (2010) Zhu Z., Hartmann L., Gammie C., 2010, ApJ, 713, 1143.
• Zhu et al. (2011) Zhu Z., Nelson R. P., Hartmann L., Espaillat C., Calvet N., 2011, ApJ, 729, 47.
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters