In population synthesis calculations (e.g., Ida & Lin 2004; Mordasini et al. 2009a), one uses distributions of (initial) protoplanetary disk masses in order to reflect the varying initial conditions for planet formation.
As discussed in Sect. 3.2, massive cold disks get self-gravitationally unstable, which would invalidate the usage of a classical α model with one constant α across the disk. Ida & Lin (2004) and Mordasini et al. (2009a) have therefore cut off the uppermost part of the observed disk mass distribution with the argument that they cannot be stable. The cut off was usually done at the (often assumed) disk mass stability limit of about a tenth of the stellar mass, without a direct calculation of this stability limit. Observationally, high disk masses can be mimicked by residual dust in the remains of the envelope from which the star formed, contributing to the observed flux (Andrews & Williams 2005).
In this Appendix, we use our upgraded disk model to determine the maximum stable disk mass for both irradiated and non-irradiated disk and make some more general remarks about the stability of irradiated α disks against the development of spiral waves and clumping.
For simplicity, we use for the initial disk profile (Eq. (7)) a power-law exponent γ = 0.9 and a characteristic radius Rc = 30 AU. Both values correspond to the maxima of the distributions observed by Andrews et al. (2010). With these values, the initial gas surface density is (A.1)It is clear that this procedure faces difficulty because it implicitly assumes that we can use the observed dust grains to trace the mass and evolution of the gas. Especially grain growth could make this assumption partially invalid (e.g., Andrews & Williams 2007).
Toomre parameter QToomre (red solid line) and τcoolΩ (blue dotted line) as a function of semimajor axis a at a moment shortly after the beginning of disk evolution. The horizontal lines correspond to the critical values of QToomre,crit = 1.7 and τcoolΩ = 3, respectively. The left panel is for a disk with an initial mass Md(t = 0) = 0.024 M⊙, while the right one corresponds to a higher mass, Md(t = 0) = 0.11 M⊙.
|Open with DEXTER|
When integrating Eq. (7) from r = 0 to infinity using the parameters mentioned before, one finds a total initial disk mass (cf. Miguel et al. 2011) of (A.2)Thus, the MMSN of Weidenschilling (1977, 2005) corresponds to Σ0 ≈ 200 g/cm2, while for Hayashi’s (1981) MMSN of 0.013 M⊙, Σ0 ≈ 108 g/cm2. It is likely that the true initial mass of the solar nebula was a few times larger than these minimal masses (Weidenschilling 2005).
We note that we do not allow in the model accretion rates larger than 3 × 10-7 M⊙/yr in order to avoid convection in the vertical direction, but if necessary reduce locally the initial Σ in order not to exceed this limit. Therefore, the initial gas mass for disks with Σ0 larger than ~1000 g/cm2 is somewhat smaller than predicted by Eq. (A.2). For Σ0 = 2000 g/cm2, the initial mass is, for example, about 20% smaller than given by Eq. (A.2). We note further that in contrast to earlier models, we can now include the total disk mass inside the computational domain and not only the fraction contained in the innermost 30 AU (Mordasini et al. 2009a).
Figure A.1 shows QToomre (solid red lines) and τcoolΩ (blue dotted lines) as a function of distance from the star for two disks (see also Bell et al. 1997) at t ≈ 0 (the moment we start the disk evolution). The left-hand panel shows a 1 − 2 × MMSN disk (Md(t = 0) = 0.024 M⊙, Σ0 = 200 g/cm2), while the right one shows a more massive disk with Md(t = 0) = 0.11 M⊙ (Σ0 = 1000 g/cm2). Both examples are calculated including the effect of stellar irradiation on the disk temperature structure, M∗ = 1 M⊙, α = 7 × 10-3 and the initial surface density profile discussed above.
We see that for both disks, QToomre reaches a minimum at about 25 AU. For the MMSN disk, the minimum value is about 12, while for the more massive disk, it is 2.5. The critical value of QToomre,crit ≈ 1.7 is not reached. Regarding the cooling, the particular structure of the curve is due to opacity transitions. We see that in the inner parts of the disk and outside 15 − 25 AU, there are parts of the disk that could in principle cool quickly enough, so that τcoolΩ < 3. The cooling timescale is the time in which the disk would cool if we suddenly switched off the heating mechanisms. This is however meaningless for gravitational fragmentation due to the too high values of Q everywhere in the disks. Therefore, fragmentation does not occur, either.
We have studied the overall minimal QToomre,min occurring in disks as a function of time, radius, and Σ0 in order to determine the maximum Σ0 that can be used. In these calculations, we again assume α = 7 × 10-3 and M∗ = 1 M⊙. It is found that the minimal QToomre always occur very shortly after the start of the disk evolution (within some 104 years) and that afterwards QToomre always increase. Thus, the minimal QToomre are a direct consequence of the initial conditions (as expected), and the disks evolve towards stability.
The orbital distance where the overall minimal QToomre,min occur correspond to about 20 AU for disks calculated without stellar irradiation and, as also seen in the two examples in Fig. A.1, to about 25 AU for disks with irradiation, independent of the disk mass and as long as Md(t = 0) ≲ 0.1 − 0.2 M⊙. For more massive disks, the distance of the minimal Q increases. This result can be understood at least for the irradiated disks if we assume that in the outer parts the temperature structure is given by the stellar irradiation only. We can consider two cases: first, an optically thin disk with Tmid ∝ r − 1/2. In this case, the orbital distance RToomre,min,t where the minimal QToomre,min occurs is given for γ < 7/4 as (using Eqs. (7), (10) and Tmid ∝ r − 1/2) (A.3)This corresponds to ratios RToomre,min,t/Rc = 1/4, 3/4, 0.791, 0.886 for a disk surface density exponent (Eq. (7)) γ = 3/2, 1, 0.9, and 1/2.
Second, we can consider a passively irradiated disk without viscous dissipation, where for orbital distances between about 0.4 and 84 AU, Tmid ∝ r − 3/7 (Chiang & Goldreich 1997). In this case, one finds an orbital distance RToomre,min,p where the minimal QToomre,min occurs at (A.4)This corresponds to ratios RToomre,min,p/Rc = 0.18, 0.71, 0.76, 0.87 for γ = 3/2, 1, 0.9, and 1/2. We see that the results for an optically thin and a passively irradiated disk are similar. This is due to the fact that the temperature depends in both cases in a similar way on the distance.
For Rc = 30 AU and γ = 0.9, this leads for a passively irradiated disk (which is the situation assumed here, see Fouchet et al. 2012) to a RToomre,min,p = 22.8 AU, close to the result ( ≈ 25 AU) seen in the simulations. The somewhat larger value in the simulations is likely due to residual viscous heating.
In Fig. A.2 we plot the overall minimum QToomre,min (i.e., the lowest value reached at any distance and time) as a function of Σ0 (or equivalently the initial disk mass), for disks with (red solid) and without (blue dotted) irradiation. In both cases, Q follows a simple power-law that scales as , as expected from Eq. (10).
One sees that Q is higher in disks with irradiation, which is expected as these disks are hotter and thus more stable, in particular in the outer regions, where Q becomes small and where viscous heating is not important.
Overall minimal Toomre parameter QToomre,min as a function of initial gas surface density at 5.2 AU, Σ0 which is a direct proxy of the initial disk mass (Eq. (A.2)). The red solid line is for disks with irradiation, the dotted blue one is for disks where viscosity is the only heating source. The dashed line gives the critical value of 1.7.
|Open with DEXTER|
From Fig. A.2 we find a maximum initial gas surface density Σ0,max (measured at 5.2 AU) leading to stable disks (QToomre,crit = 1.7) for disks without irradiation of 790 g/cm2, corresponding to an initial disk mass of 0.091 M⊙, while for disks with irradiation, the maximum allowed value is 1510 g/cm2, corresponding to an initial mass of 0.16 M⊙.
Gas surface density Σ in a protoplanetary disk as a function of distance and time. The surface density is plotted in intervals of 2 × 104 years. The uppermost line shows a state shortly after the beginning of disk evolution. The lowermost line is the profile when the calculation is stopped. The two panels differ from each other by the boundary condition that is used at the inner edge of the disk at 0.1 AU, as described in the text. All other parameters are identical and given as M∗ = 1 M⊙, α = 7 × 10-3 and Σ0 = 200 g/cm2.
|Open with DEXTER|
At disk masses higher than ~10 − 25% the mass of the central object, the possibility of global gravitational instability must be considered, too (e.g., Harsono et al. 2011).
The particular case of global instability to m = 1 modes (one-armed spirals) has been studied numerically and analytically in Adams et al. (1989) and Shu et al. (1990). If this mode is amplified by the SLING mechanism (which is due to the fact that for a one-armed spiral, the star is displaced through the conservation of the center of mass), there is a finite threshold for the instability to set in, which is a function of the Toomre criterion. For the case that QToomre = 1 at the outer disk edge, it is possible to derive analytically a critical disk mass, below which the disk is gravitationally stable against all modes. In a linear stability analysis, Shu et al. (1990) find for this situation a critical “maximum-mass solar nebula” at Mdisk/(Mdisk + M∗) = 3/(4π), corresponding to Mdisk/M∗ = 0.31.
The analytical analysis of Shu et al. (1990) suffers from a number of limitations (cf. Noh et al. 1992 for a discussion), and might not be applicable in the strict sense in the context here. There are two reasons for this: Our disks have a smooth outer edge, so that they might not reflect the density waves, and in addition, for our disks, the minimal Q does not occur at the outer disk edge, but further in (Eq. (A.4)). If we compare the critical disk mass of Shu et al. (1990), Mdisk/M∗ = 0.31, we find, despite the limitations with the most massive disk stable to local disturbances for the model here (Mdisk/M∗ = 0.16), we find that the disks should also be globally stable.
Since the works of Adams et al. (1989) and Shu et al. (1990), many studies have revisited the question of (global) disk stability. While a direct connection to the SLING mechanism can in general not be made (Nelson et al. 1998), it is nevertheless found that at constant Q, with increasing Mdisk/M∗, the character of the instability changes. At low Mdisk, thin, multi-armed structures develop that are characterized by high-order patterns (m ≳ 5), whereas at high Mdisk, global, low order (m = 2 − 5) instabilities dominate (Nelson et al. 1998). The transition between the two regimes occurs at approximately Mdisk/M∗ = 0.2 − 0.4. These results are confirmed at much higher numerical resolution by Harsono et al. (2011) and Lodato & Rice (2004). The latter authors find that local effects dominate as long as the disk mass is less than 0.25 M∗ and the disk aspect ratio H/R ≲ 0.1. Our massive disks are characterized by H/R between 0.06 to 0.09 in their outer parts, so that they fulfill this criterion.
These results indicate that the transition from local to global instabilities might rather be a gradual one and not set in at a specific disk mass, as originally advocated by Shu et al. (1990). The disk mass where global effects become important seems to be at about 0.25 M∗ (Nelson et al. 1998; Lodato & Rice 2004; Harsono et al. 2011), which is comparable to the original criterion by Shu et al. (1990). These results indicate that our most massive, Toomre-stable disk with 0.16 M⊙ should also be stable to global modes. We must, however, keep in mind that the numerical simulations mentioned here did not use exactly the same surface density and temperature profile. Therefore, in order to get a firm conclusion, dedicated hydrodynamic simulations are necessary.
In this Appendix, we illustrate the characteristic evolution of the protoplanetary disks under the combined action of viscosity and photoevaporation, as found in the update model, in order to check our results against previous studies. Figure B.1 shows the gas surface density as a function of time for a disk with α = 7 × 10-3, Σ0 = 200 g/cm2, and a mean external photoevaporation rate over the lifetime of the disk of about 7 × 10-9 M⊙/yr. The effects of irradiation on the thermal structure of the disk is taken into account as described in Fouchet et al. (2012). These initial conditions result in a disk lifetime of 2.0 Myr (remaining disk mass 10-5 M⊙). In the figure, a line is plotted all 2 × 104 years, and the minimal allowed surface density is set to 10-3 g/cm2.
The two panels of the figure differ only by the inner boundary condition at Rmin = 0.1 AU. In the left-hand panel, we use the same boundary conditions as in our earlier models, which means that the flux through the innermost cell instantaneously adopts to the flux coming from further out. This is equivalent to the statement that the disk structure would continue to some smaller radius inside 0.1 AU, which is not modeled. In the right-hand panel, Σ is forced to fall to zero at 0.1 AU. Physically, we can associate the two situations with different sizes of the magnetospheric cavity, i.e., with a weak and a strong magnetic field of the host star, respectively (e.g., Bouvier et al. 2007). A detailed description of the structure of the disk close to the host star that takes into account magnetic fields will be presented in an upcoming work (Cabral et al., in prep.).
While the inner boundary condition has a very strong impact on the migration behavior of low-mass planets close to the star (Benitez-Llambay et al. 2011), it has otherwise no effect on the characteristic evolution of the disk.
Figure B.1 shows that the evolution of the disk can approximatively be separated in four phases:
The initial exterior radius specified by the initial conditions getsquickly reduced by external photoevaporation to a radius wheremass removal due to photoevaporation, and the viscousspreading of the disk are in a quasi-equilibrium, as described byAdams et al. (2004). In the specific case, the radiusdecreases from initially about 200 AU to ~100 AU. In the inner part, the disk very quickly evolves from the initial profile towards equilibrium.
In the second, dominant phase a quasi self-similar evolution of the disk occurs. The inner part of the disk (r ≲ 10 AU) is in near equilibrium (i.e., the mass accretion rate is nearly constant as a function of radius), and the slope of the gas surface density γ is approximately −0.9 (but varying between −0.4 and −1.5 due to opacity transitions). The outer radius is slowly moving inwards, from about 100 AU to 60 AU. For more massive disk, this equilibrium radius is further out.
Once the surface density has fallen to about 0.01−0.1 g/cm2 at ~10 AU, a gap opens somewhat outside of Rg,II. The evolution of the disk now speeds up, which corresponds to the so-called “two-timescale” behavior (Clarke et al. 2001).
Quickly afterwards, the total disk mass has fallen to 10-5 M⊙, where we stop the calculation. We note that the evolution at very small Σ is not important for our purpose of planet-formation modeling. Therefore, we currently do not include the effect of the direct radiation field, which would clear the disk quickly from inside out once the gap has opened (Alexander & Armitage 2009).
Such an evolution is very similar to the findings of previous studies (Matsuyama et al. 2003; Clarke et al. 2001). We note that the evolution can become more complex if there is additionally a planet accreting significant amounts of gas. This is the case if the planetary core becomes massive enough to trigger gas runaway accretion, so that a giant planet forms. Such a planet then effectively acts as a sink cell in the disk (see Paper I), so that another gap would form at its position towards the end of disk evolution, even if we neglect the tidal gap formation.
With the inclusion of a both detailed model for the photoevaporation and the calculation of the luminosity of the planets in the gas runaway accretion phase (where forming giant planets can be bright, cf. Paper I), we are able to address new observational constraints. As shown by Fouchet et al. (2012), we can use the disk structure to calculate the spectral energy distribution, in which we can now include the contributions from the star, the disk, and the growing planet. As with an accreting star, we expect two contributions from a giant planet undergoing rapid gas accretion: a contribution in the infrared coming from the internal luminosity, and a hard component from the accretion shock. This will be addressed in a dedicated study (Mordasini et al., in prep.). With upcoming observational facilities like the Atacama Large Millimeter/submillimeter Array (ALMA), observing planet formation as it happens (Wolf et al. 2002; Klahr & Kley 2006) will become possible and put a whole new class of constraints on formation models.
© ESO, 2012