The Astrophysical Journal, 551:L167L170, 2001 April 20
2001. The American Astronomical Society. All rights reserved. Printed in U.S.A.
FORMATION OF PLANETARY-MASS OBJECTS BY PROTOSTELLAR COLLAPSE AND FRAGMENTATION Alan P. Boss
Department of Terrestrial Magnetism, Carnegie Institution of Washington, 5241 Broad Branch Road, NW, Washington, DC 20015-1305 Received 2001 February 9; accepted 2001 March 8; published 2001 April 17
ABSTRACT Searches for very low mass objects in young star clusters have uncovered evidence for free-oating objects with inferred masses possibly as low as 515 Jupiter masses (MJup), similar to the masses of several extrasolar planets. We show here that the process that forms single and multiple protostars, namely, collapse and fragmentation of molecular clouds, might be able to produce self-gravitating objects with initial masses less than 1 MJup, provided that magnetic eld tension effects are important and can be represented approximately by diluting the gravitational eld. If such fragments can be ejected from an unstable quadruple protostar system, prior to gaining signicantly more mass, protostellar collapse might then be able to explain the formation of free-oating objects with masses below 13 MJup. These objects might then be best termed subbrown dwarf stars. Subject headings: binaries: general hydrodynamics ISM: clouds ISM: kinematics and dynamics MHD stars: formation

1. INTRODUCTION

Searches for brown dwarf stars in the eld (Delfosse et al. 1999; Kirkpatrick et al. 2000; Leggett et al. 2000) and young star clusters (Bejar, Zapatero Osorio, & Rebolo 1999; Ardila, Martn, & Basri 2000; Luhman et al. 2000; Najita, Teide, & Carr 2000) have discovered hundreds of candidates with estimated masses below the hydrogen-burning limit of 75 MJup. Freeoating objects have even been found (Lucas & Roche 2000; Zapatero Osorio et al. 2000) with inferred masses possibly below the deuterium-burning limit of 13 MJup. On the basis of the suggestion (Saumon et al. 1996) that objects incapable of burning deuterium be dened as planets, free-oating objects with inferred masses in the range from 5 to 13 MJup have been termed planets (Lucas & Roche 2000; Zapatero Osorio et al. 2000). Radial velocity surveys (Marcy, Butler, & Vogt 2000; Udry et al. 2000) have detected over 50 likely planetary companions to solar-type stars, with minimum masses in the range from 0.25 to 15 MJup. Evidently the least massive, isolated objects found in young star clusters could be less massive than the most massive planetary companions, blurring the distinction between stars and planets. These observations thus raise an important theoretical question: Can very low mass free-oating objects be formed directly in star-forming regions, or must they form in planetary systems and later be ejected through dynamical interactions? Until recently there was believed to be a distinct gap in mass between the most massive planet (Jupiter) and the least massive brown dwarf star, a mass gap that would simplify the categorization of these objects. This expected gap was based in part on theoretical estimates of the minimum mass of a star that might be produced during the gravitational collapse of dense molecular cloud cores to form protostellar systems. Analytical Jeans mass estimates (Low & Lynden-Bell 1976; Rees 1976; Silk 1977) of opacity-limited fragmentation suggested that the minimum stellar mass might be in the range from 3 to 7 MJup. Numerical hydrodynamical calculations (Boss 1988, 1993) of protostellar collapse and fragmentation produced a minimum protostellar fragment mass from 3 to 10 MJup, in very good agreement with the analytical estimates. However, both the analytical and numerical estimates ignored the effects of magnetic elds, which might be expected to discourage fragmentation and so to increase the minimum protostellar mass. L167
These estimates also only referred to the initial masses of selfgravitating clumps, which will continue to accrete mass from the infalling cloud envelope. Hence, these estimates seemed to be quite rm lower bounds on stellar masses, with the true minimum protostellar mass presumably being considerably higher. Theoretical estimates (Lin & Papaloizou 1980) of the maximum mass of a gas giant planet suggested that a gap would open once a planet reached about a Jupiter mass, preventing further growth of the planet. Planets with masses between 1 and 10 MJup thus were not thought to be common, neatly separating planets from brown dwarfs. However, observations of planetary-mass companions to solar-type stars (Marcy et al. 2000; Udry et al. 2000) and free-oating objects in young star clusters (Lucas & Roche 2000; Zapatero Osorio et al. 2000) do not show any evidence for a universal mass gap in this rangethere seem to be objects in both settings with masses in the range from 1 to 10 MJup. Recent work on gas giant planet formation by core accretion suggests that gap formation need not completely stop further growth (Bryden et al. 1999; Kley 2000), while the alternative mechanism of forming gas giant protoplanets rapidly through a disk instability (Boss 2000a) may be able to form protoplanets with masses in the range from 1 to 10 MJup. However, the question remains: Can the star formation process produce gravitationally bound objects as small in mass as 5 MJup? We present here the rst numerical hydrodynamical calculation showing the formation of planetary-mass objects during a well-resolved calculation of the self-gravitational collapse of a dense molecular cloud core. A dynamically unstable, quadruple protostar system forms, with initial fragment masses well below a Jupiter mass, low enough to yield a considerable margin between the initial fragment masses and observations of free-oating planetary-mass objects (Lucas & Roche 2000; Zapatero Osorio et al. 2000). The new result depends on the approximate inclusion of magnetic eld effects, leading to a central rebound and cooling that allows smaller mass objects to form by fragmentation and still be gravitationally bound than in the previous nonmagnetic calculations (Boss 1988, 1993).

