A&A 417, 151-158 (2004)
DOI: 10.1051/0004-6361:20034036

An alternative look at the snowline in protoplanetary disks

K. Kornet1 - M. Rózyczka1 - T. F. Stepinski2


1 - Nicolaus Copernicus Astronomical Center, Bartycka 18, Warsaw, 00-716, Poland
2 - Lunar and Planetary Institute, 3600 Bay Area Blvd., Houston, TX 77058, USA

Received 2 July 2003 / Accepted 12 November 2003

Abstract
We have calculated an evolution of protoplanetary disk from an extensive set of initial conditions using a time-dependent model capable of simultaneously keeping track of the global evolution of gas and water-ice. A number of simplifications and idealizations allows for an embodiment of gas-particle coupling, coagulation, sedimentation, and evaporation/condensation processes. We have shown that, when the evolution of ice is explicitly included, the location of the snowline has to be calculated directly as the inner edge of the region where ice is present and not as the radius where disk's temperature equals the evaporation temperature of water-ice. The final location of the snowline is set by an interplay between all involved processes and is farther from the star than implied by the location of the evaporation temperature radius. The evolution process naturally leads to an order of magnitude enhancement in surface density of icy material.

Key words: accretion, accretion disks - solar system: formation

1 Introduction

Water-ice in the protoplanetary disk can exist beyond a minimum distance from the star called the snowline. Although the snowline is defined directly by the presence of ice, the standard method of calculating its location is indirect and uses disks temperatures. Thus, in the present literature, the snowline is the location where the disks temperature is equal to the sublimation/condensation temperature of water-ice. The total surface density of all solids,  $\Sigma_{\rm s}$, increases rapidly outside the snowline because water-ice, the most abundant species of solid, becomes available and its contribution dominates the value of $\Sigma_{\rm s}$. In the core accretion - gas capture scenario of giant planet formation high values of  $\Sigma_{\rm s}$ are necessary to produce solid cores on time scales consistent with the presence of a gaseous nebula. Thus, at least in such a scenario, the importance of the snowline derives from a notion that it marks the inner edge of the giant planet formation zone.

In the context of the Solar System formation, the so-called "minimum-mass solar nebula'' model has been used frequently for quantitative analysis. In this steady-state model the temperature is derived from radiative equilibrium with the solar radiation field (Hayashi 1981) and the snowline is located at a distance  $R_{\rm sl}=2.7$ AU from the Sun. Just the presence of ice exterior to  $R_{\rm sl}$ increases solid abundance by a factor of 4 yielding $\Sigma_{\rm s}=6.8$ g cm-2 at the snowline and $\Sigma_{\rm s}=2.7$ g cm-2 at r=5 AU (present, and, by assumption, the original location of Jupiter). However, even with addition of ice, the values of  $\Sigma_{\rm s}$ are still too low for a rapid formation of giant planets cores. Lissauer (1987) calculated that surface density of solids 5-10 times greater than this given by the Hayashi model are required to grow the Jovian core fast enough to accrete the gaseous envelope before the solar nebula is dispersed. Stevenson & Lunine (1988) proposed a mechanism for a further enhancement of the abundance of water-ice just exterior to the snowline by diffusive redistribution of water vapor through the snowline. They calculated that such a mechanism can enhance the value of  $\Sigma_{\rm s}$ by 1 to 2 orders of magnitude, making the rapid formation of a Jovian core possible at the location of the snowline. However, diffusive redistribution cannot explain core formation for other giant planets in the Solar System.

The discovery of $\sim$100 giant extrasolar planets (for a recent review see Bodenheimer & Lin 2002) posed formidable theoretical problems for the traditional scenario of giant planet formation because most exoplanets were found interior to the snowline as calculated from the "minimum mass'' model, and a significant fraction of those (the so-called hot Jupiters) - closer than 0.1 AU, where both high temperature and insufficient amount of matter should have prevented their formation. It appears that post-formation orbital evolution such as migration of a planet under the influence of tidal torque from the disk matter (see e.g. Lin et al. 2000) or gravitational close encounters between planets leading to significant changes of their orbits (see e.g. Rasio & Ford 1996) need to be invoked to explain the orbital properties of hot Jupiters. This notwithstanding, a revisiting of the snowline concept is also needed to understand the original formation of giant exoplanets, and to explain orbital properties of those exoplanets that seem to escape a significant post-formation orbital evolution.

