A&A 442, 359-364 (2005)
F. Marzari1 - H. Scholl2 - P. Tricarico3
1 - Dipartimento di Fisica, University of Padova, Via Marzolo 8, 35131 Padova, Italy
2 - Observatoire de la Côte d'Azur, BP 4229, 06304 Nice Cedex 4, France
3 - Department of Physics, Washington State University, PO Box 642814, Pullman, WA, 99164-2814, USA
Received 31 March 2005 / Accepted 3 July 2005
We investigate the dynamical stability of the 3/1 resonant planets b and c in the extrasolar planetary system around 55 Cancri by applying Laskar's frequency map analysis. We find that the region with low diffusion speed extends to high eccentricity for both planets and that the observed system is deeply embedded in this region. The dynamics while in resonance and the influence of planet e on the resonant couple are also analysed. Using a simple model for the capture of the planets in resonance during migration we show that evolutionary tracks from smaller to the presently observed eccentricities lie in the stable region.
Key words: celestial mechanics - planetary systems - methods: N-body simulations
The extrasolar planetary system around the main-sequence star 55 Cancri (= Cancri) is presently the most crowded known system with four planets. Three planets reside in a comparatively small region between 0.038 and 0.24 AU. The fourth planet is further out at 5.26 AU. The masses of the four planets range from Neptune-like to Jupiter-like masses and their orbital eccentricities range from almost zero to 0.44 (McArthur et al. 2004). The central star is a G8V main sequence star with a high metallicity ([Fe/H] = 0.27, McArthur et al. 2004) and an estimated age of about 5 Gyr. The inner planet, 55 Cnc e, is the least massive known extra-solar planet orbiting a Sun-like star with a mass comparable to that of Neptune. This planet at 0.038 AU is not yet fully confirmed and doubts about its location have been raised by Wisdom (2005) very recently, who suggets that this planet has a much larger semimajor axis of nearly 0.8 AU and a mass of about 1.8 Neptune masses.
The large eccentricities of the planetary orbits, their vicinity to the central star and the age of the central star suggest that the planetary system has dynamically evolved. It cannot be excluded that the orbits of the innermost planets are still evolving due to tidal interactions with the central star. The second and third planet, termed planets b and c, respectively, are known to be locked in a 3/1 mean motion resonance (Zhou et al. 2004; Ji et al. 2003). Their respective semimajor axes are 0.115 and 0.240 AU corresponding to orbital periods of about 14.67 and 43.93 days. The eccentricity of planet b oscillates between almost zero and 0.22 while it oscillates between 0.25 and 0.44 for planet c. Possible close encounters between the two planets due to these eccentricity oscillations are avoided by the 3/1 resonance locking which is essential for the stability of the system. When moving one of the planets slightly out of the resonance, while keeping their observed eccentricities, close approaches would occur within timescales of 104 years according to our results. A "Jumping Jupiter'' phase (Weidenschilling & Marzari 1996; Rasio & Ford 1996) would start: The planets have their orbits frequently altered by close encounters until a planet is ejected from the system. Since the inner planets are so close to the central star, a planet may also fall into the star before ejection. However, the system may be stable for lower eccentricities, in particular for planet c. We found stability for eccentricities smaller than 0.25.
It is, therefore, plausible to assume that at least one of the planets, either b or c, migrated and entered the 3/1 resonance. It is well known that two bodies on converging orbits can be captured in a mean motion resonance (e.g. Murray & Dermott 1999). The masses of the two planets, their orbital eccentricities and their migration velocities are the three main parameters for capture. The 2/1 and 3/1 resonances are known to be the most efficient ones for capture (Murray & Dermott 1999). If, after capture, the dynamical system changes adiabatically, both planets remain locked in the resonance, increasing their orbital eccentricities during migration. The mass ratio between the two planets determines the outcome. Usually, the more massive planet right after resonance capture has a slower eccentricity growth. Planet c, the less massive planet, has a higher eccentricity, which supports the capture hypothesis. After an initial phase of eccentricity increase, eccentricities may evolve in different ways depending on the underlying migration scenario, on the mass ratio of the two resonant planets and on the behaviour of angular variables which are relevant for resonant motion. The eccentricity of one of the planets may decrease while the eccentricity of the other planet continues to increase. If the driving mechanism for migration does not disappear, one of the planets is removed either by ejection on a hyperbolic orbit, by falling in the central star or even due to a close encounter among the two planets (e.g. Moorhead & Adams 2005). If the origin of the resonant planets b and c is due to migration, it is plausible to assume that the driving mechanism disappeared leaving the two planets on their present orbits. Such a scenario was proposed to explain the resonant planetary system GJ 876 (Lee & Peale 2002; Kley et al. 2005).
Three major mechanisms are known to cause orbital planetary migration: planet-disk interaction (Goldreich & Tremaine 1980; Lin & Papaloizou 1986; Ward 1997; Tanaka et al. 2002; Kley 2003), tidal interaction between planet and central star (e.g. Rasio et al. 1996), and close encounters of a planet with planetesimals (Malhotra 1993; Murray et al. 1998). The first mechanism is the most obvious of the three to lead to resonance capture since it may result in an inward migration of a planet driven towards another planet as investigated, for instance, by Kley (2003). He showed for the 55 CnC system, by fully viscous hydrodynamical simulations and by N-body simulations, how an outside disk drives planet c towards planet b with a subsequent capture in the 3/1 resonance. He also modeled the capture of planets in a 2/1 resonance in the HD 82943 and GJ 876 systems. For the latter system, Lee & Peale (2002) investigated migration and resonance capture driven by disk-planet interaction using N-body simulations. Their model also includes possible eccentricity damping due to the action of the disk on the planets.
Resonance capture due to tidal interactions between planets and a central star is not so obvious. One would expect that the closer the planet is to the star the faster it migrates inwards, producing diverging orbits. If migration is driven by tidal forces, the two planets locked in resonance either originated very close to the resonance or were captured in resonance due to a different mechanism. The orbital evolution of migrating bodies in resonance driven by tidal forces in the frame of the Darwin-Mignard model was, for instance, investigated by Ferraz-Mello et al. (2003).
The third migration mechanism based on planetesimal scattering may result in converging planetary orbits. However, since in this case migration is driven by stochastic kicks of the semimajor axes, planets may be easily moved out of resonance locking.
Two conditions must be fulfilled for resonance capture of planets b and c independently of the underlying migration mechanism: Firstly, before capture, the orbital eccentricities of both planets b and c must be smaller than the present ones and, secondly, their orbits must converge. Since only after capture do their eccentricities increase to reach their present values, the dynamical system must evolve adiabatically. Therefore, it is important to know for this scenario whether or not there is a path in phase space from small to large planetary eccentricities embedded in a stable region.
In order to test the hypothesis that the present system is a result of migration and resonance capture, it is necessary to demonstrate the long-term stability of the present system and the stability of the system during migration while locked in the 3/1 resonance. In the present paper we concentrate on the long term stability and show for a few migrating systems that their evolutionary tracks reside in the most stable region. A further more detailed paper is devoted to capture and migration.
We apply Laskar's Frequency Map Analysis (hereinafter FMA) in the framework of a purely gravitational model which allows a fast exploration of the stability of a large number of orbits with very different starting values. This method determines the main frequencies of the system and computes their diffusion rates. Orbits with low diffusion rates are the most stable ones while those with fast diffusion rates are chaotic. The major advantage of the FMA is the comparatively short timespan of numerical orbital integration required to determine the stability properties of the system. Moreover, by computing the basic frequencies and amplitudes of the angular elements of the planets, the FMA gives additional information about the free eccentricities of the planets as well as about libration and circulation periods of relevant angular variables. This allows a rich statistical exploration of phase space without an excessive computational effort.
A total of about 10 000 systems with different values for orbital parameters of planets b and c are integrated numerically in the framework of a 5-Body problem including the central star and the four known planets, and a 4-Body problem in which planet e is not considered. Using the symplectic integrator SYMBA (Duncan et al. 1998) we cover a period of 105 yr. This time interval is long enough to measure with the FMA the most important secular frequencies of both planets b and c. A short timestep of 0.05 days is adopted in the numerical integration in order to account for the short orbital periods of the planets and for their high eccentricities. The initial semimajor axes, eccentricities and orbital angles of planet b and c are randomly sampled around their nominal values as given in Table 3 of McArthur et al. (2004). The authors provide different planetary masses, derived from radial velocity solely and also from astrometry. We have tried both mass sets without finding any significant difference in our main results concerning the size and shape of the most stable regions. All the orbits of the planets are assumed to be coplanar.
We retain only those simulations where at least one of the following two classical resonant arguments (e.g. Zhou et al. 2004) , librates over the whole timespan. The third resonant argument is not independent as it is a linear combination of the two other arguments. /2 (Murray & Dermott 1999).
The FMA technique (Laskar et al. 1992; Laskar 1993a,b; Sidlichovský & Nesvorný 1997) is based on a refined Fourier analysis of the angular variables of a dynamical system. The variation with time of the fundamental frequencies gives a measure of the diffusion of the trajectories in the phase space. It provides a complete dynamical map of the stability properties of a system and it also illustrates the global dynamics of the system. The FMA has been previously applied to study dynamical stability of minor bodies in the solar system (Nesvorný & Ferraz-Mello 1997; Melita & Brunini 2001; Marzari et al. 2002; Marzari et al. 2003a,b) and extrasolar planetary systems (for - Andromedae see Robutel & Laskar 2001).
We analyse the following
An average libration amplitude D and the related center of libration of the nominal critical argument of the 3/1 resonance is estimated as the mean of the maximum libration amplitude over short sub-windows ( yr) of the whole integration timespan. Similarly, we compute the libration amplitude of the apsidal resonance over a sub-window of yr. The choice of the sub-windows is made after a detailed analysis of the libration period of the two resonances in our simulations.
The diffusion speed of the resonant system in phase space is measured as the negative logarithm of the standard deviation of on running windows of 5000 years over the entire integration timespan. The circulation or libration period of is less than 1000 years. Thus, the running window is long enough for a precise computation of the frequency. Why do we measure the diffusion speed of the signal and not the one of the more conventional signals or ? The reason is that both and might oscillate on a short timescale with small amplitudes for those systems where the planets are not only in the 3/1 mean motion resonance but also in apsidal resonance. As a consequence, both and may not be good indicators of the diffusion speed of the dynamical system since their variation may be related to the resonance and not to chaotic diffusion.
After the FMA analysis, each orbit is characterized by the values of , , D, , , , .
|Figure 1: Diffusion maps showing the stability properties of the resonance in the plane. In the upper plot planet e is included in the system while in the lower plot only planets b,c and a are considered. Different gray levels represent values of ranging from 1 to 3. Filled circles mark the systems that are in apsidal libration. The empty circle represents the nominal planetary system.|
|Open with DEXTER|
|Figure 2: Diffusion maps in the plane. As in Fig. 1, the upper plot is for dynamical systems with planet e included while the lower plot does not include planet e. Symbols are as in Fig. 1.|
|Open with DEXTER|
The FMA measures the diffusion rate of a resonant system in phase space. The corresponding variables are the free eccentricities of the two planets and the libration amplitude of either or as outlined above. In Fig. 1 we illustrate diffusion portraits of the resonance in the space [ , ] for systems with (upper plot) and without planet e (lower plot). Indices b and c refer to planets b and c, respectively. Different gray levels represent values of the diffusion speed ranging from 1 to 3 (see the scale to the right of the figure). Light gray shading corresponds to large values of (>3) and to a low diffusion rate. It implies high stability and longest dynamical lifetime. The dark regions have small values of (around 1), a fast diffusion rate, and the orbits in these regions are highly chaotic. The nominal system is marked by a circle. It is situated deeply in the most stable region that extends in particular to the lower left corner which corresponds to small eccentricities of both planets. For free eccentricities of almost zero for both planets, the stable region is very small and the allowed range for resonant motion is and . This restricts the free eccentricities of both planets at capture. Moreover, the two planets must have comparable eccentricities, otherwise resonant motion is not possible. The path towards larger eccentricities is definitively quite narrow.
The influence of planet e on the stability of the resonance can be seen by comparing the upper map in Fig. 1 with the lower map obtained for a dynamical system including only planets b, c, and a. In the case with all four planets, the stability region has approximately a triangular shape. The resonance locking appears to be more sensitive to the values of free eccentricity of planet b. In fact, systems with values of larger than 0.35 are chaotic and become quickly unstable. Any value of the outer planet's free eccentricity appears to be permitted up to 0.6. This cut for large eccentricities of planet b is due to the presence of planet e. In the lower plot of Fig. 1 the stable region extends even to high values of . This is not unexpected since planet e is rather close to planet b and it has a maximum eccentricity of about 0.22. Once out of the resonance, large eccentricities lead rapidly to close encounters between planets b and c.
The small dots in the diffusion maps mark those systems where apsidal libration occurs simultaneously with the 3/1 mean motion resonance. Apsidal libration is not a necessary condition for stability. Apsidal libration occurs preferentially at smaller free eccentricities of planet b. This means that at capture, the apsides of the two planets are almost aligned. After capture, and after the initial phase of eccentricity increase, the original libration may switch to circulation without destabilizing the 3/1 resonance, or it may be preserved.
Figure 2 illustrates the diffusion maps in the (D, ) space, where D is the amplitude of the resonant argument . It shows clearly two stable regions, one with libration amplitudes lower than , and a smaller one roughly encompassed between . It is interesting to note that for very small eccentricities near zero of planet c, only librators with large amplitudes can occur. When planet e is taken away from the system (Fig. 2, lower plot) the stable region is more extended and additional non-chaotic orbits appear at large libration amplitudes and large eccentricities of planet b. The sharp cut in the libration amplitude is not an artifact of our numerical scheme, but we do not have a simple explanation for it.
The application of the FMA method yields important clues on dynamical features of the 3/1 planetary resonance, in particular of the resonance arguments . As shown in Fig. 3, the libration center of the resonance strongly depends on the libration amplitude D. For large libration amplitudes there is a single libration center located at . As D decreases, a bifurcation point is met at and the libration center splits in two branches approaching and , respectively. In Fig. 4 we illustrate two extreme cases of this behaviour. In the upper plot the critical argument librates around with an approximate libration amplitude of . The bottom plot shows the behaviour before the bifurcation point where the libration occurs around . Note that in the top case the libration period is well defined and it is close to 140 yr. In the bottom case the critical angle evolution is determined by two frequencies, one very short with a period of about 4 yr and a second one around 160 yr. The longer period is very close to the circulation (or libration in case of apsidal resonance) period of .
|Figure 3: Center of libration as a function of the corresponding libration amplitude D of the critical angle of the 3/1 resonance.|
|Open with DEXTER|
|Figure 4: Evolution with time of the critical angle of the 3/1 resonance for two different cases in our sample.|
|Open with DEXTER|
Combining Figs. 3 and 2, we can predict the behaviour of the libration argument when the planets are captured at small eccentricities in the resonance. In the beginning, the resonance argument librates around with very large amplitude. While eccentricities increase and the bifurcation point in Fig. 3 is reached, the libration center may switch to the upper or lower part in Fig. 3 with smaller libration amplitudes. With increasing eccentricities, the libration amplitude decreases.
|Figure 5: Migration tracks for planets b and c prior (dashed line) and after (solid line) capture in resonance. The tracks are superimposed on the diffusion map of Figs. 1 and 2, lower plots.|
|Open with DEXTER|
In Fig. 5 we superimpose three typical paths of migrating planets onto the diffusion maps of Figs. 1 and 2, upper plot. The region around the nominal system is expanded. There are no significant differences in the two diffusion maps of Fig. 1 for the values of and considered, so we took the upper plot where planet e is included. For the diffusion map in the (D, ) space, the case where the inner planet is included (Fig. 2 upper plot) has a less extended stability area. The choice of this map can be considered an extreme case.
In these examples of planetary migration, the inward drift of planet c is artificially induced by imposing a drag force , where denotes the mass of planet c, its velocity vector. is the timescale over which the semimajor axis of planet c decays. No drag force is applied to planet b. This drag force was used by Chiang (2003) to model the action of a viscous disk on a planet. Migration starts outside the resonance with eccentricities lower than 0.1. The migration tracks in the figure have timescales of the order of 1 Myr. Dashed lines correspond to the period before capture. Once trapped (solid line), eccentricities are pumped up. Depending on the initial conditions outside the resonance, different evolutionary tracks are obtained. For each planet, we compute approximate values of free eccentricities over short time intervals.
All migration tracks lie within the stable region that supports the scenario where the present eccentricities are due to migration while locked in resonance. After capture, the resonance critical argument librates around with a large amplitude. While the eccentricity increases, initially the libration amplitude is reduced until the bifurcation point shown in Fig. 3 is reached. As the amplitude continues to decrease, the libration center falls in either of the two branches illustrated in Fig. 3. The three migration tracks end up very close to the nominal system (the circle in Fig. 5) both in free eccentricities and libration amplitude. All the systems after the migration phase lie on the upper branch of Fig. 3 as the nominal system.
As outlined in the introduction, eccentricities of migrating planets may, after a first increasing phase, decrease. For the particular mass ratio of planets b and c and the underlying type II migration scenario, the eccentricity of planet c increases reaching at least a value of 0.5. The eccentricity of planet b may, after the initial growth phase, decrease very slightly. A small change in the initial conditions may lead a migrating system into a different branch of libration.
A more detailed analysis of trapping with different underlying migration scenarios will be published in a forthcoming paper.
Our results concerning the dynamics of the 3/1 resonance in the 55 CnC planetary system can be summarized as follows: