Free Access
Volume 528, April 2011
Article Number A86
Number of page(s) 7
Section Planets and planetary systems
Published online 04 March 2011

© ESO, 2011

1. Introduction

In the late stages of the Solar System’s formation, the four Jovian planets had accumulated in the highly dissipative solar nebular gas and formed in well-separated, nearly circular, and coplanar orbits (Pollack 1996). Since the accretion process of the solid cores of the giant planets is not likely to have been 100% efficient, there was a residual planetesimal disk of mass 10–100   M left in the outer Solar System (Murray-Clay & Chiang 2006). By either accreting or gravitationally scattering the small bodies away, the four Jovian planets experienced orbital migrations caused by energy and angular momentum exchange. This scenario was first brought to light in the accretion models of Uranus and Neptune by Fernández & IP (1984).

The Kuiper belt (KB) with most of its observed features is believed to bear witness to the orbital migration of Jovian planets. Malhotra (1993, 1995) have proposed that Pluto’s eccentric and resonance-locked orbit could be the result of planetary migration. As Neptune migrated outwards, its sweeping 2:3 mean motion resonance (MMR) captured Pluto and excited its eccentricity. This planetary migration/resonance capture mechanism is not just for Pluto, but can also explain the orbital configurations of the Kuiper belt objects (KBOs) located in 4:5, 3:4, 2:3, 3:5, 4:7, and 1:2 MMRs with Neptune. Besides this, the observation shows that the classical KBOs (CKBOs) between 42 and 48 AU from the Sun can be classified into two major components: the dynamically cold objects with small inclinations and the hot ones with inclinations even greater than 30° (Brown 2001). These two populations have different physical properties, such as their sizes and colors (Trujillo & Brown 2002; Levison & Stern 2001; Peixinho et al. 2008).

Generally speaking, the hot and cold CKBOs could be considered to originate in different regions. Gomes (2003) suggested that the hot population may come from the inner region of the primordial planetesimal disk. While Neptune migrated outwards, a lot of planetesimals around 25AU were scattered to the current KB and attained high inclinations by close encounters with Neptune. In a subsequent paper by Gomes et al. (2004), it was shown that an external edge of the disk at  ~35 AU is also required to keep Neptune from migrating farther out than  ~30 AU for reasonable disk masses. In addition, a portion of planetesimals formed within 35 AU were captured by Neptune’s 1:2 MMR and pushed outwards. Owing to the stochastic effect in Neptune’s outward migration, some of them leaked out of the 1:2 MMR and formed the cold CKBOs. The mass of the trapped 1:2 MMR particles is a key factor in producing these low-eccentricity bodies in the KB region (Levison & Morbidelli 2003).

Moreover, other evidence of the unadulterated imprint of planetary migration can be seen in the Solar System. Since Uranus and Neptune could not form in situ on a reasonable timescale (Safronov & Zvjagina 1969; Lissauer & Stewart 1993; Lissauer et al. 1995; Levison & Stewart 2001), Thommes et al. (1999) and Tsiganis et al. (2005) suggest that these two planets may form much closer to the Sun. Since they were gravitationally scattered outwards after Jupiter and Saturn had accreted a large amount of nebular gas, Uranus and Neptune migrated outwards and evolved into their present orbital architectures. In addition, the planetary migration mechanism could also be responsible for the origins of the main asteroid belt (Minton & Malhotra 2009), trojan asteroids of Jupiter and Neptune (Morbidelli et al. 2005; Li et al. 2007; Nesvorný & Vokrouhlický 2009), irregular satellites (Nesvorný et al. 2007), the Oort cloud (Gomes 1998), and the late heavy bombardment period of the terrestrial planets (Gomes et al. 2005). Therefore, studying the orbital migration of Jovian planets could offer a lot of valuable clues to the formation and primordial evolution of our Solar System.

Following the pioneering work by Fernández & IP (1984), many studies have been done for better understanding of the planetary migration process. In a recent paper, Kirsh et al. (2009) drew a picture in which a single planet migrates in a massive planetesimal disk. The planetesimals in the inner encounter zone, with z-component of the angular momentum smaller than that of the planet (Hp), would drive the planet inwards via close encounters, and vice versa for the planetesimals in the outer encounter zone where H > Hp. Due to the shorter timescale of the scattering process in the inner encounter zone, more mass would be scattered outwards and gain angular momentum from the planet. With the imbalance of such a scattering bias, the planet tends to migrate inwards on average. In Gomes et al. (2004), it is shown that the situation can change when the mass density of the outer disk is much higher than in the inner one. Furthermore, if the four Jovian planets act together, since Uranus depletes a large amount of material in Neptune’s inner encounter zone, the mass in Neptune’s outer encounter zone would tend to be much higher. This mass asymmetry could overwhelm the above scattering bias and lead to the outward migration of Neptune. The same reason can also be applied to Uranus and Saturn, which on average acquire orbital angular momentum from the planetesimals and migrate outwards, while as the ultimate source of the angular momentum and energy, Jupiter moves sunward.