Two important developments call for the re-evaluation of the snowline concept. First, it becomes clear from astronomical observations (for a compact review see Beckwith 1994) that protoplanetary disks, of which a solar nebula is thought to be a particular example, are active, ever-changing entities with limited life-spans that cannot be modeled successfully by a steady-state, phenomenological model like the "minimum mass'' model. Instead, numerical and analytical models based on the concept of time-dependent accretion disk theory have been developed (Ruden & Pollack 1991; Sterzik & Morfill 1994; Stepinski 1998). In such models different initial conditions lead to disks of different character and evolutionary history, offering the potential to explain observed diversity in exoplanets properties. In an evolving disk, the location of the sublimation/condensation radius changes, and it has to be considered as a function of time  $R_{\rm evap}(t)$. This function is further parameterized by initial conditions.

Second, it has been recognized that the solid particles evolve differently to the gas (see e.g. Weidenschilling & Cuzzi 1993 and references therein). The differences between global evolution of gaseous and solid components of the protoplanetary disk have been calculated by Stepinski & Valageas (1996, 1997). They have shown that in the evolving protoplanetary disk, the sub-disk of water-ice solids decouples over time from the gaseous disk leading to significant departures from a constant solids-to-gas mass ratio The ice sub-disk evolves towards a swarm of icy planetesimals characterized by the value of  $\Sigma_{\rm s}$considerably higher than at the beginning of the evolution. Within such a model the standard method of calculating the location of the snowline is not useful. Because the gas and the solids decouple, the snowline location cannot be based on the temperature of the gas. Instead, the snowline  $R_{\rm sl}(t)$ should be calculated directly from its original definition as the inner edge of the region where water-ice is present.

Recently, the location of the snowline has been recalculated. Sasselov & Lecar (2000) considered passive and low accretion rate models of a protoplanetary disk instead of the "minimum mass'' model to calculate the location of the snowline using a standard technique of equating  $R_{\rm sl}$ with  $R_{\rm evap}$. They obtained $R_{\rm sl}\approx 1$ AU, significantly closer to a star than a "distant'' snowline at 2.7 AU as predicted by the "minimum mass'' model. The focus of the Sasselov & Lecar work was on using a model of the disk that is more consistent with observations than the "minimum mass'' model. Podolak (2003) used a steady-state model of an accreting disk to study how the location of the snowline depends on accretion rate, ice grain size, and contamination of ice by other materials. He defines the snowline as the location where the rate of ice grain evaporation equals the rate of ice condensation from water vapor. This is a more accurate, but still an indirect way of calculating the location of the snowline. Podolak finds that in the midplane, where most grains are to be found, the snowline is nearly independent of grain size and composition, but dependent on accretion rate.

These efforts do not address evolutionary character of protoplanetary disk, and, in particular, the decoupling of solids from the gas. However, the work by Stepinski & Valageas (1997), recently confirmed by more computationally sophisticated calculations by Weidenschilling (2003) indicates that the redistribution of solids is a very important process in the evolution of the disk. Due to coagulation and sedimentation the solid particles grow on a time scale that is short in comparison to the lifetime of the disk. As they grow, they start to drift inward, toward the star, and as they grow to $\sim$1 km sizes they settle into a fixed Keplerian orbit. This ability of the solids to stop migrating inwards allows for a seemingly unlikely situation in which  $R_{\rm sl}$ does not coincide with  $R_{\rm evap}$. The snowline has to be calculated directly from its definition by keeping track of an evolving distribution of water-ice. Calculating the snowline indirectly from gas properties would result in wrongly placing the snowline too close to the star.

In the present communication we demonstrate the difference between  $R_{\rm sl}$ and  $R_{\rm evap}$ for a broad sample of disk models analyzed by Kornet et al. (2001, 2002). We identify the radial drift as the major factor responsible for the final location of the snowline, and we show how redistribution of solids can enhance the abundance of solid material leading naturally to an emergence of a giant planet formation zone. The basic assumptions of our disk model and the basic numerical methods employed to solve the equations governing its evolution are briefly introduced in Sect. 2. The results of calculations are presented and analyzed in Sect. 3, while in Sect. 4 we discuss our approach and results in relation to ideas and results reported in the literature.

2 Method of calculation