FORMATION OF PLANETARY-MASS OBJECTS

2. NUMERICAL METHODS

Vol. 551
The model was calculated with a three-dimensional numerical code (Boss & Myhill 1992) that computes nite-difference solutions of the equations of hydrodynamics, gravitation, and radiation transport for a gas and dust cloud. The code is secondorder accurate in both space and time. Radiative transfer is handled in the Eddington and diffusion approximations, including detailed equations of state and dust grain opacities. Poissons equation for the clouds gravitational potential is solved by a spherical harmonic expansion, including terms up to l, m p 16. The computational grid consists of a spherical coordinate grid with Nr p 200, Nv p 22 for p/2 v 0 (symmetry through the midplane is assumed for p v 1 p/2), and Nf p 64 for p f 0 (symmetry through the rotational axis is assumed for 2p f 1 p). The radial grid contracts to follow the collapsing inner regions of the cloud and to provide sufcient spatial resolution to ensure satisfaction of the four Jeans conditions for a spherical coordinate grid (Truelove et al. 1997; Boss et al. 2000). The f-grid is uniformly spaced, while the v-grid is compressed toward the midplane, where the minimum grid spacing is 03. Boss (2000b) describes calculations with substantially higher azimuthal resolution than the present calculation, which serve to demonstrate the robustness of the results.

3. INITIAL CONDITIONS

The model begins from initial conditions chosen to represent as closely as possible observed magnetically supported, dense cloud cores (Myers et al. 1991; Ward-Thompson, Motte, & Andre 1999; Crutcher 1999). The clouds initial radial density prole is Gaussian with a central density (r 0 p 2 # 1018 g cm3) 20 times higher than that at the boundary. The cloud is initially oblate with an axis ratio of 2 : 1 and with a random pattern of noise added to the density at the 10% level. The cloud radius is R p 1.0 # cm and the initial temperature is 10 K, yielding a mass of 2.1 M, and an initial ratio of thermal to gravitational energy a i p 0.30. An initial solid-body rotation rate of 1014 rad s1 (Goodman et al. 1993) around the short axis of the cloud leads to an initial ratio of rotational to gravitational energy bi p 1.1 # 104. A number of other numerical models have been calculated with different initial rotation rates, with different initial central to boundary density ratios, and with initially prolate rather than oblate shapes (A. P. Boss 2001, in preparation), yielding results quite similar to those presented below except that the initially prolate clouds tend to produce binary rather than multiple protostar systems.

4. MAGNETIC FIELD APPROXIMATIONS
Mouschovias 1995) can be crudely approximated (Boss 1997, 1999, 2000b) by specifying that B 0 p B 0 (t) p B 0i (1 t/tAD ), where B 0i is a constant equal to 200 mG, t is the model time, and the ambipolar diffusion timescale tAD is 10t f f. The free-fall time t f f p t f f (r 0 ) p 1.486 # s. With B 0i p 200 mG at n 0 p 6 # cm3, the model starts with a ratio of magnetic to gravitational energy of gi p 0.43. With the addition of the magnetic pressure, the cloud is initially in a quasi-equilibrium state such that the cloud rotates for several free-fall times without collapsing or expanding signicantly (Boss 1997). Loss of magnetic eld support by ambipolar diffusion leads to a steady increase in the maximum density of the cloud and eventually to a full-edged dynamic collapse phase. During the collapse phase, initially straight magnetic eld lines bend inward and exert a tension force that counteracts gravity. For a thin disk with a constant mass-to-ux ratio, the acceleration due to magnetic tension is proportional to the gravitational acceleration, greatly simplifying its inclusion in theoretical models (Shu & Li 1997; Basu 1997; Nakamura & Hanawa 1997). Once the central density exceeds 1015 g cm3, denoting the onset of the dynamic collapse phase, magnetic tension forces are approximated in the present models by diluting the overall gravitational potential F by a factor involving a function of time F r F {1 1.1 [(tAD t) /tAD ]2}, where 1.1 is a coefcient that depends on B 0, r0, and the initial temperature. This approximation strictly applies only once a thin disk is formed. However, magnetic tension has the same qualitative effect of opposing gravity during the dynamic collapse phase immediately prior to disk formation, and so this approximation is used throughout the cloud once the maximum density exceeds 1015 g cm3. Fragmentation does not begin in these models until well after a central rebound occurs and the maximum density has increased by a factor of 104 from that at the initiation of gravity dilution, so the choice of the exact onset of gravity dilution should not be critical to the outcome of the model. The validity of this approximate treatment of magnetic eld effects needs to be veried by subsequent calculations with a three-dimensional magnetohydrodynamical code that also includes self-gravity and radiative transfer.

5. RESULTS