As we know, in numerical simulations performed so far of the evolution of planets embedded in a planetesimal disk (Hahn & Malhotra 1999; Gomes 2003; Gomes et al. 2004; Tsiganis et al. 2005), two important physical effects are not included: (1) all small planetesimals are treated as point masses, thus the possibility of collision or accretion phenomenon during the planetary migration has been neglected; (2) mutual gravitational forces among planetesimals are ignored, while the self-stirring may be able to excite planetesimals more distant in the disk into high-eccentricity orbits and cause them to enter Neptune’s scattering region. Therefore, more investigation of the general migration process is warranted.

The aim of our work is to explore both the evolution of the four Jovian planets in a more realistic planetesimal disk with self-gravity and the accretion of proto-Uranus and proto-Neptune taking place simultaneously. The rest of this paper is organized as follows. We start in Sect. 2 with a brief estimation of the growth timescale for KBO-sized planetesimals in the region of planetesimal disk that we are concerned with. Section 3 explores a model where Neptune starts its journey about 7 AU inside the present location and determines some physical parameters. Since Neptune’s mass hardly increases due to planetesimal accretion as a result of Sect. 3, we consider a more compact planetary configuration where Neptune formed even closer to the Sun in Sect. 4. Here, we discuss, in addition, how the disk self-gravity affects the migration process. In Sect. 5, we compare the distribution of surviving planetesimals gained in Sect. 4 with the observational data of the KBOs. The conclusions and discussion are given in Sect. 6.

2. Formation of KBO-sized planetesimals

As the time of writing, more than 1000 KBOs have been discovered1 beyond Neptune’s orbit. Most of them have diameters (D) larger than 100 km, while there are roughly a dozen giant objects with diameters of about 1000 km or more, such as Pluto (D ≈ 2320 km) and 2003 UB313 (D ≈ 2400 km). Previous studies show that the current belt mass between 30 and 50 AU of less than 0.1   M could not produce the KBOs in situ (Stern 1996; Stern & Colwell 1997; Kenyon & Luu 1998; Kenyon & Luu 1999). Thus, the KBOs should form in a denser region much closer to the Sun and were subsequently transported outwards to their present locations (Levison & Morbidelli 2003; Gomes 2003). Now comes the question, whether the KBOs could grow to their present sizes in the protoplanetary disk within 30 AU on a reasonable timescale?

According to the core accretion model in the gaseous disk, the timescale τgrow for a body to reach a mass M can be described by (Ida & Lin 2004) (1)where the scaling parameters fd = fg = 0.7 with respect to the minimum mass solar nebula model (Hayashi 1981), and the volatile-ice enhancement ηice = 4.2 outside the snow line at aice ≈ 2.7 AU.

thumbnail Fig. 1

The timescales τgrow for forming planetesimals with different diameters D (100, 200 and 2000 km) at various locations a. The density of these planetesimals is assumed to be 2 g cm-3 as Pluto’s. The dashed line indicates the gas depletion timescale  ~107 yr.

By the above formula (1), we can estimate the timescale τgrow for the candidates of the KBOs with different sizes D at various locations a, as shown in Fig. 1. We can find that there is no problem in forming typical KBO-sized objects between 10 and 30 AU on a timescale  ≤ 107 yr prior to severe gas depletion. Therefore, the planetesimal disk considered below may not only be the fuel for the migration of Jovian planets but also the source of the KBOs. As we know, the high-inclination CKBOs have relatively large sizes (Levison & Stern 2001), thus forming in the inner region of the planetesimal disk offers a good explanation (Fig. 1). It should be pointed out that the protoplantary disk could be more massive than that in the minimum mass solar nebula model. This tends to overestimate the growth timescale of an isolated planetesimal and hence the Pluto-sized objects with D ≈ 2000 km may also come into being on a reasonable timescale within 30 AU.

3. Migration simulations in model i

The Solar System in our numerical simulations consists of the Sun with the masses of four terrestrial planets added, four Jovian planets, and a self-gravitating planetesimal disk. In this section, the setup of the experiment is similar to the simulations in Hahn & Malhotra (1999). The initial orbital semimajor axes are chosen to be 5.4, 8.7, 16.3, and 23.2 AU for Jupiter, Saturn, Uranus, and Neptune, respectively. The planets are assumed to have reached their present sizes. The initial planetesimal disk has a mass of 30   M and extends from 10 to 30 AU, so its inner edge lies just outside the outer feeding zone (~3.5 Hill radii) of Saturn, and the outer edge is placed at Neptune’s current position. In the simulations of this paper, the disk surface density varies as ak, and the index k is set to be  − 1, as normally used in most works (Gomes et al. 2004; Tsiganis et al. 2005; Kirsh et al. 2009, etc.). Supposing that the disk is composed of N equal-mass particles, we perform three bins of simulations below with N = 1000, 2000, 3000 (considering the computational power), and individual disk planetesimals have masses of m = 0.03,0.015,0.01   M. The spacing between planetesimals is at least 7 Hill radii. The initial eccentricities and inclinations of the particles are set to be e = i = 0.001. All other angles of an orbit (longitude of ascending node Ω, longitude of perihelion ϖ, and mean longitude λ) are chosen randomly in the range 0–2π.