The protoplanetary disk is modeled as a two-component fluid consisting of gas and solids. The evolution of the gas component is described by an analytic solution to the viscous diffusion equation, which gives the surface density of the gas as a function of radius r and time t (Stepinski 1998). The viscosity is given by the usual $\alpha $ model. The temperature of the gas is calculated in the thin-disk approximation, assuming vertical thermal balance, according to Eqs. (2) through (6) in Stepinski (1998).

Initially the solids are in the form of grains, but because in our model the solids grow all the way to planetesimals our calculations handle solid objects with a wide range of sizes. We refer to them collectively as "dust'', although term "solid particles'' and "dust particles'' are also used. The model of their evolution includes gas drag effect, sedimentation, coagulation and evaporation. Only one component of dust is considered, in this paper corresponding to water-ice, which has a sublimation temperature  $T_{\rm evap} = 150$ K and a bulk density $\varrho = 1$ g cm-3. The main assumptions used in the calculation are (1) at each radius the particles are all assumed to have the same size (which, of course, varies in time), (2) all collisions between particles lead to coagulation, (3) in disk regions with temperature exceeding  $T_{\rm evap}$ all water is in the form of vapour and evolves at the same radial velocity as the gas component, (4) initially, in disk regions with temperature below  $T_{\rm evap}$, the particles have sizes equal to  $a_{\rm min}=10^{-3}$ cm, (5) the systematic radial velocity of grains is entirely determined by the effects of gas drag, (6) the evolution of solids does not affect the evolution of the gaseous disk The evolution of solids is governed by the set of two equations (see Kornet et al. 2001). The first of them is the standard continuity equation for  $\Sigma_{\rm s}$. The second equation can be interpreted as the continuity equation for size-weighted surface density of solids  $\Sigma_a(r)\equiv a(r)\Sigma_{\rm s}(r)$ with the source term accounting for the growth of solid particles due to collisional aggregation. The radius of a solid particle at a distance r from the star is denoted by a(r). The set is solved numerically. Details of our numerical method is given in Kornet et al. (2001).

The initial conditions are parameterized by the quantities m0 (the mass of the disk in units of solar mass $M_\odot$), and j0 (the total angular momentum of the disk in units of 1052 g cm2 s-1). Once m0 and j0 are specified, the analytic solution of Stepinski (1998) gives the gas surface density $\Sigma_{\rm g,0}(r)=\Sigma_{\rm g}(r,t=0)$ that serves as an intial condition for the evolution of the gaseous disk.

At t=0 the ratio $\Sigma_{\rm s}/\Sigma_{\rm g}$ and the particle radius a are independent of r (values of 0.01 and  $a_{\rm min}=10^{-3}$ cm are adopted, respectively). Thus, the initial condition for the evolution of solid disk is $\Sigma_{\rm s,0}(r)=0.01\Sigma_{\rm g,0}(r)$. We have checked that the results are not sensitive to a particular selection of  $a_{\rm min}$ by recalculating some models with  $a_{\rm min} =
10^{-4}$ cm.


  \begin{figure}
\par\includegraphics[width=7.5cm,clip]{0036fig1.ps}
\end{figure} Figure 1: The evolution of solids in an  $\alpha = 0.1$ disk with  m0 = 0.02 and j0 = 1. Top: surface density of solids at the times indicated in the frame (solids that have sublimated but not yet accreted onto the star are not taken into account). Bottom: locations of sublimation limit (dashed) and snowline (solid) as a function of time The snowline is defined as the inner edge of the region where ice grains are present.
Open with DEXTER

3 Results

Our results are divided into three grids of disk models. Each grid groups models with the same constant value of the viscosity parameter $\alpha $ and contains complete evolutionary information for disks starting from 99 different initial conditions specified by  $0.02 \le
m_0 \le 0.2~M_\odot$ and  $1\le j_0 \le 25$. The grids are characterized by  $\alpha=10^{- 3}$, 10-2 and 10-1, respectively. The models are evolved until either (1) the outer edge of the dust disk falls within 0.1 AU, in which case all dust is assumed to have accreted onto the star, or (2) the total elapsed time is 107 yr. In the second case, it usually occurs that the surface density distribution of solids, $\Sigma_{\rm s}(r)$, converges to a stationary configuration well before 107 yr.


  \begin{figure}
\par\includegraphics[width=16cm,clip]{aa0036f2.eps}
\end{figure} Figure 2: The decimal logarithm of surface density enhancement  $\xi _{\rm s}$, defined by Eq. (1) for models with  $\alpha = 0.001$ ( left) and  $\alpha = 0.1$ ( right). The values of  $\log \xi_{\rm s}$ are plotted on the plane of m0 and j0 (the initial values of mass and angular momentum of the disk). Within the light-gray region at the lower right corner of each frame the solid component has completely accreted onto the star. The step-like shape of the boundary of that region results from insufficient resolution of the grid (the locations of our models are indicated by dots).
Open with DEXTER

An example of such convergence is illustrated in Fig. 1, where, for a disk with  m0 = 0.02, j0 = 1 and  $\alpha = 0.1$, the surface density of solids is shown at several selected times as a function of the distance r from the star. Initially the outer edge of the disk  $R_{\rm out}$ and the sublimation limit  $R_{\rm evap}$ are located, respectively, at $\sim$50 AU and $\sim$4 AU. The surface density of solids drops from about 0.5 g cm-2 at  $R_{\rm evap}$ to about 0.1 g cm-2 at  $R_{\rm out}$. Interior to  $R_{\rm evap}$ the solids remain in the form of vapour, and as such they are tightly coupled to the gas; however at  $r > R_{\rm evap}$ the solid particles begin to drift with respect to the gas. While drifting, they grow at a rate dependent on r (more slowly at larger r because of lower densities and collision rates). By t = 500 yr enough solid material from the outer disk has drifted across the sublimation limit to cause a noticeable depletion of solids between r = 4 AU and r = 10 AU. At t = 1000 yr the particles arriving from the outer disk to the sublimation limit are large enough to fall onto the descending branch of the drift velocity curve (shown, for example, in Fig. 1 of Weidenschilling 1997), and as a result a maximum of  $\Sigma_{\rm s}$ develops immediately exterior to  $R_{\rm evap}$. Due to large viscosity, the gaseous component of the disk evolves so quickly and cools so rapidly that at $\sim$ $5\times10^3$ yr  $R_{\rm evap}$ begins to decrease, and the inner edge of the dust disk moves closer to the star. At about the same time the solid particles that have arrived at the sublimation limit are so large that they practically stop drifting. As their motion is essentially Keplerian, one may say that they form a "Keplerian barrier'' blocking the flow of solids across  $R_{\rm evap}$. Beyond $\sim$15 AU $\Sigma_{\rm s}$ is now markedly lower than in the initial model, as the solids originally residing there have moved closer to  $R_{\rm evap}$. An outline of the final dust disk is already visible in a region between r = 3.5 AU and r = 15 AU. The evolution is completed at $t = {\sim}3\times10^4$ yr, when all solids that have not evaporated are collected in that region. $\Sigma_{\rm s}$ is now almost uniform, reaching a value $\sim$4 times larger than at the same locations in the initial model. The typical particle radius a is 300 m. At $t = {\sim}9\times10^4$ yr a=1 km - i.e. the solids have evolved into planetesimals. At still later times the planetesimals would tend to grow even further, but the model is no longer valid because it does not include gravitational effects. The snowline (i.e. the inner edge of the planetesimals swarm) stabilizes at  $R_{\rm sl} = 3.5$ AU, while the sublimation limit continues to shrink, and at  $t\sim3\times10^6$ yr it is equal to $\sim$0.5 AU.

Note that the movement of the sublimation limit is always slower than the velocity of the gas, so the vapor does not leak to the region beyond  $R_{\rm evap}$. In addition, unlike in the model proposed by Stevenson & Lunine (1988), in our model there is no significant outward diffusion of water vapor through the sublimation/condensation radius. This is because the vapor is embedded in a dynamically dominant gas and diffuses together with it. In viscous accretion disk the net result of diffusion of adjunct rings of gas is to produce an inward mass flow at times and radii relevant to the sublimation radius. Thus, there is a negligible leak of vapor to the region  $r > R_{\rm evap}$. Overall, in our model, the disk develops a zone where  $T<T_{\rm evap}$, but ice does not exist, as it has been stopped at larger radii by particle growth and an associated decay of particle radial drift.

In the above example the final planetesimals swarm forms a well-defined ring $\sim$ $3.5~{\rm AU} \le r \le {\sim}15$ AU, with  $\Sigma_{\rm s}$ in excess of 2 g cm-2. Whether such a swarm leads to the formation of giant planets cores is unclear as  $\Sigma_{\rm s}$ is lower than required (see Lissauer 1987), but the surface density is relatively constant throughout the entire zone in contrast to the steep decline ( $\Sigma_{\rm s} \sim r^{-3/2}$) in the "minimum mass'' model. Similar structures with different sizes and locations are observed in all models in which the solids are present at the end of the evolution. In all models the ratio of final-to-initial  $\Sigma_{\rm s}$ is greater than 1. The latter property of the final planetesimal swarms is illustrated in Fig. 2. To obtain Fig. 2, the maximum value of the enhancement ratio

 \begin{displaymath}\xi_{\rm s} = \frac{\Sigma_{\rm s}^{\rm max}(t=10^7)}{\Sigma_{\rm s}(R^{\rm max},t=0)}
\end{displaymath} (1)

was calculated for each model, and the distribution of  $\xi _{\rm s}$ was plotted on the (m0, j0) plane. In the above formula  $\Sigma_{\rm s}^{\rm max}$ is the maximum density found in a given disk, and  $R^{\rm max}$ is the location at which that maximum was found. The distributions are shown for  $\alpha=10^{- 3}$ and  $\alpha = 10^{-1}$(the plot for  $\alpha = 10^{-2}$ is omitted because it does not qualitatively differ from the two cases displayed in Fig. 2). In each frame the gray area indicates the region occupied by disks in which all solids are accreted onto the star. The step-like shape of its boundary results from the poor resolution of the grid (the locations of our models are indicated by dots). The enhancement ratio  $\xi _{\rm s}$ increases as j0decreases and m0 increases, i.e. as the initial models become more compact and dense. This is because in denser disks the solid particles grow to "Keplerian'' sizes more quickly, and once the Keplerian barrier is formed at  $R_{\rm evap}$ the solids that still remain at  $r > R_{\rm evap}$are saved from sublimation. At later times the distribution of solids can only evolve toward a more compact configuration, i.e. $\xi _{\rm s}$ can only grow. The earlier the Keplerian barrier forms, the more solids are locked in the outer disk, the further proceeds the "compaction'', and the larger may be the final value of  $\xi _{\rm s}$.

Note that for every m0 a critical value  $j_{\rm0,c}$ can be chosen such that disks with  $j_0<j_{\rm0,c}$ are entirely void of solids. The transit from disks with maximum  $\xi _{\rm s}$ to disks in which all solids accrete onto the star is very abrupt. This is because the Keplerian barrier is not formed at all when the sublimation limit falls initially so close to  $R_{\rm out}$ that there is no time for the particles to grow to Keplerian sizes before they arrive at  $R_{\rm evap}$.

One can also see that in more viscous disks  $\xi _{\rm s}$ is larger. This is because the particles in hotter disks gain larger drift velocities. From Eq. (6) in Kornet et al. (2001) we can estimate the maximum inward drift particle velocity

               $\displaystyle V_{\rm dr}^{\rm max}$ $\textstyle \approx$ $\displaystyle \frac{1}{\rho}\left(\frac{\partial P}{\partial
r}\right)\frac{r}{V_{\rm k}}$  
  $\textstyle \approx$ $\displaystyle \frac{1}{\rho} \left(\frac{{\rm d} P}{{\rm d}
\rho}\right) \left(\frac{\partial \rho}{\partial
r}\right) \frac{r}{V_{\rm k}}$  
  $\textstyle \approx$ $\displaystyle \left(\frac{\partial \ln \rho}{\partial \ln r}\right) \frac{C_{\rm s}^2}{V_{\rm k}}$  

where $\rho$ is gas density, P - gas pressure, and $C_{\rm s}$ - sound velocity. As a result, the outer radius of the planetesimal swarm in hotter disks is smaller, and  $\xi _{\rm s}$ may reach larger values through more "compaction''.

Our results indicate that the most important event in the evolution of the disk occurs when the first "Keplerian'' particles appear at  $R_{\rm evap}$. Let  $R_{\rm evap}^{\rm K}$ be the value of  $R_{\rm evap}$ at the moment the particles arriving at the sublimation limit are so large that their drift time

 \begin{displaymath}t_{\rm dr} = \frac{R_{\rm evap}}{v_{\rm d}\left(a(R_{\rm evap})\right)}
\end{displaymath} (2)

is equal to 106 yr. In the above formula  $v_{\rm d}(a)$ is the drift velocity of particles with a radius a, and a(r) is the radius of particles at a distance r from the star (recall that in our model at each r the particles are all assumed to have the same size). Figure 3 shows how well $R_{\rm evap}^{\rm K}$ correlates with the inner radius of the final planetesimal swarm, i.e. with the final location of the snowline, $R_{\rm sl}^{\rm f}$. To explain the correlation note that as the disk cools, $R_{\rm evap}$ decreases at the rate  $\dot
R_{\rm evap}$.

As long as particles at $R_{\rm evap}$ have velocities  $v_{\rm d}(a(R_{\rm evap}))>\dot R_{\rm evap}$ they move inward across  $R_{\rm evap}$and are destroyed maintaining  $R_{\rm sl}$ at  $R_{\rm evap}$. However, at later times, particle at  $R_{\rm evap}$ becomes larger, decouple from the gas and settle into Keplerian orbits. This leads to $v_{\rm d}(a(R_{\rm evap}=R^{\rm K}_{\rm evap}))\ll \dot R_{\rm evap}$, and while  $R_{\rm evap}$continues to decrease with time, the solids lack inward velocity to follow it. The finite snowline  $R_{\rm sl}^{\rm f}$ forms approximately at the location where the inner movement of  $R_{\rm sl}$ stalls, that is at  $R_{\rm evap}^{\rm K}$. In disk's outer zone ( $r>R_{\rm sl}^{\rm f}$) the inward motion of solids continues for some additional time leading to further "compaction'' of solids.


  \begin{figure}
\par\includegraphics[width=8cm,clip]{0036fig3.ps}
\end{figure} Figure 3: Correlation between the final location of the snowline  $R_{\rm sl}^{\rm f}$ and the distance  $R_{\rm evap}^{\rm K}$ at which "Keplerian particles'' appear at the sublimation limit. Points with error bars indicate the location of our disk models on the ( $R_{\rm sl}$, $R_{\rm evap}^{\rm K}$) plane. The length of the error bars is equal to the local grid spacing in r In the lower frame  $R_{\rm sl}^{\rm f}$ is smaller because disks with smaller $\alpha $ are cooler.
Open with DEXTER

Because many exoplanets have been found relatively close to their parent stars (see Sect. 1), it is natural to identify models leading to the formation of giant planets as close to the star as possible. Note, however, that these are not necessarly the models with the smallest value of  $R_{\rm sl}^{\rm f}$. For given values of $\alpha $ and m0, the models with maximum value of j0 are extended and cool and would evolve to have smallest possible values of  $R_{\rm sl}^{\rm f}$. However, due to their low densities, these are not the disks that are expected to produce planetary cores, thus their snowlines, although located close to the star, are not of great interest. To identify disk models forming giant planets at the closest possible distances to the star would require integration of our models with the models of giant planet formation, a task that is beyond the scope of the present paper. Kornet et al. (2002) have probed this issue working with models similar to those investigated here, but with high-temperature silicates instead of water-ice being the sole species of solids. They have found that, for given values of $\alpha $ and  $m_0>m_0^{\rm crit}$, the disks with minimum (but still swarm-producing) values of j0 are the best candidates to form a close-to-the-star giant planet at  $r \simeq 2$ AU. With water-ice as the sole species of solids the results are qualitatively the same, except the inner edge for giant planet formations is at  $r \simeq 4.4$ AU and $\simeq$8.3 AU for a disk characterized by  $\alpha=10^{- 3}$ and 10-1 respectively.

4 Discussion

In this paper we have concentrated on the issue of the snowline in a protoplanetary disk in which solids, including water-ice, transform to ever larger particles through hierarchical coagulation. In such a disk the solids, initially in the form of dust grains, grow all the way to planetesimals acquiring a significant inward radial velocities in the process. However, these velocities are lost when the size of solid particles approaches planetesimal size, preventing their penetration into the inner disk regardless of its temperature. This evolution of solids occurs on a time scale that is short in comparison with the disk's lifetime. As a result solids decouple from the gas and the location of the snowline cannot be calculated from properties of the gaseous disk, such as the gas temperature corresponding to ice evaporation. Instead, the snowline must be calculated from its direct definition as the minimum radius at which ice exists.

Using an example of a disk evolving from a set of particular initial conditions we have demonstrated that although the snowline ( $R_{\rm sl}$) coincides with the location of an ice sublimation/condensation temperature ( $R_{\rm evap}$) during early stages of disk evolution, at later times  $R_{\rm sl}$ remains fixed at  $R_{\rm sl}^{\rm f}$ whereas  $R_{\rm evap}$decreases as the gaseous disk cools. Thus, at these later times, there exists a zone, $R_{\rm evap}<r<R_{\rm sl}^{\rm f}$, free of ice despite having gas temperatures  $T<T_{\rm evap}$.

