ABSTRACT
Recent imaging campaigns indicate the likely existence of massive planets (∼1–10 MJ) on ≳1000 year orbits about a few percent of stars. Such objects are not easily explained in most current planet formation models. In this Letter, we use ensembles of 100 N-body simulations to evaluate the potential for planet scattering during relaxation of dynamically active systems to explain the population of giant planets with projected separations up to a few 100 AU. We find that such a mechanism could indeed be at play, and that statistical samples of long-period planets could place interesting constraints on early stage planet formation scenarios. Results from direct imaging and microlensing surveys are complementary probes of this dynamical relaxation process.
Export citation and abstract BibTeX RIS
1. INTRODUCTION
The ongoing discovery of exoplanets is providing a wealth of opportunities to test models of planetary system formation, and has already revealed that our own solar system is not necessarily typical (see Udry & Santos (2007) for a review). Due to the natural biases of radial velocity and transit detection methods, a major focus has been on short-period planets of all ages. Direct imaging efforts are biased oppositely, toward the detection of long-period, young, and massive planets (see, e.g., Beuzit et al. (2007) for a review). Microlensing searches are also sensitive to longer period objects, independent of system age (see, e.g., Gaudi (2008) for a review). Recent direct imaging efforts have yielded constraints on the likelihood of certain classes of long-period gas-giant planets (Masciadri et al. 2005; Biller et al. 2007; Kasper et al. 2007; Nielsen et al. 2008), together with a number of increasingly robust claims of detection of long-period planets, or substellar companions (e.g., Neuhäuser et al. 2005; Mohanty et al. 2007; Oppenheimer et al. 2008; Lafrenière et al. 2008; and in particular Kalas et al. 2008; Marois et al. 2008). These objects have projected separations of ∼20 to a few 100 AU. The result of Lafrenière et al. (2008, hereafter LJvK08) is interesting, since it claims detection of an 8+4−1 MJ object at apparent separation of 330 AU from the parent star (1RXS J160929.1−210524)—an approximately solar-mass star in the young (∼5 Myr) Upper Scorpius association. Taken at face value this represents one positive detection out of ∼80 potential planet hosting stars, or a detection rate of ∼1%—although clearly this is a very crude evaluation since the detectability of companions is a complex function varying from star to star.
Most protostellar/protoplanetary disks are approximately 100 AU in radius, and core accretion models predict giant planet formation close to the water–ice snow line (less than 10 AU), where coalescence timescales are close to, or within, the disk lifetime against photoevaporation (e.g., Hollenbach et al. 2000). Gravitational instability models (Boss 1997) might have a higher probability of forming massive planets at large radii—but this would require quite unusually massive and extended disks. As pointed out by LJvK08, this suggests that migration and/or scattering mechanisms might be necessary to explain a planet at projected separation of 330 AU. Another possibility is that highly mass-imbalanced wide binary systems could form, with stellar and substellar companions.
In this Letter, we present the results of numerical simulations of ensembles of 100 planetary systems that capture the effect of strong dynamical evolution through planet–planet scattering. We present a number of proof-of-principle predictions for objects with large semimajor orbital axes (greater than 100 AU), in particular those in young systems. Our methodology stems primarily from the work of Juric & Tremaine (2008), and other investigations of dynamical effects (e.g., Papaloizou & Terquem 2001; Adams & Laughlin 2003, 2006; Ida & Lin 2004; Veras & Armitage 2004; Debes & Sigurdsson 2006; Chatterjee et al. 2008; Ford & Rasio 2008; Raymond et al. 2008; Thommes et al. 2008b). Juric & Tremaine (2008) suggest that, for the purposes of studying dynamical evolution, planet formation may be simplified into two stages. The initial episode (Stage 1) involves the assembly and evolution of a circumstellar and proto-planetary disk over ∼106–107 years. During Stage 1, planetary embryos are formed, destroyed, and interact strongly with their disk surroundings. Stage 2 takes place once the major part of the gaseous protoplanetary disk has dissipated. It involves the dynamical evolution of any planets that remain, including collisions, mergers, scattering, and significant planet ejection from the systems over the following tens to hundreds of millions of years. Using large ensembles of statistically identical initial conditions for planetary systems, Juric & Tremaine (2008) demonstrated that, in quite general circumstances, aspects of the observed distribution of giant exoplanet orbital eccentricities can be reproduced. This result is not strongly dependent on initial conditions as long as systems are dynamically active, or partially active. It also demonstrates the utility of studying the statistical dynamical properties of a population of planetary systems, rather than trying to identify a physical process that prescribes orbital characteristics in any single system.
Planets in dynamically active systems experience strong mutual interactions and frequent encounters, and undergo relaxation. This typically occurs when the spacing relative to the mutual Hill radii (e.g., Henon & Petit 1986) is small. Partially active systems contain only a subset of objects undergoing strong and frequent interactions, and while inactive systems still evolve, there are typically few or no major encounters in this case. Both active and partially active systems tend to lose planets as a function of time (instability timescales are given by Chatterjee et al. (2008)). This occurs through collisions (and mergers), capture by the parent star, or very often through complete ejection from the system. Thus, the Stage 2 dynamical paradigm predicts both a large population of ejected planets and a systematic excess of planets in very young systems compared to older systems.
2. SIMULATIONS AND METHODOLOGY
We employ the publicly available N-body integration routines in the MERCURY6 package (Chambers 1999) to follow the dynamical evolution of statistically identical ensembles of planetary systems. MERCURY6 also allows for close encounters and fully inelastic (no fragmentation) star–planet or planet–planet collisional mergers. We have found that a high-accuracy Burlisch–Stoer integrator typically provides relative energy errors of less than 10−5 up to integration times of 107 years for the systems we are modeling. For simplicity, we therefore use the Burlisch–Stoer BS2 integrator in MERCURY6 and restrict our integrations to 107 years—commensurate with the ages of the young planetary systems we hope to place constraints upon. Our methodology is otherwise similar to that of Juric & Tremaine (2008). We note that for dynamically active systems all but the final ∼10% of planet ejection/loss events have occurred by 107 years. Our simulations do not include the effects of tidal dissipation between planetary objects and the parent star. Dissipation due to star–planet tides can act to dampen orbital eccentricities, and is an important physical process for planets with short periods. We choose to circumvent these issues by ignoring planets with semimajor axes of less than 0.1 AU in our final results. This is well motivated in that tidal damping is expected to operate on timescales ranging from ∼109 to ∼1012 years for objects with ∼20 day orbital periods around solar-mass stars (e.g., Matsumura et al. 2008, and references therein).
We adopt a (somewhat arbitrary) physical ejection radius of 2000 AU, but note that for semimajor axes beyond ∼1000 AU, secular effects from Galactic tides may become important (see, e.g., Higuchi et al. 2007). It is possible that these tides would cause some of the planets we label as ejected to remain on bound orbits at very large separations. We neglect this possibility altogether in the present study.
2.1. Initial Populations of Objects
We perform two sets of runs, each with statistical ensembles of 100 systems. We populate the systems with 10 planets and adopt initial properties such that the systems are dynamically active. Masses are drawn randomly from a uniform distribution in log(M) and are limited to a range of 0.1–10 MJ. Initial orbital inclinations (i) for all runs are drawn randomly from a Rayleigh (Schwarzschild) distribution—where the peak, or mode, is equal to the distribution width parameter σ, and we set σi = 0.3 degrees. Orbital eccentricities are also drawn from a Rayleigh distribution with σe = 0.3. In one run, that we denote as R30, the initial semimajor axes are drawn randomly from a uniform distribution in log(distance) between 0.1 AU and 30 AU. In the second run, denoted as R100, the outer truncation is increased to 100 AU. All systems contain a central stellar mass of 1 M☉.
The R30 and R100 runs illustrate as simply as possible the dependence on initial conditions of the dynamical relaxation process. The two leading scenarios for planet formation, core accretion and direct gravitational collapse, should lead to significantly different initial orbital outcomes for young planets (e.g., Lissauer 1993; Boss 1997; Papaloizou & Terquem 2006). Each of these scenarios faces challenges (e.g., Matzner & Levin 2005; Rafikov 2005; Dong et al. 2008) and they could in principle both operate in individual systems or in the protoplanetary population as a whole (e.g., Ribas & Miralda-Escude 2007). Given these complications, what is important in the present context is that the late dynamical architecture of exoplanetary systems could potentially allow us to discriminate between competing formation scenarios. Specifically, formation via core accretion is favored in the vicinity of the snow line in the protoplanetary disk, typically well within a few tens of AU, while gravitational collapse is favored at larger orbital distances. This is the primary motivation for the two types of runs, R30 and R100, performed here.
2.2. Monte Carlo "Snapshots"
In order to model apparent separations of planet–star systems from our simulation outputs, we must take into account the orbital characteristics and potential orientation of systems relative to a hypothetical observer. We adopt a simple algorithmic approach. First, the semimajor axis, eccentricity, and orbital period are used to construct a probability distribution for the instantaneous radial star–planet distance, assuming an observation made at a random time during the orbital period. We then assume that the orbital plane is randomly orientated with respect to the observer and project the star–planet radial separation onto the sky to obtain an apparent separation.
We can draw an arbitrary number of "snapshot" projections for all remaining bound planets in our simulations. Since this step is not computationally intensive, the actual number of snapshots is simply determined by the need to adequately populate the mass-projected separation plane. Here we generate 1000 such projections per individual planet, irrespective of any other bound planet in the system. For simplicity, we then generate a cumulative distribution of projected separations for all individual planets, ignoring entirely conditional probabilities associated with multiple planetary systems. This is sufficient for our present purpose, focused on the apparently most distant planet, but a more careful analysis would clearly require conditional probabilities for projected separations to be accounted for.
Uncertainties in the separation distributions are estimated from the dominant Poisson statistics of the original simulations with a total of only 10 × 100 statistically identical planets, as opposed to the more numerous projection realizations. For instance, if N uniquely tagged planets, out of the total 1000 original ones, contribute to the statistics beyond some fixed projected separation, an uncertainty 
is applied to the apparent separation distribution. Plotted error bars/envelopes are 1σ.
3. RESULTS
In Figure 1, we plot the mean number of planets per system as a function of time for the two runs R30 and R100. The more closely packed planets in R30 result in a faster decay, and by 107 years, when both runs are much more dynamically relaxed, the mean number of remaining planets in R30 is slightly lower than that in R100, with ∼2 planets per system in both cases. These results are in good agreement with those of Juric & Tremaine (2008), who demonstrated that integrating to 108 years (by employing a symplectic integrator) results in only moderate further evolution, owing to the loss already of an average 80% of the initial population by 107 years.
Figure 1. Mean number of remaining planets per system is plotted as a function of time for runs R30 and R100. Most collisions (planet–planet and planet–star) and ejections occur within the first 105–106 simulated years.
Download figure:
Standard image High-resolution imageBy 107 years in run R30, 25% of the initial population have collided with the central star, 12% undergo collision and merger with other planets, and 44% are ejected (the majority by scattering into an unbound e > 1 orbit, some by passing beyond 2000 AU, see below). For R100, 21% collide with the central star, 11% undergo collision and merger, and 44% are ejected. This indicates that the primary difference in final planet numbers at 107 years between R30 and R100 is almost entirely due to the higher rate of collision with the central star in R30.
In Figure 2, the eccentricity versus semimajor axis is plotted for all objects at several timeslices for the R100 run to illustrate the orbital evolution of the population. In the first ∼106 years, a significant number of objects are scattered onto high eccentricity and large semimajor axis orbits. The great majority (greater than 99.2%) of objects that reach eccentricities e > 0.9 become unbound (e > 1) due to scattering close to periastron.
Figure 2. Orbital eccentricities and semimajor axes of objects plotted at a number of times for the R100 run. Starting at the upper left and moving right, plots show the distribution at 0, 103, and 104 years. Lower panels, left to right, are for 105, 106, and 107 years. Of the highest eccentricity objects (e > 0.9) at any given time, less than ∼0.8% remain on bound orbits by 107 years.
Download figure:
Standard image High-resolution imageIn Figure 3, we plot the cumulative number of planets beyond a given apparent separation. In both R30 and R100 runs, some outward diffusion of the planet population is clearly seen by comparison to the initial distribution. Cumulative distributions after 106 and 107 years are shown, and only modest evolution is seen at these late times. In both runs, there is some reduction in the fraction of objects between a few AU and ∼100 AU, from 106 to 107 years.
Figure 3. Percentage of the remaining planet population with apparent separation greater than a given value is shown for the two models explored. Left panel: R100 ensemble results. Black curve and green error bars correspond to the population after 106 years, blue curve and gray error bars correspond to the population after 107 years. Errors are computed from the Poisson uncertainty in the underlying number of planets used to seed the Monte Carlo projected separation distributions. The red dashed curve shows the initial conditions. Right panel: R30 ensemble results. In both panels, the consequences of an outward diffusion of planets from their initial distribution are clearly seen. Planets in the tail of the projected separation distribution typically have eccentricities e≳ 0.4 (see Figure 2).
Download figure:
Standard image High-resolution image
For comparison with the putative planet detection of LJvK08, we find that the population of planets with 330 AU apparent separation or greater in the R100 run is 0.9(±0.6)% with a variation of only ∼0.2% between 106 and 107 years. In the R30 run, the corresponding numbers are 0.6(±0.3)% at 106 years and 0.2(±0.2)% at 107 years. We note that these results are broadly consistent with Veras et al. (2008) who find that ∼0.5% of dynamically evolved (scattered) systems at late times (greater than 107 years) could still harbor a long-term stable planet with e > 0.8 and a > 1000 AU.
For a minimum apparent separation of 100 AU, at 106 years, these fractions increase to 7.1(±2.5)% in the R100 runs and 2.3(±1.0)% in the R30 runs. Later, at 107 years, there is little evolution in the R100 results, but a drop to 1.4(±0.6)% in the R30 results.
Although these numbers are not suitable for a direct comparison to current imaging survey constraints, we nonetheless note that the constraint of Nielsen et al. (2008) that less than 20% of systems should harbor >4 MJ planets between 20 and 100 AU is broadly consistent with the results of our simulations. The R30 runs indicate that less than ∼16% of planets of all masses would be expected at larger than 20 AU apparent separations (at 107 years). The corresponding R100 result is less than ∼30%.
The results shown in Figure 3 do not discriminate between planets of different masses as a function of projected separation. This is a particularly important issue for direct imaging surveys, with yields depending on a combination of ages, masses, and projected separations. We believe that a larger number of statistical realizations than performed here is needed to address this point reliably and postpone a more detailed investigation to future work.
4. DISCUSSION
Dynamical relaxation of young planetary systems can produce a population of long-period planets within 106–107 years of planet formation.
The results presented here are subject to a number of limitations. Although the outcome of dynamical relaxation of active systems in terms of eccentricities is not strongly dependent on the initial conditions set during Stage 1 planet formation (Juric & Tremaine 2008), the precise relationship between other population characteristics at ∼107 years in Stage 2 and initial parameters needs further exploration. For example, the differences between our R30 and R100 runs indicate that information on semimajor axis distributions is indeed retained in the population of objects at later times. The population of long-period planets could therefore potentially provide important constraints on the outcome of Stage 1 planet formation.
Furthermore, to make meaningful comparisons with observational surveys, one must allow for at least two biases: (1) not all stars surveyed will necessarily form a population of giant planets during Stage 1. (2) Not all systems with giant planets will be dynamically active at the start of Stage 2. For instance, Juric & Tremaine (2008) estimate on the basis of the observed distribution of exoplanet eccentricities that perhaps 75% of currently known planetary systems have been dynamically active. Less dynamically active systems, which were not directly addressed here, may provide better records of the planet formation outcome but at the expense of a smaller range of projected separations (with respect to the dynamically relaxed apparent separations shown in Figure 3).
The detailed strategy of an imaging survey is also important for making meaningful comparisons with theoretical expectations. If the planet formation phase (Stage 1) lasts between 106 and 107 years, the timescale for disk dispersal/evaporation (e.g., Hollenbach et al. 2000), then a "blind" survey of young stars—where the target objects span a wide range of ages up to a few 107 years—may not be well suited to probe the earliest Stage 2 phase. Indeed, the fast (Stage 2) relaxation phase (Figure 1), relative to a longer duration Stage 1 gaseous phase, suggests that it may be statistically difficult to identify the systematically larger number of planets expected in younger, dynamically active protoplanetary systems. On the other hand, Thommes et al. (2008a) present evolutionary scenarios based on mean-motion resonance locking and remnant planetesimal disk evolution in which the transition from Stage 1 to the Stage 2 relaxation occurs much more gradually. Finding a larger number of planets in younger systems would provide strong observational evidence in support of the dynamical relaxation scenario.
We conclude that the detection of massive planets at apparent separations of a few 100 AU need not be a challenge to the core accretion model of planet formation if dynamical relaxation occurs through extensive planet–planet scattering. Furthermore, the number of very long period planets resulting from scattering during the first few Myr following disk dispersal may be consistent with current upper limits, and claimed detections, from all existing planet imaging campaigns.
These are encouraging but very preliminary results. They suggest that a more systematic exploration of dynamical relaxation outcomes for long-period planets, as a function of initial conditions and system age, could provide a valuable additional means to interpret the existing and future yields of direct imaging surveys. Furthermore, since microlensing surveys probe the same projected planetary separations, but typically at much later system ages, they provide very useful complementary information on the late stages of the dynamical relaxation process. Together such data would allow consistency checks to be performed on the systematic trends expected between early and late stages of orbital evolution in this promising model.
The referee is thanked for comments that helped improve the Letter. C.S. acknowledges the support of Columbia University's Initiatives in Science and Engineering, and NASA Astrobiology: Exobiology grant NNG05GO79G.