We numerically simulate the evolution of the above systems using the code NBODY4 (Aarseth 2003), which considers the total gravitational forces due to the Sun, planets, and the planetesimals. NBODY4 is an Hermite individual time step code and runs with the special-purpose GRAPE-type supercomputer hardware, and it can strongly increase the speed of the force calculation. This procedure does not rely on softening of the force, close two-body encounters are handled by the Kustaanheimo-Stiefel regularization method, and strong interactions involving three or more bodies are selected for chain regularization.

thumbnail Fig. 2

The temporal evolution of Neptune’s semimajor axis for the planetesimal disk with a total mass of 30   M and a surface density decaying as a-1 between 10 and 30 AU. The planetesimal disks are assumed to be composed of N = 1000 (green), 2000 (red) and 3000 (blue) equal-mass particles, respectively.

We first focus on Neptune’s orbital migration since it is correlated with a series of puzzles such as the formation of the KB (Malhotra 1993, 1995). Figure 2 shows the typical temporal evolution of Neptune’s semimajor axis as it scatters the surrounding disk particles during the 50 Myr runs with different N. At the very beginning, Neptune experiences a short period of inward migration as in the scenario of a single planet embedded in the planetesimal disk. After Uranus gradually clears the planetesimals with the angular momentum H < HNeptune inside Neptune’s orbits, Neptune tends to move outwards instead.

The final semimajor axis of Neptune may serve as a very important diagnostic for our models. For the case of N = 1000, Neptune’s outward motion rapidly slows down, and it can only reach a distance as far as  ~27 AU that is well within the outer edge of the disk. This outcome is very similar to what is obtained in Hahn & Malhotra (1999). When we increase the number N of planetesimals to 2000 and 3000, we observe that Neptune could maintain its migration a few AU farther and then stops at  ~30 AU. For slightly different initial condition for the disk in each N bins, the results are equivalent to Fig. 2. We may understand the different behavior of Neptune’s migration related to the resolution of the disk according to the following. For a fixed disk mass, a smaller number of planetesimals induces more pronounced and stochastic density fluctuations in the disk, therefore affecting the orbital distribution of the planetesimals and resulting in slow migration of planets (Murray-Clay & Chiang 2006). While the disk consists more planetesimals, the stochastic effects could be effectively averaged in space and time, and the final location of Neptune could converge (Gomes et al. 2004). Therefore, we deem that the results obtained in the cases of N = 2000,3000 are more reliable.

thumbnail Fig. 3

Evolution of semimajor axes of the four Jovian planets embedded in the planetesimal disk as in Fig. 2, but only for the case of N = 3000.

On the basis of the above results, we do not go on to discuss the low-resolution case of N = 1000. Figure 3 summarizes the variations in the semimajor axes of the four planets with respect to time for the case of N = 3000. As might be expected, Neptune, Uranus, and Saturn migrate outwards, while Jupiter moves slightly sunward. The final state resembles the current configuration of the outer Solar System very well, and a similar result can also be observed in the case of N = 2000. However, although the case of N = 3000 may make the planetary migration more realistic, a huge amount of CPU time (due to the many more close encounters that delay the code) is needed to fulfill this task. As a result, we think that the disk composed of 2000 planetesimals (m = 0.015   M, around three time Pluto’s mass) is the most suitable model, because it is not only economical but can also be used to describe the evolution history of planets.

The observed distribution of the resonant KBOs provides a strong constraint on the migration timescale of Jovian planets. Ida et al. (2000) propose that Neptune would be able to capture planetesimals into its 2:3 MMR if it migrates on a timescale  ≥ (106) yr. Although the region outside the location of Neptune’s 1:2 MMR is empty initially, many planetesimals would be scattered outwards to this virgin territory during the planet migration period. To prevent the capture of planetesimals into the 1:2 MMR, Neptune’s migration timescale should also be  ≤ (107) yr (Ida et al. 2000). In Fig. 3 we can see that the migration timescale of planets is about 107 yr, which is consistent with the above arguments and with the results presented in some other previous works (Hahn & Malhotra 1999; Li et al. 2006).