This counter-intuitive result can be understood as follows. In an accretion disk with the overall inward mass flux, a radial zone is constantly replenished by material from beyond its outer radius. At the beginning of the disk's evolution the "no-ice'' zone described in the previous paragraph was too hot to support ice. When the temperature in this zone dropped below the condensation level, the zone no longer consisted of its original vapor-rich material, instead it consists of material carried from an outer disk that has been depleted of ice by decoupling of solids from the gas.

Because we are considering an evolving disk, the evaporation radius is a function of time, $R_{\rm evap}(t)$, and calculating a unique snowline using the traditional method of solving equation  $T(r,t)=T_{\rm evap}$ is not possible. The so defined "snowline'' would be a function of time and could be located at arbitrarily small radii. However, $R_{\rm sl}$ calculated directly from ice presence converges to a single value  $R_{\rm sl}^{\rm f}$ and constitutes a reasonable estimate for an inner edge of giant planet formation zone. The location of  $R_{\rm sl}^{\rm f}$ depends on the disk's initial conditions and the value of $\alpha $. We have demonstrated that the value of  $R_{\rm sl}^{\rm f}$ is set by the value of  $R_{\rm evap}^{\rm K}$, an evaporation radius at the moment when solids there are so large that they settle into Keplerian orbits. Thus, the snowline is determined by a complicated interplay between coagulation, sedimentation and gas properties.

We have shown that in the potential giant planet formation zone, $r>R_{\rm sl}^{\rm f}$, there is a significant enhancement of density of solids above its initial value. This enrichment is due to drift of solids from an outer disk. All solids are "compacted'' into a ring between  $R_{\rm sl}^{\rm f}$ and some outer radius that, however, is significantly smaller than the outer radius of the gaseous disk. Thus, in our model, the relatively high surface density of ice can be achieved naturally, as a result of solids decoupling from the gas.

The above conclusions are based on a simplified description of the processes governing the evolution of protoplanetary disks, and one may wonder how robust they are. Below we critically assess several assumptions on which our model is based.

The most radical is the assumption about the size distribution of solid particles, which at each distance from the star are assumed all to have the same diameter. This assumption is responsible for formation of a "totally impenetrable'' barrier at  $R_{\rm sl}^{\rm f}$. However, detailed work by Morfill (1985) and Weidenschilling (1997) showed that the size distribution of solids in a protoplanetry disk quickly converges to a stage in which most of the mass is concentrated in a narrow range of sizes approching the maximum size. Thus, our approximation can be regarded as a reasonable idealization. The emergence of an impenetrable barrier in our models is corroborated by the behavior of solids in the two-dimensional disk model by Weidenschilling (2003) where the incoming mass tends to pile up at a distance where large bodies form, producing a sharp transition in both surface density and mean size. We expect that, in reality, the barrier at  $R_{\rm sl}^{\rm f}$ would leak some ice particles inward, but they would have a negligible mass.

Our second assumption concerns the 100% efficiency of coagulation. In reality, one hardly expects that collisions between solid ice particles always lead to sticking without fragmentation. Interestingly, our calculations suggest that if solids in the protoplanetary disk indeed transform themselves into larger sizes by hierarchical coagulation, and the efficiency of coagulation is low, the disk would be depleted of most of its solids, diminishing the opportunity for planet formation.

To illustrate how the evolution of solids depends on the efficiency of coagulation we have considered again the particular model described in Sect. 3, but with an additional parameter  $0\le\epsilon\le1$ regulating the efficiency of coagulation. To vary this efficiency, we multiply by $\epsilon $ the source function for particle growth by coagulation, f, (see Kornet et al. 2001). The model shown in Fig. 1 was obtained for  $\epsilon=1$. Decreasing $\epsilon $ to 0.5 and 0.25 caused the final outer radius of the icy planetesimals swarm to shrink from 15 AU to 7.5 AU and 2.75 AU, respectively. The final inner radius moved from 2.65 AU to 2.50 AU and 2.33 AU, respectively. For  $\epsilon=0.2$ all solids were accreted onto the star. This demonstrates how low coagulation efficiency inhibits development of giant planets formation zone, at least within a paradigm of hierarchical coagulation.


  \begin{figure}