Magnetic pressure effects are represented by adding the magnetic eld pressure B 2/8p to the gas pressure, which is an exact approximation for straight magnetic eld lines and high conductivity, two very good approximations for dense cloud cores. Observationally (Crutcher 1999) and theoretically (Mouschovias 1991) the magnetic eld strength B depends on the gas density in dense cloud cores as B p B 0 (r/r 0 ) k, where k and B 0 is the eld strength at a reference density r 0. Recent work (Desch & Mouschovias 2001) has shown that ambipolar diffusion, rather than Ohmic dissipation, is responsible for magnetic decoupling, when the magnetic eld ceases to have a signicant effect on the predominantly neutral species in the collapsing cloud. The effects of ambipolar diffusion (Ciolek &

The cloud is initially supported against collapse primarily by magnetic pressure but begins to undergo dynamic collapse after 2.45tff as a result of ambipolar diffusion. The cloud then becomes optically thick and undergoes compressional heating at the center. The combination of magnetic tension forces (i.e., gravity dilution) and thermal and magnetic pressure (centrifugal forces are minimal in this slowing rotating cloud) leads to a reversal of the collapse at the center and to an outward rebound at 2.61tff. The temperature begins to rise strongly in the region between the central rebound and the infalling outer cloud, and eventually an asymmetric, off-axis density maximum forms and begins to fragment into rst two and then four clumps. Figure 1 shows that by a time of 2.615tff the cloud has fragmented into a quadruple protostar system. In the absence of magnetic tension forces, a central rebound and fragmentation do not occur for this cloud (Boss 2000b). The four fragments shown in Figure 1 have ratios of thermal to gravitational energy of 0.3 and so are gravitationally bound. The fragments have masses of either 0.3 or 0.5 MJup, where a fragment is dened to include all adjacent regions with a density greater than 0.1 that of the maximum density in the fragment.

No. 2, 2001

Fig. 1.(a) Midplane density contours, (b) temperature contours, and (c) velocity vectors at the center of an initially oblate, magnetic cloud core that has collapsed, rebounded, and fragmented into a quadruple protostar system with fragment masses less than a Jupiter mass. All three plots are shown at a time t p 2.615tff and depict a region with a radius of 4.7 # 1013 cm (3.1 AU), containing 65 radial grid cells and effectively 128 azimuthal grid cells with the assumed symmetry. The density maximum in (a) is 6.3 # 1010 g cm3, and the temperature maximum in (b) is 100 K. In both (a) and (b), contours represent changes by factors of 1.3. Density and temperature minima occur in the central regions of the cloud, where the velocity eld is outwardly directed. Density maxima occur at the centers of the four clumps, while the maximum temperature occurs in an annulus outside the four clumps. The maximum velocity in (c) is 9.0 # 104 cm s1. Only every fourth cell in radius and every fourth cell in azimuth are shown in (c) for clarity.
These fragment masses are about a factor of 10 or more times smaller than the lowest mass fragments obtained in nonmagnetic collapses (Boss 1988, 1993). The average temperature of the fragments is 31 K, considerably less than would be expected for a fragment with an average density of 1010 g cm3. On the basis of nonmagnetic collapse calculations, a temperature of 200 K would be expected (Boss 1988) instead. Lower temperatures occur in the new model because the magnetically mediated rebound leads to decompressional cooling in the central regions and to lower fragment temperatures. This in turn permits

the formation of lower mass fragments. The Jeans mass depends on temperature to the 3/2 power, so a decrease in temperature by a factor of 200/31 p 6.5 results in a reduction in the Jeans mass by a factor of 16, consistent with the lower fragment masses. The fragment masses are in fact comparable to the Jeans mass. At the time shown in Figure 1, the calculation obeys the four Jeans conditions for a spherical coordinate grid (Truelove et al. 1997; Boss et al. 2000). However, when the calculation is continued beyond the point depicted in Figure 1, so that the
Jeans conditions are no longer satised and there is a danger of spurious numerical fragmentation occurring, the fragments shown in Figure 1 continue to contract and evolve, while no further fragmentation occurs. The fragments become increasingly distinct and contract to maximum densities a factor of 1000 times higher than shown in Figure 1, implying their continued survival, even when faced with the effects of magnetic pressure and gravity dilution representing magnetic tension.
6. DISCUSSION AND CONCLUSION
While the subsequent evolution of the quadruple system shown in Figure 1 is highly uncertain, and may involve mergers and accretion of substantial gas, if three or four fragments survive long enough, it is possible that one or more will be ejected during subsequent close mutual encounters because of the nonhierarchical nature of the system. In order to escape the cloud, a fragment must gain enough kinetic energy (per unit mass of the fragment) during a close three-body encounter to exit the gravitational well of the entire cloud. We can approximate the former by the depth of the gravitational well close to one of the fragments, GMf /R f , where G is the gravitational constant, Mf is the fragment mass, and R f is the fragment radius. In order for escape to occur, this energy per unit mass must exceed the clouds gravitational binding energy, GMc /R c, where Mc is the cloud mass and R c is the cloud radius. Escape thus requires that the ratio Ef /Ec p (Mf /Mc )(R c /R f ) 1 1. For the fragments shown in Figure 1, Mf 0.5 MJup, Mc 2100 MJup, R f 5 # cm, and R c p cm, so the ratio Ef /Ec 5, sufcient for escape to occur. However, this estimate errs on the side of escape, because the clouds potential well will deepen beyond GMc /R c as collapse proceeds. One can more directly assess the depth of the clouds potential well through investigating the gravitational potential generated by the codes Poisson solver. At the time shown in Figure 1, the minima in the gravitational potential occur at the locations of the fragments, while the po-

tential elsewhere near the center of the cloud is about 5% smaller in absolute value. This small difference is probably insufcient to lead to ejection. However, the fragments have the advantage that their collapse time is much shorter than that of the cloud as a whole, owing to their much higher mean density, and so they should be able to out-contract the overall cloud. Shortly after the time shown in Figure 1, in fact, the maximum fragment density has increased by a factor of about 30, and the gravitational potential minimum at a fragment is about 20% greater than at the clouds center, making escape more probable. Thus, by the time that dynamical interactions between the fragments shown in Figure 1 could lead to a threebody encounter, the fragments may have contracted to smaller radii and to much higher densities, making the escape of a single fragment easier to achieve. Clearly, this conjecture remains to be veried by further work. The quadruple protostar system shown in Figure 1 could disintegrate and eject a protostar that will become a free-oating object with a mass less than 13 MJup. The primary process of protostellar collapse might then explain the formation of planetary-mass, free-oating objects, obviating any need to invoke the secondary process that later produces planetary companions to explain the formation of very low mass free-oating objects. Because these objects would then form through the same process that forms normal main-sequence (dwarf) stars and brown dwarf stars, free-oating objects with masses in the range of 513 MJup might then be termed subbrown dwarf stars, in analogy with subdwarf stars, which are less luminous than dwarf stars. The numerical calculations were performed on the Carnegie Alpha Cluster, which, along with this work, is partially supported by the National Science Foundation under grants AST 99-83530 and MRI 99-76645. I thank Gotthard Saghi-Szabo for his management of the Cluster and Richard Durisen for a perceptive review of the manuscript.
REFERENCES Ardila, D., Martn, E., & Basri, G. 2000, AJ, 120, 479 Basu, S. 1997, ApJ, 485, 240 Bejar, V. J. S., Zapatero Osorio, M. R., & Rebolo, R. 1999, ApJ, 521, 671 Boss, A. P. 1988, ApJ, 331, 370. 1993, in The Realm of Interacting Binary Stars, ed. J. Sahade, G. McCluskey, & Y. Kondo (Dordrecht: Kluwer), 355. 1997, ApJ, 483, 309. 1999, ApJ, 520, 744. 2000a, ApJ, 536, L101. 2000b, ApJ, 545, L61 Boss, A. P., Fisher, R. T., Klein, R. I., & McKee, C. F. 2000, ApJ, 528, 325 Boss, A. P., & Myhill, E. A. 1992, ApJS, 83, 311 Bryden, G., Chen, X., Lin, D. N. C., Nelson, R. P., & Papaloizou, J. C. B. 1999, ApJ, 514, 344 Ciolek, G. E., & Mouschovias, T. Ch. 1995, ApJ, 454, 194 Crutcher, R. M. 1999, ApJ, 520, 706 Delfosse, X., Tinney, C. G., Forveille, T., Epchtein, N., Borsenberger, J., Fouque, P., Kimeswenger, S., & Tiphene, D. 1999, A&AS, 135, 41 ` Desch, S. J., & Mouschovias, T. Ch. 2001, ApJ, 550, 314 Goodman, A. A., Benson, P. J., Fuller, G. A., & Myers, P. C. 1993, ApJ, 406, 528 Kirkpatrick, J. D., et al. 2000, AJ, 120, 447 Kley, W. 2000, MNRAS, 313, L47 Leggett, S. K., et al. 2000, ApJ, 536, L35 Lin, D. N. C., & Papaloizou, J. 1980, MNRAS, 191, 37 Low, C., & Lynden-Bell, D. 1976, MNRAS, 176, 367 Lucas, P. W., & Roche, P. F. 2000, MNRAS, 314, 858 Luhman, K. L., Rieke, G. H., Young, E. T., Cotera, A. S., Chen, H., Rieke, M. J., Schneider, G., & Thompson, R. I. 2000, ApJ, 540, 1016 Marcy, G. W., Butler, R. P., & Vogt, S. S. 2000, ApJ, 536, L43 Mouschovias, T. Ch. 1991, ApJ, 373, 169 Myers, P. C., Fuller, G. A., Goodman, A. A., & Benson, P. J. 1991, ApJ, 376, 561 Najita, J. R., Teide, G. P., & Carr, J. S. 2000, ApJ, 541, 977 Nakamura, F., & Hanawa, T. 1997, ApJ, 480, 701 Rees, M. 1976, MNRAS, 176, 483 Saumon, D., Hubbard, W. B., Burrows, A., Guillot, T., Lunine, J. I., & Chabrier, G. 1996, ApJ, 460, 993 Shu, F. H., & Li, Z.-Y. 1997, ApJ, 475, 251 Silk, J. 1977, ApJ, 214, 152 Truelove, J. K., Klein, R. I., McKee, C. F., Holliman, J. H., II, Howell, L. H., & Greenough, J. A. 1997, ApJ, 489, L179 Udry, S., et al. 2000, A&A, 356, 590 Ward-Thompson, D., Motte, F., & Andre, P. 1999, MNRAS, 305, 143 Zapatero Osorio, M. R., Bejar, V. J. S., Martin, E. L., Rebolo, R., Barrado y Navascues, D., Bailer-Jones, C. A. L., & Mundt, R. 2000, Science, 290, 103

The Astrophysical Journal, 536:L101L104, 2000 June 20
2000. American Astronomical Society. All rights reserved. Printed in U.S.A.
POSSIBLE RAPID GAS GIANT PLANET FORMATION IN THE SOLAR NEBULA AND OTHER PROTOPLANETARY DISKS Alan P. Boss
Department of Terrestrial Magnetism, Carnegie Institution of Washington, 5241 Broad Branch Road, NW, Washington, DC 20015-1305; boss@dtm.ciw.edu Received 2000 March 27; accepted 2000 May 5; published 2000 June 13
ABSTRACT Gas giant planets have been detected in orbit around an increasing number of nearby stars. Two theories have been advanced for the formation of such planets: core accretion and disk instability. Core accretion, the generally accepted mechanism, requires several million years or more to form a gas giant planet in a protoplanetary disk like the solar nebula. Disk instability, on the other hand, can form a gas giant protoplanet in a few hundred years. However, disk instability has previously been thought to be important only in relatively massive disks. New three-dimensional, locally isothermal, hydrodynamical models without velocity damping show that a disk instability can form Jupiter-mass clumps, even in a disk with a mass (0.091 M, within 20 AU) low enough to be in the range inferred for the solar nebula. The clumps form with initially eccentric orbits, and their survival will depend on their ability to contract to higher densities before they can be tidally disrupted at successive periastrons. Because the disk mass in these models is comparable to that apparently required for the core accretion mechanism to operate, the models imply that disk instability could obviate the core accretion mechanism in the solar nebula and elsewhere. Subject headings: accretion, accretion disks hydrodynamics planetary systems solar system: formation

1. INTRODUCTION

Spectroscopic searches have discovered about three dozen very low mass companions to nearby solar-type stars (Marcy, Cochran, & Mayor 2000). These companions are widely assumed to be gas giant planets similar to Jupiter, an identication that has been strengthened by the discovery of the rst transiting extrasolar planet (Charbonneau et al. 2000; Henry et al. 2000) and the resulting determination of a planetary radius and density close to that expected for a hot Jupiter (Burrows et al. 2000). However, most (about 3 ) of these objects are more mas4 sive than Jupiter, with likely masses in the range 110 Jupiter masses (MJup). The rst conrmed system of such objects, around the star u Andromedae, has a combined mass roughly 10 times that of our solar systems planets (Butler et al. 1999). Theories of gas giant planet formation derived to explain the formation of Jupiter may not be adequate to account for these considerably more massive objects. Conventional wisdom suggests that Jupiter and Saturn formed by the core accretion mechanism (Mizuno 1980; Pollack 1984), which is a two-step process. First, collisional accumulation of icy and rocky planetesimals leads to the runaway growth (Lissauer 1987) of a 10 M solid core. Second, as this core grows, it acquires a growing atmosphere of nebular gas, which eventually becomes unstable to collapse, leading to a phase of rapid accretion of hydrogen and helium gas onto the protoplanets envelope. The rst step requires on the order of a half-million years, while the second step requires several million years or more according to the most recent models (Pollack et al. 1996). Core accretion became the favored mechanism for several reasons. The critical core mass required for gas accretion is not strongly dependent on orbital distance in the solar nebula (Mizuno 1980; Pollack 1984), offering an explanation for the then-inferred similarity in core masses for the giant planets (Stevenson 1982). By specifying the amount of nebular gas accreted during the rapid accretion phase, core accretion can lead to arbitrary ratios of core mass to envelope L101

mass and hence to nonsolar bulk compositions. Finally, core accretion begins with the collisional accumulation of solids, the uncontested mechanism of terrestrial planet formation (Wetherill 1990, 1996). However, core accretion is not without problems, some of which have only recently emerged. The perpetual problem has been the timescale, which is comparable to or longer than estimates of the lifetime of planet-forming disks. Such estimates range from a few million years or less in low-mass starforming regions (Wolk & Walter 1996; Jayawardhana et al. 1999) to much less than one million years in high-mass starforming regions (Bally et al. 1998), where most stars are born. If the disk gas has already been dissipated, a massive gaseous envelope cannot be accreted. Other problems are more recent. New models of the interior of Jupiter suggest a considerably smaller core mass or even no core at all (Guillot, Gautier, & Hubbard 1997). If the core masses are not large enough to initiate sustained gas accretion, the core accretion mechanism will fail. Furthermore, a 10 M core is expected to migrate inward toward the protosun in 104 yr through gravitational interactions with the nebula, well before it can accrete a signicant gaseous envelope (Goldreich & Tremaine 1980; Ward 1997a, 1997b; Tanaka & Ida 1999; Miyoshi et al. 1999). If a core manages to avoid this fate by accreting a sufciently massive gaseous envelope, it will open a gap in the disk, slowing the further growth of the planet (Bryden et al. 1999). The only alternative mechanism appears to be disk instability, where a gravitationally unstable disk fragments directly into self-gravitating clumps of gas and dust that can contract and become giant gaseous protoplanets. Disk instability was rst proposed decades ago (Kuiper 1951; Cameron 1978), but was largely discarded because disk instability seemed to be unable to form large solid cores, as solids are expected to dissolve in the hot envelope rather than settle to the center of a gas giant planet (Slattery, DeCampli, & Cameron 1980; Stevenson 1982). However, a solid core could form by coagulation

PLANET FORMATION IN SOLAR NEBULA
TABLE 1 Effects of Increasing Spatial Resolution on Disk Instability Models Model 64d. 64f. 128f. 256f. 256pf. 512pf. Damping Yes No No No No No Nf NYlm 6.3 7.9 1.6 2.0 1.6 4.3 rmax # # # # # # 108 107

Vol. 536

of dust grains and their sedimentation to a protoplanets center in 103 yr (Boss 1997, 1998a), well before the protoplanet contracts to planetary densities and temperatures, which requires 105 yr (Bodenheimer et al. 1980). A 1 MJup protoplanet of solar composition (Z = 0.02) could then quickly form a 6 M solid core, a value that lies in the middle of the presently inferred range for Jupiter (Guillot et al. 1997). The disk instability process itself is quite fast, as it occurs on orbital timescales, so that clump formation and dust grain sedimentation proceed nearly simultaneously, on a timescale of 103 yr. A disk instability is able to form multiple-MJup planets (Boss 1998a), because the more massive (and cooler) a disk is, the more likely it is to undergo the instability. Because the massive clumps form rapidly and only afterward open up disk gaps, there is no problem with growing to the large masses inferred for many extrasolar planets. Similarly, the disk instability mechanism avoids entirely the danger of signicant orbital migration prior to reaching 1 MJup. Once a gap forms around a giant protoplanet formed by either mechanism, the protoplanet will be subject to orbital migration driven by the subsequent evolution of the disk, a likely means for explaining the short orbital periods of the hot Jupiters (Lin & Papaloizou 1986; Lin, Bodenheimer, & Richardson 1996; Trilling et al. 1998). The disk instability mechanism also has problems. Gravitationally unstable disks are thought to evolve as a result of gravitational torques that transport mass inward and angular momentum outward and thereby avoid forming long-lived clumps (Cassen et al. 1981; Papaloizou & Savonije 1991; Laughlin & Bodenheimer 1994). A marginally unstable disk may then require a trigger to induce the instability, such as episodic mass accretion onto the disk (Boss 1997) or a close encounter with another star (Bofn et al. 1998). The instability may require a fairly massive disk (Cassen et al. 1981; Laughlin & Bodenheimer 1994; Boss 1998a), considerably more massive than the minimum mass necessary to make the solar system. As a result, disk instability may not be able to make gas giant planets as low in mass as Saturn or gas giant planets with nonsolar bulk compositions unless their gaseous envelopes are preferentially lost by overow through their Roche lobes during inward orbital migration. Here we demonstrate that disk instability can operate in a disk with a mass comparable to that inferred for the solar nebula and indeed seemingly required for core accretion to succeed, removing one of the major problems with disk instability.

2. HYDRODYNAMICAL MODELS

Previous three-dimensional hydrodynamics models (Boss 1997, 1998a) showed that the instability could proceed in a disk with a mass of 0.14 M, inside a radius of 10 AU; however, the total mass for a disk extending beyond 10 AU could be
signicantly higher. These models had a 1 M, central protostar, so that the ratio of the disk mass to the star mass was at least Md /Ms = 0.14. Semianalytical studies of a specic eccentric disk instability mechanism suggested that instability could occur only in disks with Md /Ms 1 0.19 (Shu et al. 1990). Several numerical studies have found clumping to occur when Md /Ms 1 (Cassen et al. 1981; Laughlin & Bodenheimer 1994). One study of two-dimensional (thin) disks found clumping to occur in smoothed particle hydrodynamics calculations of disks with Md /Ms = 0.05 and 0.1, but was unable to conrm the result with nite differences calculations; however, for Md /Ms 0.2, both numerical methods yielded qualitatively similar behavior (Nelson et al. 1998). Estimates of a lower bound on the mass of the solar nebula fall in the range 0.010.1 M, (Weidenschilling 1977), so that the solar nebula must have had a ratio Md /Ms 0.010.1. The numerical models presented here have a disk mass of 0.091 M, within a radius of 20 AU and a 1 M, central protostar, so that Md /Ms = 0.091. The nite differences code used here solves the threedimensional equations of hydrodynamics and the Poisson equation on a spherical coordinate grid. The code is the same as that previously employed (Boss 1997, 1998a; with the exception noted below) and has been shown to be accurate to second order in space and time (Boss & Myhill 1992). The number of grid points in each direction is Nr = 101, Nv = 23 in p/2 v 0, and Nf = 64, 128, 256, or 512 (Table 1). The radial grid is uniformly spaced between either 1 or 4 AU and 20 AU, with boundary conditions at both 1 (or 4) AU and 20 AU chosen to absorb radial velocity perturbations. The v grid is compressed into the midplane to ensure adequate vertical resolution (Dv = 03 at the midplane). The f grid is uniform. The. central protostar is allowed to wobble in response to the growth of nonaxisymmetry in the disk, thereby preserving the location of the center of mass of the star/disk system. The number of terms in the spherical harmonic expansion for the gravitational potential of the disk is varied from NYlm = 16 to 32 to 48. The disk initially has a surface density prole with j r1/2r1 in the inner disk, steepening to j r3/2 in the outer disk, similar to that thought to be appropriate for the core accretion mechanism (Lissauer 1987). The surface density at 5 AU falls within the likely bounds for the solar nebula (Weidenschilling 1977) and is about 50% higher than in the standard core accretion model (Pollack et al. 1996). The initial disk density distribution is seeded with a mixture of nonaxisymmetric cos (mf) density perturbations (m = 1 , 2, 3, and 4 modes with amplitude am = 0.01) and with random noise at a somewhat lower amplitude. The three-dimensional models start with the thermal structure of the corresponding axisymmetric disk as calculated by a two-dimensional radiative hydrodynamics code (Boss 1996). Because the axisymmetric model extended to only 10 AU, whereas the new models extend to 20 AU, the temperature prole from 10 to 20 AU was taken to be constant at 40 K. This appears to be a reasonable temperature for the outer solar nebula, based on a variety of cosmochemical constraints (Boss 1998b). With this temperature prole, the disk initially has a (Toomre 1964) Qmin 1.3; marginally unstable disks have Qmin-values less than 1.5 (Papaloizou & Savonije 1991; Nelson et al. 1998; Boss 1998a). The disk is assumed to evolve as a locally isothermal disk, meaning that the initial radial temperature prole is held xed (Boss 1997, 1998a). Previous models by Boss (1997, 1998a) damped the disks translational velocities whenever vr or vv ! 0 in order to

No. 2, 2000

maintain a stable inner disk without suffering a severe time step penalty. Marginally unstable disk models calculated by Pickett et al. (1998, 2000a), however, did not include this damping and did not lead to the formation of the long-lived clumps found by Boss (1997, 1998a). The present study includes models with and without velocity damping in order to learn the effects of this artice with the present code (see Pickett et al. 2000b).

3. RESULTS

Table 1 summarizes the six models and the approximate maximum density encountered during the evolution. As the spatial resolution is increased, the maximum density encountered increases, as expected for calculations that are converging toward a solution involving clump formation. Model 64d was calculated in much the same manner as previous models (Boss 1997, 1998a), with limited f spatial resolution and velocity damping. After 3000 yr, the result was similar: two multiple-MJup clumps formed and orbited apparently stably on circular orbits. Model 64f was identical to model 64d except for having the velocity damping removed beyond 5 AU. Model 64f became nonaxisymmetric within 300 yr and still produced multiple-MJup clumps, but the clumps were on eccentric orbits and disappeared after at most a single revolution period. Thin laments persisted throughout the evolution, and new clumps continued to form out of them and then disappear. Model 128f, with doubled f resolution, behaved much like model 64f, but the clumps reached maximum densities about a factor of 2 higher. Evidently, velocity damping had the unintended effects of preserving the clumps indenitely and of slowing the evolution, as has also been found by Pickett et al. (2000b). For the subsequent models, velocity damping was removed completely and the inner boundary was moved out to 4 AU. Model 256f again doubled the f resolution, but the maximum density achieved by the clumps increased only slightly, implying convergence in the sense of the hydrodynamical grid. However, because clump formation is driven by self-gravity, it is important to also increase the resolution of the Poisson solver. Model 256pf had doubled NYlm = 32 compared to model 256f and led to almost a factor of 10 increase in the maximum clump density. After this model had formed a few clumps, it was continued from that point with Nf = 512 and NYlm = 48. As shown in Figure 1, this model, 512pf, led to an even higher maximum clump density, about 27 times higher than in model 256pf and about 7 times higher than that encountered in the initial velocity-damped model 64d. The clump shown in Figure 1 survived for two orbital periods (60 yr). As Nf and NYlm are increased, the clumps become better dened and reach much higher densities but still disappear if the calculation is run indenitely. This suggests that even much higher spatial resolution is needed to follow the clumps faithfully and that the formation of the clumps is not the result of insufcient spatial resolution. At the time shown in Figure 1, the most prominent clump has a mass of 5.2 MJup, compared to a Jeans mass of 0.02 MJup at the average density of the clump, showing how tightly gravitationally bound the clump is. The clump also has an average spherical radius of 0.2 AU, considerably smaller than the critical tidal radius of 1.2 AU for its mass, orbital distance, and the central protostar mass, implying stability to tidal disruption (Boss 1998a). The clumps free-fall time at its maximum density is about 0.1 yr, a small fraction of its orbital

Fig. 1.Midplane density contours for the highest spatial resolution model (512pf) at 374 yr. Contours represent factors of 2 change in density. Diameter of region shown is 40 AU; inner hole has a diameter of 8 AU. Filaments, clumps, and holes are evident. The maximum density of 4.3 # 107 g cm3, which is particularly well-dened, occurs in the clump near 12 oclock, survives for two orbits, and is likely to contract toward planetary densities in a calculation with even higher spatial resolution.
period of about 30 yr, implying that the clump should be able to contract toward protoplanetary densities if it were allowed to by the spatial resolution of the grid. The clumps orbital eccentricity is e 0.3, suggesting that disk instability could lead to protoplanets on initially noncircular orbits.

4. CONCLUSIONS

The models show that strongly self-gravitating clumps can form in a marginally gravitationally unstable disk even without velocity damping. While the clumps appear to be capable of forming gas giant protoplanets, their subsequent orbital evolution and survival are highly uncertain. The models suggest that a phase of disk instability is likely to be a relatively chaotic environment for planet formation. If several clumps survive and reach planetary densities, the resulting gas giant planets may suffer mutual close encounters, resulting in even higher orbital eccentricities. Disk instability leading to the formation of gas giant protoplanets thus seems to be possible in a disk with a mass (0.1 M,) comparable to the upper end of the range inferred for the solar nebula and other protoplanetary disks (e.g., Beckwith et al. 1990). The surface density at 5 AU in the model is comparable to that favored in standard core accretion models (Pollack et al. 1996). However, in order to avoid loss of the gaseous envelope by a hydrodynamic blowoff, core accretion may require that the planet be embedded in a disk with a density at 5 AU of 109 g cm3 (Wuchterl, Guillot, & Lissauer 2000), a density that is at least 5 times higher than that of the disk instability models presented here and comparable to previous disk instability models (Boss 1997, 1998a). Because the timescale for disk instability is roughly times shorter than that for core accretion, clearly if it can occur, a

disk instability will circumvent core accretion in the solar nebula or in disks orbiting around other protostars. However, it would be premature to attempt to decide which of these mechanisms is superior at forming gas giant planets, based on our present knowledge. For example, the models presented here assume locally isothermal thermodynamics, an assumption that errs strongly in favor of clumping and merits considerable further scrutiny (Boss 1997; Pickett et al. 1998, 2000a; Nelson, Benz, & Ruzmaikina 2000). The survival of the clumps found here is likely to depend to a large extent on a proper treatment
of their thermodynamical evolution, which will be the subject of a future paper. I thank Patrick Cassen, Stephen Kortenkamp, and Brian Pickett for discussions and the referee for improvements to the manuscript. Supported in part by the NASA Planetary Geology and Geophysics Program under grant NAG5-3873. Calculations were performed on the Carnegie Alpha Cluster, which is supported in part by NSF MRI grant AST 99-76645.
REFERENCES Bally, J., Testi, L., Sargent, A., & Carlstrom, J. 1998, AJ, 116, 854 Beckwith, S. V. W., Sargent, A. I., Chini, R. S., & Gusten, R. 1990, AJ, 99, 924 Bodenheimer, P., Grossman, A. S., DeCampli, W. M., Marcy, G., & Pollack, J. B. 1980, Icarus, 41, 293 Bofn, H. M. J., Watkins, S. J., Bhattal, A. S., Francis, N., & Whitworth, A. P. 1998, MNRAS, 300, 1189 Boss, A. P. 1996, ApJ, 469, 906. 1997, Science, 276, 1836. 1998a, ApJ, 503, 923. 1998b, Annu. Rev. Earth Planet. Sci., 26, 53 Boss, A. P., & Myhill, E. A. 1992, ApJS, 83, 311 Bryden, G., Chen, X., Lin, D. N. C., Nelson, R. P., & Papaloizou, J. C. B. 1999, ApJ, 514, 344 Burrows, A., Guillot, T., Hubbard, W. B., Marley, M. S., Saumon, D., Lunine, J. I., & Sudarsky, D. 2000, ApJ, 534, L97 Butler, R. P., Marcy, G. W., Fischer, D. A., Brown, T. W., Contos, A. R., Korzennik, S. G., Nisenson, P., & Noyes, R. W. 1999, ApJ, 526, 916 Cameron, A. G. W. 1978, Moon Planets, 18, 5 Cassen, P. M., Smith, B. F., Miller, R. H., & Reynolds, R. T. 1981, Icarus, 48, 377 Charbonneau, D., Brown, T. M., Latham, D. W., & Mayor, M. 2000, ApJ, 529, L45 Goldreich, P., & Tremaine, S. 1980, ApJ, 241, 425 Guillot, T., Gautier, D., & Hubbard, W. B. 1997, Icarus, 130, 534 Henry, G. W., Marcy, G. W., Butler, R. P., & Vogt, S. S. 2000, ApJ, 529, L41 Jayawardhana, R., Hartmann, L., Fazio, G., Fisher, R. S., Telesco, C. M., & Pina, R. K. 1999, ApJ, 521, L129 Kuiper, G. P. 1951, Proc. Natl. Acad. Sci., 37, 1 Laughlin, G., & Bodenheimer, P. 1994, ApJ, 436, 335 Lin, D. N. C., Bodenheimer, P., & Richardson, D. C. 1996, Nature, 380, 606 Lin, D. N. C., & Papaloizou, J. C. 1986, ApJ, 309, 846 Lissauer, J. J. 1987, Icarus, 69, 249 Marcy, G. W., Cochran, W. D., & Mayor, M. 2000, in Protostars & Planets IV, ed. V. Mannings, A. P. Boss, & S. S. Russell (Tucson: Univ. Arizona Press), 1285 Miyoshi, K., Takeuchi, T., Tanaka, H., & Ida, S. 1999, ApJ, 516, 451 Mizuno, H. 1980, Prog. Theor. Phys., 64, 544 Nelson, A. F., Benz, W., Adams, F. C., & Arnett, W. D. 1998, ApJ, 502, 342 Nelson, A. F., Benz, W., & Ruzmaikina, T. V. 2000, ApJ, 529, 357 Papaloizou, J. C., & Savonije, G. J. 1991, MNRAS, 248, 353 Pickett, B. K., Cassen, P., Durisen, R. H., & Link, R. 1998, ApJ, 504, 468. 2000a, ApJ, 529, 1034 Pickett, B. K., Durisen, R. H., Cassen, P., & Meija, A. C. 2000b, in Poster Proc., IAU Symp. 200, Birth and Evolution of Binary Stars, ed. R. Mathieu & H. Zinnecker (Potsdam: Astrophys. Inst.), 190 Pollack, J. B. 1984, ARA&A, 22, 389 Pollack, J. B., Hubickyj, O., Bodenheimer, P., Lissauer, J. J., Podolak, M., & Greenzweig, Y. 1996, Icarus, 124, 62 Shu, F. H., Tremaine, S., Adams, F. C., & Ruden, S. P. 1990, ApJ, 358, 495 Slattery, W. L., DeCampli, W. M., & Cameron, A. G. W. 1980, Moon Planets, 23, 381 Stevenson, D. J. 1982, Planet. Space Sci., 30, 755 Tanaka, H., & Ida, S. 1999, Icarus, 139, 350 Toomre, A. 1964, ApJ, 139, 1217 Trilling, D. E., Benz, W., Guillot, T., Lunine, J. I., Hubbard, W. B., & Burrows, A. 1998, ApJ, 500, 428 Ward, W. R. 1997a, ApJ, 482, L211. 1997b, Icarus, 126, 261 Weidenschilling, S. J. 1977, Ap&SS, 51, 153 Wetherill, G. W. 1990, Annu. Rev. Earth Planet. Sci., 18, 205. 1996, Icarus, 119, 219 Wolk, S. J., & Walter, F. M. 1996, AJ, 111, 2066 Wuchterl, G., Guillot, T., & Lissauer, J. J. 2000, in Protostars & Planets IV, ed. V. Mannings, A. P. Boss, & S. S. Russell (Tucson: Univ. Arizona Press), 1081