Moreover, we have examined the macro-accretion process of Uranus and Neptune in our numerical simulations. Occurring simultaneously with the migration driven by gravitational scattering of the planetesimals, additional solid material may be injected into the accretion zones of Uranus and Neptune. In the code NBODY4, the coalescence between planets and planetesimals is computed carefully through their physical collisions. In all runs, Uranus and Neptune merely attract a few planetesimals, and their masses increase by about 0.1   M. We find that the pairwise collision events almost happened in the first several million years, but that they possess an extremely low rate. These outcomes are qualitatively different from those of Fernández & Ip (1984), which has some limitations, such as allowing bodies to collide if they were on mutual crossing orbits by using Öpik’s two-body formulation (Öpik 1951), and neglecting the perturbation of the Sun during close encounters, etc (see Brunini & Melita 2002 for a detailed discussion). To confirm that the accretion process is much less efficient, we made some short runs starting with proto-Uranus and proto-Neptune having one fifth of their present sizes, just as Fernández & Ip (1984) did. As expected, the results of these runs do not show any noticeable increase in the accretion efficiency. Although there are not enough planetesimals modeling the disk to reveal the realistic migration/accretion process with a completely self-consistent calculation, we do not imagine there would be any significant accretion by Uranus or Neptune in this final stage. In fact, Matter et al. (2009) evaluated the enrichment of the giant planets on the basis of the “Nice model”, and given certain probability of impact, they find that the amount of solid material accreted by these two ice giants is even less than 0.1   M. Therefore, we believe that Uranus and Neptune should form even closer to the Sun (Shu et al. 1993; Ford & Chiang 2007), maybe well inside 18 AU, and were subsequently transported outwards (Thommes et al. 1999; Tsiganis et al. 2005). Consequently, a more compact model is investigated in the next section.

4. Migration simulations in model ii

We now investigate the planetary migration process in a more compact model. Keeping the initial positions of Jupiter and Saturn the same as in Sect. 3, we put Uranus and Neptune much closer to the Sun at distances of 15.5 AU and 17.8 AU, respectively. The planetesimal disk of 30   M is still modeled by 2000 equal mass particles distributed as a-1 from 10 to 30 AU. Since Neptune and Uranus were initially so close, they would strongly perturb each other so their eccentricities would be stirred. Due to the dynamical friction, these two ice giants would evolve into near circular orbits at the expense of the surrounding planetesimals, which obtained more eccentric and inclined orbits. This may roughly resemble the heating process of the planetesimal disk, just after Neptune and Uranus were violently thrown outwards from the inner region by Jupiter or Saturn (Thommes et al. 1999; Tsiganis et al. 2005). Hereafter, we do not discuss the formation history of Uranus and Neptune’s cores in detail (Wuchterl et al. 2000), but only focus on their migration behavior, which started from such closely packed orbits.

thumbnail Fig. 4

Evolution of semimajor axes of the four Jovian planets due to a planetesimal disk of 30   M modeled by 2000 equal mass particles distributed as a-1 from 10 to 30 AU. The initial positions of Uranus and Neptune are 15.5 AU and 17.8 AU, respectively.

Figure 4 shows the profile of the planet migration from the more compact configuration. In this case, Neptune starts more than 5 AU closer to the Sun than in Sect. 3. It undergoes a rapid outward migration in less than 2 × 107 yr, and then its motion slows down. The location of Neptune is found to increase quasi-asymptotically to  ~26 AU, which is well within the outer edge of the planetesimal disk. The result is qualitatively equivalent to the phenomenon called “damped migration” observed in Gomes et al. (2004) for a disk without self-gravity. It should be noted that, at the end of integrations, planetesimals with a < 36 AU (i.e., within 10 AU beyond Neptune’s orbit) have nearly been depleted. Therefore, the outward motion of Neptune has ceased because of the empty encounter zone.

To check whether Neptune could persist its migration until it reaches the current position, we extended the outer edge of the disk outside 30 AU. Keeping the surface density and each particle’s mass constant, we displaced the outer edge of the disk at 35 AU so that 500 more equal-mass planetesimals with a total mass of 7.5   M were added between 30 and 35 AU. With this extended planetesimal disk, we then performed five new runs (hereafter RUN35) with slightly varied initial conditions. We keep up the integration until Neptune’s migration has nearly stopped, i.e., the semimajor axis quasi-asymptotically converges at some certain value.

thumbnail Fig. 5

As in Fig. 4, but the outer edge of the planetesimal disk is extended to 35 AU and 500 extra equal-mass planetesimals with a total mass of 7.5   M distributed as a-1 are added.

Three runs of RUN35 were unsuccessfully completed when the system became unstable, e.g., Uranus is ejected from the system by Saturn, whose eccentricity is excited as a result of the 2:5 MMR crossing with Jupiter, while in the other two runs, the system could evolve to a configuration that is very similar to the present outer Solar System, and an example of the variations in the planets’ semimajor axes is shown in Fig. 5. Neptune first experiences a rapid radial displacement to about 25 AU in about 2 × 107 yr, and then slowly migrates, tending towards its current position near 30 AU. This outcome is consistent with our perspective but differs markedly from the simulations by Gomes et al. (2004), in which the planetesimal disk (up to 50 AU) has the same surface density variation but is a bit different from our disk’s mass. We infer that the qualitative difference stems from including mutual gravitational interactions among planetesimals themselves in our model. It is reasonable to imagine that, since tens of Earth masses of material have been scattered outwards, the planetesimals beyond Neptune’s orbit would shift to more eccentric orbits. Besides this, the local self-stirring would also make positive contributions to the eccentricity excitation (Kenyon et al. 2008). Consequently, some particles that are a few AU farther out in eccentric orbits could evolve down to perihelia as low as in Neptune’s outer encounter zone, and then they may possibly suffer the perihelion instability and become the fresh replenishment to transfer more angular momentum to Neptune and displace it farther away. As a matter of fact, each RUN35 is quite “expensive”, and more than several months of CPU time is needed, so we did not make additional simulations, but find that the main features in Fig. 5 are robust.

Indeed, a straightforward way to illustrate the difference between our RUN35 and Gomes’ model (Gomes et al. 2004) is to perform some extra runs without the self-gravity in our disk model. In a comparison of runs with initial conditions the same as RUN35, for convenience, we employed the MERCURY6 N-body integrator (Chambers 1999) to track the orbital evolution of the Jovian planets in a non-self-gravity planetesimal disk. In all these experiments, just as in Gomes et al. (2004), Neptune’s outward migration motion would stop inside 25 AU and the part of disk a few AU beyond Neptune’s orbit preserves its initial mass and cold status. Once a self-gravitating planetesimal disk is taken into account, some of these cold particles could suffer eccentricity excitation and enter the outer encounter zone, and Neptune’s outward migration would continue.

Figure 5 shows that Neptune expands its orbit slightly out to 30 AU at 2 × 108 yr, suggesting that the parameters in our model still need to be adjusted delicately. However, it is worth noting that the planetesimal disk inside 36 AU have almost been depleted by the end of the integration (Fig. 6), implying that the displacement of Neptune’s position should have stopped.

5. Implication for the Kuiper belt

As Neptune migrated by exchange of energy and angular momentum with planetesimals, some of these small objects were scattered outwards and entered the KB region. Figure 6 shows the distributions of semimajor axes, eccentricities, and inclinations for the surviving objects between 30–100 AU after 2 × 108 yr integration for the example of RUN35 shown in Fig. 5. Since Neptune stops at  ~31 AU in this run, while its current location is near 30 AU, for comparison with the real KB region, planetesimals’ semimajor axes are all scaled by a factor of 30/31. The solid curve in the upper panel stands for the Neptune-crossing orbit limit where planetesimals’ perihelia q = aNeptune.

thumbnail Fig. 6

The distributions of semimajor axes, eccentricities and inclinations for the surviving planetesimals between 30–100 AU after 2 × 108 yr integration for the migration case shown in Fig. 5.

In Fig. 6, we can see some apparent features of the surviving planetesimals between 30–50 AU. (1) The region between 30–36 AU is nearly unoccupied. (2) Some objects are concentrated in the 2:3 MMR with Neptune, while only a small proportion of surviving particles fall into the 1:2 resonance. (3) There are objects with low eccentricities e < 0.1 in the inner classical KB region from 42–45 AU, and a deficiency of objects with e < 0.1 in the outer region from 45–48 AU. (4) Dynamically “cold” (small i) and “hot” (large i) populations coexist in the classical KB region. (5) A drop occurs in number density near Neptune’s 1:2 MMR, and objects beyond all move on highly eccentric orbits (e > 0.3) as scattered disk objects. These outcomes show that the overall picture looks like the KB. However, there are not enough bodies to do proper statistical comparisons with the real KB as in Levison et al. (2008), the resolution of primordial disk in our simulations is too low to see whether self-gravity plays a key role, and the low resolution may also inevitably lead to an exaggerated role for close encounters, which are probably critical for the formation of the remnant belt. These situations need to be studied in the future.

As mentioned in the introduction, the CKBOs show a correlation between their colors and inclinations (Trujillo & Brown 2002). Later, Peixinho et al. (2008) re-examined this correlation and proposed a color break at i ≈ 12°: hot CKBOs with i > 12° are homogeneously blue, whereas cold CKBOs with i < 12° exhibit red colors. Distinct distributions of color and inclination suggest that the hot population would consist of planetesimals that formed closer to the Sun than the cold one. Our simulation’s results show that most high-inclination objects originate in the region from 18 AU to 24 AU, while objects with small inclinations come from the location outside 24 AU. Such bi-modal inclination distribution of CKBOs is consistent with their observed properties. Besides, as discussed in Sect. 2, the high-inclination CKBOs formed in the inner region of the planetesimal disk would grow larger on the same timescale.