\par\includegraphics[width=8cm,clip]{0036fig4.ps}
\end{figure} Figure 4: The evolution of solids in an $\alpha = 0.001$ disk with m0 = 0.2and j0 = 9. Top: surface density of solids at the times indicated in the frame. Bottom: final distributions of solids obtained for the indicated vaules of the coagulation efficiency parameter $\epsilon $.
Open with DEXTER

The inferred location of icy planetesimals zone in the early Solar System is between 3 AU and $\sim$50 AU. The estimation of an inner edge (i.e. the snowline) is based on dynamical considerations; a closer snowline could produce Jupiter at close location, which would prevent formation of terrestrial planets. The estimate of the outer radius is based on an observation that the Kuiper Belt seems to be truncated at $\sim$50 AU (Weidenschilling 2003). Of all of our models, the model with $\alpha=10^{- 3}$, $m_0 = 0.2~M_\odot$, j0 = 9 and  $\epsilon=1$ produces the icy planetesimal swarm closest to the one inferred for the Solar System. The evolution of solids surface density in such a model is shown in Fig. 4, the final swarm extends from 3 AU to $\sim$25 AU with surface density ranging from $\sim$8 to $\sim$12 g cm-2. Note that the solid surface density at 5 AU is above 10 g cm-2, which was the standard value used by Pollack et al. (1996) to form Jupiter in less than 10 Myr.

Our model provides a reasonable scenario for the formation of icy planetesimal swarms in the Solar System. Note, however, that a good fit requires  $\epsilon\approx 1$. For $\epsilon=0.5$ and  $\epsilon=0.25$ the outer radius of the swarm decreases to 15 AU and 10 AU, respectively, not far enough to account for the Solar System. Indeed, to account for the entire Kuiper belt, models with $\epsilon $increasing with the distance from the star, up to values greater than unity in the outer disk are needed. Note that given the numerous approximations on which the form of f is based, $\epsilon>1$ is not as incongruous as it appears. For example, in low-viscosity disks ( $\alpha<{\sim}0.01$) the relative particle velocities induced by turbulence become smaller than the differential drift velocities (Weidenschilling 1997). The latter are not included in f, and, as a result, in such a regime our approach underestimates the particle growth rate which could be countered by increasing the value of $\epsilon $. Also, in our models the coagulation rate is underestimated by the absence of other species of solids. Their presence would increase the solid density and thus coagulation rate. This also can be qualitatively countered by increasing the value of $\epsilon $.

These factors notwithstanding, a high coagulation efficiency is needed in our model to produce an icy planetesimal swarm leading to a Solar-like planetary system. However, the high efficiency of coagulation is questionable. Indeed, a so-called "meter-sized particle barrier'', wherein particles growing to meter-size rubble achieve high relative velocities and their collisions lead to fragmentation, is discussed in the literature (see, for example Weidenschilling & Cuzzi 1993). Thus, the demand for a high coagulation efficiency in our model, to produce planets, can be view as an argument against a hierarchical coagulation as the sole agent of solid aggregation. Perhaps some form of a collective process to transform solids from dust directly into planetesimals, thus bypassing the interim stages and avoiding large drift velocities that lead to the loss of solids, has to take place for the disk to create planets. Such a mechanism, in the form of gravitational instability of the dust layer, was originally proposed by Goldreich & Ward (1973), but rejected when it became clear that the presence of turbulence would prevent dust from settling into a thin enough layer to become gravitationally unstable. Later, a variant of the gravitational instability mechanism was proposed wherein the dust density was enhanced in gaseous vortices to produce clumps dense enough to be collapsed into planetesimals by gravitational instability (Tanga et al. 1996). Barring the existence of some collective process, the paradigm of hierarchical coagulation requires a high efficiency of coagulation for the formation of giant planets. This requirement is independent of the details of our model and is rooted in the high drift velocities of intermediate-size particles.

Finally, we stress that the models discussed in the present communication are not meant to be quantitative. Instead, they should be viewed as an illustration of the basic processes involved, which, we believe, is qualitatively correct despite its simplicity. Future improvements of the model should allow for a range of particle sizes at each position in the disk and for interactions between different types of particles; a more detailed description of the gas component is also a desirable option.

Acknowledgements
This work was supported in part by the Polish Committee for Scientific Research through the grant No. 2P03D01419. K.K. and M.R. also benefited from the European Commission RTN grant No. HPRN-CT-2002-00308. This is Lunar and Planetary Institute Contribution No. 1188.

References



Copyright ESO 2004