In RUN35, the total number of objects that survived in the region 30–50 AU is about 3% corresponding to the initial disk that contains 30   M material, indicating the total residual mass is about 1 M. As shown in Fig. 6, the non resonant objects above the solid line (q < aN) would experience close encounters with Neptune and leave the KB, and the 40–42 AU region is particularly unstable due to the presence of ν8 secular resonance. Inferred from Gomes (2004), only about 20% of objects (i.e.,  ~0.2   M) may persist over the Solar System’s age. Nevertheless, this estimated mass of the remnant belt is still a bit higher than the one we observed (Chiang & Brown 1999; Trujillo et al. 2001; Gladman et al. 2001). Unfortunately, considering the total systematic error in the integration, a complete exploration of a time evolution that is as long as the Solar System’s age does not seem reliable for our code.

6. Conclusions and discussion

The existence of the KB, the main asteroid belt, trojan asteroids of Jupiter and Neptune, the inefficiency in planet formation, and the Oort cloud, all suggest that a residual planetesimal disk remained in the outer Solar System following the birth of Jovian planets. Thanks to the gravitational scattering of these small bodies, Neptune, Uranus, and Saturn gained orbital angular momentum and migrated outwards, while as the ultimate source of the angular momentum and energy, Jupiter moved sunward. In this paper we have investigated the orbital evolution of the four Jovian planets embedded in a self-gravitating planetesimal disk, and the simultaneous accretion of planetesimals by proto-Uranus and proto-Neptune. Different initial conditions were discussed in our work and several important conclusions can be summarized as follows.

  • (1)

    The accretion of small planetesimals by proto-Uranus and proto-Neptune through physical collisions was very inefficient by about 0.1   M, in contrast to those previously found through the Öpik’s formalism (Fernández & IP 1984). Thus, we believe that Uranus and Neptune had already reached their present sizes before the planetary migration era, most likely in a region inside 18 AU.

  • (2)

    Considering mutual gravitational forces among the planetesimals themselves, a fraction of the outer disk objects that cannot get close to Neptune will be dynamically heated to eccentric orbits by the planet-scattered planetesimals and local disk self-stirring. When these distant planetesimals enter the outer encounter zone of Neptune, they would transfer more angular momentum to Neptune and drive it farther away. If the disk self-gravity is ignored, the radial displacement of Neptune is smaller, and the outer part of the disk discussed before will preserve its original status and mass.

  • (3)

    In the preferred simulations where we start Neptune from 17.8 AU, it could migrate outwards to the present location near 30 AU in a self-gravitating planetesimal disk, which has a mass of 35   M and is distributed as a-1 from 10 to 35 AU. At the end of integration, Neptune’s outward migration has ceased because the region inside 36 AU has nearly been depleted, and the distribution of the surviving planetesimals is analogous to the KB. However, under the same initial conditions, but without the disk self-gravity, Neptune would stop at  ~25 AU, since its outer encounter zone is empty and cannot be replenished by new planetesimals in the region a few AU away.

In view of the above, we conclude that the planetesimal disk’s self-gravity should play a non-negligible role in the formation and primordial evolution of the outer Solar System; however, there are several questions about the primordial disk adopted in our work. The main one is that the planetesimals are too big: even for N = 3000 the particles are Lunar-sized. This means that close encounters between planets and planetesimals would be far more severe than in a more realistic disk, thus the stochastic effect in Neptune’s orbital migration could be much more prominent; besides, close encounters between planetesimals are also too strong, and they produce dynamical “heating” to a greater extent and rate. All these would lead to the replenishment of Neptune’s outer feeding zone, as well as the final distribution of the remnant “KB”. The second is that we start the disk at 10 AU, while the region inside the orbit of the outermost planet (Neptune) is more likely to be empty when the giant planets are fully formed. Embedded in the disk with an inner edge beyond Neptune, the Jovian planets can also evolve from some certain configurations to their current stable orbits, and similar results for the KB could be produced (Tsiganis et al. 2005; Morbidelli et al. 2007). The last one is that the planet migration driven by planetesimal scattering depends on the disk density, while this parameter is kept constant in our numerical experiments. All these issues make it worth studying the consequences of the different set-ups of disk in the future. However, we argue that our key findings for including the disk self-gravity would not alter.

Theoretically speaking, a planet embedded in a self-gravitating planetesimal disk can launch spiral density waves at its Lindblad resonances (Goldreich & Tremaine 1980). Since the disk inside Neptune’s orbit is strongly perturbed by Uranus and rapidly depleted, we are not concerned about the inner Lindblad resonance. As a result, the disk torques exerted on Neptune due to outer Lindblad resonances would tend to inhibit Neptune’s outward migration. However, this delay effect is not observed in any of our simulations, which is expected. As the planetesimals are somewhat large in our disk model, the disk is heated quickly (the root-mean-squared eccentricity  ⟨ e2 ⟩ 1/2 ≥ 0.25) in less than 1 Myr via self-gravity. Thus, the disk becomes too stirred up to sustain density waves at resonances. Indeed, the role of density waves could be better studied by simulating a disk composed of much more smaller planetesimals, but such a calculation exceeds our current computational resources. However, the resolution problem as mentioned above is crucial, which implies that some of the results might be different for higher resolutions.


This work was supported by the National Natural Science Foundation of China (Nos. 10833001, 11003008, 11073012 and 11078001) and the National Basic Research Program of China (No. 2007CB814800). We are grateful to Dr. Sverre Aarseth for supplying the numerical code and improving the manuscript. The authors would also like to express their thanks to the referee for his valuable comments.


  1. Aarseth, S. J. 2003, Gravitational N-Body Simulations (Cambridge: Cambridge Univ. Press) [Google Scholar]
  2. Brown, M. 2001, AJ, 121, 2804 [Google Scholar]
  3. Brunini, A., & Melita, M. D. 2002, MNRAS, 330, 184 [NASA ADS] [CrossRef] [Google Scholar]
  4. Chambers, J. E. 1999, MNRAS, 304, 793 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  5. Chiang, E. I., & Brown, M. E. 1999, AJ, 118, 1411 [NASA ADS] [CrossRef] [Google Scholar]
  6. Fernández, J. A., & Ip, W.-H. 1984, Icarus, 58, 109 [NASA ADS] [CrossRef] [Google Scholar]
  7. Ford, E. B., & Chiang E. I. 2007, ApJ, 661, 602 [NASA ADS] [CrossRef] [Google Scholar]
  8. Gladman, B., Kavelaars, J. J., Petit, J.-M., et al. 2001, AJ, 122, 1051 [NASA ADS] [CrossRef] [Google Scholar]
  9. Goldreich, P., & Tremaine, S. 1980, ApJ, 241, 425 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  10. Gomes, R. S. 1998, AJ, 116, 2590 [NASA ADS] [CrossRef] [Google Scholar]
  11. Gomes, R. S. 2003, Icarus, 161, 404 [Google Scholar]
  12. Gomes, R., Levison, H. F., Tsiganis, K., & Morbidelli, A. 2005, Nature, 435, 466 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  13. Gomes, R., Morbidelli, A., & Levison H. F. 2004, Icarus, 170, 492 [NASA ADS] [CrossRef] [Google Scholar]
  14. Hahn, J. M., & Malhotra, R. 1999, AJ, 117, 3041 [NASA ADS] [CrossRef] [Google Scholar]
  15. Hayashi, C. 1981, Prog. Theor. Phys. Suppl., 70, 35 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  16. Ida, S., & Lin, D. N. C. 2004, ApJ, 604, 388 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  17. Ida, S., Bryden, G., Lin, D. N. C., & Tanaka, H. 2000, ApJ, 534, 428 [NASA ADS] [CrossRef] [Google Scholar]
  18. Kenyon, S. J., & Luu, J. X. 1998, AJ, 115, 2136 [NASA ADS] [CrossRef] [Google Scholar]
  19. Kenyon, S. J., & Luu, J. X. 1999, AJ, 118, 1101 [NASA ADS] [CrossRef] [Google Scholar]
  20. Kenyon, S. J., Bromley, B. C., OBrien, D. P., & Davis, D. R., 2008, in The Solar System Beyond Neptune, ed. M. A. Barucci, H. Boehnhardt, D. P. Cruikshank, & A. Morbidelli (Tucson: Univ. Arizona Press), 293 [Google Scholar]
  21. Kirsh, D. R., Duncan, M., Brasser, R., & Levison, H. F. 2009, Icarus, 199, 197 [NASA ADS] [CrossRef] [Google Scholar]
  22. Levison, H. F., & Stern, S. A. 2001, AJ, 121, 1730 [NASA ADS] [CrossRef] [Google Scholar]
  23. Levison, H. F., & Stewart, G. 2001, Icarus, 153, 1960 [Google Scholar]
  24. Levison, H. F., & Morbidelli, A. 2003, Nature, 426, 419 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  25. Levison, H. F., Morbidelli, A., VanLaerhoven C., Gomes R., & Tsiganis, K. 2008, Icarus, 196, 258 [NASA ADS] [CrossRef] [Google Scholar]
  26. Lissauer, J. J., & Stewart, G. R. 1993, in Protostars and Planets III, ed. E. H. Levy, & J. I. Lunine (Tucson: Univ. Arizona Press), 1061 [Google Scholar]
  27. Lissauer, J. J., Pollack, J. B., Wetherill, G. W., & Stevenson, D. J. 1995, in Neptune and Triton, ed. D. P. Cruikshank, M. S. Matthews, & A. M. Schumann (Tucson: Univ. Arizona Press), 37 [Google Scholar]
  28. Li, J., Zhou, L.-Y., & Sun, Y.-S. 2006, ChJAA, 6, 588 [Google Scholar]
  29. Li, J., Zhou, L.-Y., & Sun, Y.-S. 2007, A&A, 464, 775 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Malhotra, R. 1993, Nature, 365, 819 [CrossRef] [Google Scholar]
  31. Malhotra, R. 1995, AJ, 110, 420 [NASA ADS] [CrossRef] [Google Scholar]
  32. Matter, A., Guillot, T., & Morbidelli, A., 2009, Planet. Space Sci., 57, 816 [NASA ADS] [CrossRef] [Google Scholar]
  33. Minton, D. A., & Malhotra, R. 2009, Nature, 457, 1109 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  34. Morbidelli, A., Levison, H. F., Tsiganis, K., & Gomes, R. 2005, Nature, 435, 462 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  35. Morbidelli, A., Tsiganis, K., Crida, A., Levison, H. F., & Gomes, R. 2007, AJ, 134, 1790 [NASA ADS] [CrossRef] [Google Scholar]
  36. Murray-Clay, R. A., & Chiang, E. I. 2006, ApJ, 651, 1194 [NASA ADS] [CrossRef] [Google Scholar]
  37. Nesvorný, D., & Vokrouhlický, D. 2009, AJ, 137, 5003 [NASA ADS] [CrossRef] [Google Scholar]
  38. Nesvorný, D., Vokrouhlický, D., & Morbidelli, A. 2007, AJ, 133, 1962 [NASA ADS] [CrossRef] [Google Scholar]
  39. Öpik, E. J. 1951, Proc. Roy. Irish Acad. 54, 165 [Google Scholar]
  40. Peixinho, N., Lacerda, P., & Jewitt, D. 2008, AJ, 136, 1837 [NASA ADS] [CrossRef] [Google Scholar]
  41. Pollack J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62 [NASA ADS] [CrossRef] [Google Scholar]
  42. Safronov, V. S., & Zvjagina, E. V. 1969, Icarus, 10, 109 [Google Scholar]
  43. Shu, F. H., Johnstone, D., & Hollenbach, D. 1993, Icarus, 106, 92 [NASA ADS] [CrossRef] [Google Scholar]
  44. Stern, S. A. 1996, AJ, 112, 1203 [NASA ADS] [CrossRef] [Google Scholar]
  45. Stern, S. A., & Colwell, J. E. 1997, AJ, 114, 841 [NASA ADS] [CrossRef] [Google Scholar]
  46. Thommes, E. W., Duncan, M. J., & Levison, H. F. 1999, Nature, 402, 635 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  47. Trujillo, C. A., & Brown, M. E. 2002, ApJ, 566, L125 [NASA ADS] [CrossRef] [Google Scholar]
  48. Trujillo, C. A., Jewitt, D. C., & Luu, J. X. 2001, AJ, 122, 457 [NASA ADS] [CrossRef] [Google Scholar]
  49. Tsiganis, K., Gomes, R., Morbidelli, A., & Levison, H. F. 2005, Nature, 435, 459 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  50. Wuchterl, G., Guillot, T., & Lissauer, J. J. 2000, in Protostars and Planets IV, ed. V. Mannings, A. P. Boss, & S. S. Russel (Tucson: Univ. Arizona Press), 1081 [Google Scholar]

All Figures

thumbnail Fig. 1

The timescales τgrow for forming planetesimals with different diameters D (100, 200 and 2000 km) at various locations a. The density of these planetesimals is assumed to be 2 g cm-3 as Pluto’s. The dashed line indicates the gas depletion timescale  ~107 yr.

In the text
thumbnail Fig. 2

The temporal evolution of Neptune’s semimajor axis for the planetesimal disk with a total mass of 30   M and a surface density decaying as a-1 between 10 and 30 AU. The planetesimal disks are assumed to be composed of N = 1000 (green), 2000 (red) and 3000 (blue) equal-mass particles, respectively.

In the text
thumbnail Fig. 3

Evolution of semimajor axes of the four Jovian planets embedded in the planetesimal disk as in Fig. 2, but only for the case of N = 3000.

In the text
thumbnail Fig. 4

Evolution of semimajor axes of the four Jovian planets due to a planetesimal disk of 30   M modeled by 2000 equal mass particles distributed as a-1 from 10 to 30 AU. The initial positions of Uranus and Neptune are 15.5 AU and 17.8 AU, respectively.

In the text
thumbnail Fig. 5

As in Fig. 4, but the outer edge of the planetesimal disk is extended to 35 AU and 500 extra equal-mass planetesimals with a total mass of 7.5   M distributed as a-1 are added.

In the text
thumbnail Fig. 6

The distributions of semimajor axes, eccentricities and inclinations for the surviving planetesimals between 30–100 AU after 2 × 108 yr integration for the migration case shown in Fig. 5.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.