Letter to the Editor
Strong chaos induced by close encounters with Ceres and Vesta
^{1}
ASD, IMCCECNRS UMR 8028, Observatoire de Paris, UPMC, 77 Av. DenfertRochereau, 75014 Paris, France
email laskar@imcce.fr
^{2}
Observatoire de Besançon – CNRS UMR 6213, 41bis Av. de l’Observatoire, 25000 Besançon, France
Received: 17 June 2011
Accepted: 7 July 2011
We consider the full Solar System including(1) Ceres and some of the main asteroids, (2) Pallas, (4) Vesta, (7) Iris, and (324) Bamberga. We show that close encounters amongthese small bodies induce strong chaotic behavior in their orbits and in those of many asteroids that are much more chaotic than previously thought. Even if space missions will allow very precise measurements of the positions of Ceres and Vesta, their motion will be unpredictable over 400 kyr. As a result, it will never be possible to recover the precise evolution of the Earth’s eccentricity beyond 60 Myr. Ceres and Vesta thus appear to be the main limiting factors for any precise reconstruction of the Earth orbit, which is fundamental for the astronomical calibration of geological timescales. Moreover, collisions of Ceres and Vesta are possible, with a collision probability of 0.2% per Gyr.
Key words: chaos / methods: numerical / celestial mechanics / planets and satellites: dynamical evolution and stability / minor planets, asteroids: individual: Ceres / minor planets, asteroids: individual: Vesta
© ESO, 2011
1. Introduction
(1) Ceres, (2) Pallas and (4) Vesta are the largest bodies in the asteroid belt, with respective masses of 4.758, 1.114, and 1.331 × 10^{10} M_{⊙} (Fienga et al. 2011), while the mass of the smallest planet, Mercury, is about 1660 × 10^{10} M_{⊙}, and even Pluto is already 74 × 10^{10} M_{⊙}. These asteroids^{1} induce some perturbations of several km on the position of Mars (e.g., Kuchynka et al. 2010), and have long been taken into account in the elaboration of highprecision planetary ephemerides (Standish 1998; Fienga et al. 2008, 2009), valid over a few kyr. Meanwhile, they have only recently been introduced in the more simplified models used for the computation of longterm ephemerides, valid over several Myr, and aimed at paleoclimate reconstructions and geological timescale calibrations (Laskar et al. 2004, 2011). The time of validity of these computations is limited by the exponential divergence of the planetary orbits, resulting from their chaotic secular behavior (Laskar 1989, 1990). In practice, this chaotic behavior induces an increase in the initial position error d_{0} as d = d_{0} × 10^{T/10}, where the time T is expressed in Myr (Laskar 1989, 1990, 1999). In the La2004 solution (Laskar et al. 2004), the asteroids were not included individually, and only their averaged contribution was considered. The validity of the solution was estimated to be about 40 Myr. In La2010 (Laskar et al. 2011), as in INPOP08 (Fienga et al. 2009), the fivebodies (1) Ceres, (2) Pallas, (4) Vesta, (7) Iris, and (324) Bamberga were included and considered as the planets, with all gravitational interactions, while the numerical algorithm was being improved. In a surprising way, although the solution could be extended to about 50 Myr, its numerical accuracy was not improved, and it appears that the main reason for that is the presence of the asteroids that induce some additional instabilities. In this paper, we analyze the contribution of these asteroids and the instabilities that they generate.
2. Numerical simulations
In the numerical simulations, the full model comprises the eight planets of the Solar System, Pluto, and the Moon, as well asCeres, Pallas, Vesta, Iris, and Bamberga. The model and integrator that are used are the same as in the longterm solution La2010 (Laskar et al. 2011). To evaluate the effects induced by the asteroids, we performed two different sets of integrations. In a first set of solutions, S5_{a,b}, the asteroid interactions were included, while they were omitted in the second S0_{a,b}. The stability of the trajectories were analyzed by comparing two solutions with initial conditions that differ in position in only 10^{11} AU (~1.5 m) for the planets and 10^{10} AU (~15 m) for the asteroids, which is in both cases much less than the actual uncertainty of these initial conditions.
Fig. 1 Evolution of the Earth’s eccentricity uncertainty δe in S0_{a,b} (dotted line) and S5_{a,b} (solid line). In both cases, the initial difference at t = 0 is δe = 1.4 × 10^{12}. The thin solid line, with slope −1/10 represents the secular chaotic behavior (Laskar 1989, 1990). The eccentricity solution loses its validity when its uncertainty approaches its averaged value (δe ≈ 0.03). 

Open with DEXTER 
The divergence of the two sets of solutions are presented in Fig. 1. In both cases, the differences in the Earth’s eccentricity at t = 0 is 1.4 × 10^{12}, and thisuncertainty will first grow because of the difference in initial condition, but after 10 Myr, the differences in the S5_{a,b} solutions are about ten time greater than the differences in the S0_{a,b} solutions. As a result, the time of validity of the S5_{a,b} solutions is limited to about 60 Myr, while the (unrealistic) S0_{a,b} solutions are valid over nearly 70 Myr. We note that the secular chaotic regime, with an exponential increase in the error such as 10^{T/10} where T is in Myr (Laskar 1989, 1990), dominates after about 40 Myr (Fig. 1).
Fig. 2 Relative variations in thetotal orbital energy of the Moon and planets for S0_{a} (red), S0_{b} (green), S5_{a} (blue), S5_{b} (purple). The large variations in energy of S5_{a,b} correspond to asteroid close encounters (Table 1). 

Open with DEXTER 
The relative change in the total orbital energy of the Moon and planets for all solutions S0_{a,b} and S5_{a,b} is displayed in Fig. 2, after correction from the tidal dissipation of the Moon (Laskar et al. 2004). The energy of S0_{a,b} is preserved up to about 10^{13} in relative value, but some large and rapid variations in energy are observed in S5_{a,b}. These are the results of close encounters among the largest asteroids. In Fig. 2, the occurrences of these close encounters appear to be randomly distributed, and differ very much from S5_{a} to S5_{b}. This thus provides a noncompressible error in the energy of the planets, regardless of the accuracy of the initial conditions in the planetary orbits. The rapid variations in energy of S5_{a,b} that are larger than 2.25 × 10^{13} in relative value are easily identified as the results of specific asteroids close encounters (Table 1).
Close encounters for S5a (top) and S5b (bottom).
3. Stability of asteroids
The stability of the asteroids is estimated by measuring the rate of divergence inmean longitude of S0_{a,b} and S5_{a,b}. At t = 0, the differences in longitude in the solutions a,b for the five considered asteroids are on the order of d_{0} = 10^{10}. This error then grows by two means: first, the difference in initial conditions induces a difference in mean motion of about d_{0}, and thus a linear drift d_{0} t, while the chaotic divergence behaves as d_{0}exp(t/T), where T is the Lyapunov time and 1/T the maximum Lyapunov exponent. The maximum divergence will be given by the linear drift over the [0,T] time interval, and the exponential drift on [T,t]. The full rate of divergence is then (1)To determine the Lyapunov time, we solve for T in d(t_{1}) = 1 (Eq. (1)), where t_{1} is the first time for which d(t) = 1. The Lyapunov exponents (1/T) are given for S0_{a,b} and S5_{a,b} in Table 2. We note that the Lyapunov exponents LE_{0} computed in the absence of asteroid interactions are close to the values computed in a different way by Knežević & Milani (2003) and Novaković et al. (2011). We note, however, that in our integrations S0_{a,b}, Vesta is much more regular than Ceres, unlike Novaković et al. (2011). When the asteroids are taken into account (S5_{a,b}), the behavior of all asteroids is much more chaotic, and Ceres and Vesta, in particular, become strongly chaotic, with Lyapunov times of only 28 900 and 14 282 yr (Table 2, Fig. 3). We anticipate that numerous asteroids in the main belt will behave in the same way with, as is the case for Ceres and Vesta, much more chaotic behavior than previously thought when the collective effects of the asteroids were not taken into account.
Lyapunov exponents (LE) and Lyapunov times (LT).
Fig. 3 Decimal logarithm of the differences in longitude for the asteroids in the S0_{a,b} (dotted line) and S5_{a,b} (solid line) solutions. The initial conditions differ by about 1.4 × 10^{10} in longitude. 

Open with DEXTER 
Collision probability (in 10^{3} per Gyr) for each asteroid pair (i,j).
In the region of the Vesta family, we can also look more globally at the stability of the asteroids (Fig. 4). In this experiment, the initial conditions were taken on a regular grid in (a,e) with a 0.01 step size in e and 0.001 AU step size in a for a total of 51 × 701 = 35 751 orbits, for which all the other initial elliptical elements (i,ϖ,M,Ω) are set to zero. The orbits were computed over 100 kyr with a step size of 0.01 yr. Frequency analysis was then performed in order to obtain a global view of the stability of the trajectories (Laskar 1993; Dumas & Laskar 1993). The experiment was performed with (Fig. 4a) and without (Fig. 4b) the contribution of the asteroid interactions. The presence of the asteroid interactions makes the orbits much less regular, and in particular, most of the very stable blue zone present in Fig. 4b no longer exists in Fig. 4a.
Fig. 4 Global stability of asteroid motion with a) and without b) the mutual interactions of Ceres, Pallas, Vesta, Iris, and Bamberga. The initial conditions are taken on a regular grid of 701 × 51 in (a,e). The color index, obtained through frequency analysis (Laskar 1993), is the measure of the diffusion in mean motion (and thus in semimajor axis through Kepler’s law) over two consecutive time intervals. Blue denotes very regular regions and red highly chaotic ones. 

Open with DEXTER 
4. Collision probability
The asteroid orbits present numerous close encounters. We have analyzed all close encounters among asteroids in S5_{a},S5_{b}, and in an additional solution S5_{c} similar to S5_{b}, but with the opposite initial offset. S5_{a},S5_{b},S5_{c} were computed over 121, 117, and 52 Myr, respectively and provide a total of 25 779 (resp. 10 509, 10 359, 4911) asteroid encounters closer than 0.01 AU. Because of the chaotic nature of the solutions, all these close encounters are seen as independent events and are added in a single set of random variables, with each random variable corresponding to the distances of close approach Δ_{i,j} for a given asteroid pair (i,j). We then assume that the density probability of this random variable is linear, with zero probability for the exact collision (Δ_{i,j} = 0). This is the case for the Rice density (Rice 1945), which represents the density probability of the distance to the origin of a twodimensional Gaussian random variable, which is linear in the vicinity of the origin. A linear fit is then performed for each pair (i,j) (Fig. 5) in order to provide the probability of collision, that is the probability for Δ_{i,j} to be smaller than the sum of the radii of the two bodies (Table 3). This direct method differs somewhat from previous methods that consider large sets of asteroids (Farinella & Davis 1992; Yoshikawa & Nakamura 1994; Vedder 1996, 1998). The probabilities of collision (Table 3) are on the order of 0.1% per Gyr, the highest being the probability of collision of Ceres with Vesta, which reaches 0.2% per Gyr.
5. Conclusion
Fig. 5 Probability density of the close encounters for the pairs 1–2, 1–4, 1–7, 1–324. All approaches closer than 0.01 AU are collected in bins of 0.0001 AU. The number of hits in each bin is plotted versus the distance of closest approach Δ (in AU/1000). The solid line is the linear fit to this data, which provides, after normalization, the probability density function for the close approaches (Table 3). 

Open with DEXTER 
We show here that when we consider their collective effect, the motions of the main asteroids are much more chaotic than previously thought. This is particularly true for the motions of the largest asteroids, Ceres and Vesta, which can be mistakenly considered as regular when their mutual perturbations are not taken into account. More important, the strong interactions during close encounters also affect the planetary motions, and appear as the main limiting factor for establishing very longterm planetary solution for the Earth eccentricity, beyond 60 Myr, which
would be useful for paleoclimate studies. Moreover, this limit appears to be an absolute limitation, because the horizon of predictability for the motion of Ceres and Vesta is less than 500 kyr (Fig. 3). We also expect that most of the elements of the Vesta family will be perturbed by the main asteroids and their motion is largely chaotic as for the few that are considered here (Fig. 4). It should thus be investigated whether this can induce some additional diffusion in the elements of the Vesta family that has not been considered before (Carruba et al. 2007). Finally, we have obtained a reliable estimate for the collision probability of several asteroid couples, which amounts to 0.2% per Gyr for Ceres and Vesta.
Acknowledgments
This work was supported by GTSnext, ANRASTCM, INSUCNRS, PNPCNRS, and CS, Paris Observatory. J.L. thanks P. Michel for discussions.
References
 Carruba, V., Roig, F., Michtchenko, T. A., FerrazMello, S., & Nesvorný, D. 2007, A&A, 465, 315 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Drummond, J., & Christou, J. 2008, Icarus, 197, 480 [NASA ADS] [CrossRef] [Google Scholar]
 Dumas, H. S., & Laskar, J. 1993, Phys. Rev. Lett., 70, 2975 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Farinella, P., & Davis, D. R. 1992, Icarus, 97, 111 [NASA ADS] [CrossRef] [Google Scholar]
 Fienga, A., Manche, H., Laskar, J., & Gastineau, M. 2008, A&A, 477, 315 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Fienga, A., Laskar, J., Morley, T., et al. 2009, A&A, 507, 1675 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Fienga, A., Laskar, J., Kuchynka, P., et al. 2011, Celestial Mechanics, submitted [Google Scholar]
 Knežević, Z., & Milani, A. 2003, A&A, 403, 1165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kuchynka, P., Laskar, J., Fienga, A., & Manche, H. 2010, A&A, 514, 96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Laskar, J. 1989, Nature, 338, 237 [NASA ADS] [CrossRef] [Google Scholar]
 Laskar, J. 1990, Icarus, 88, 266 [NASA ADS] [CrossRef] [Google Scholar]
 Laskar, J. 1993, Phys. D, 67, 257 [NASA ADS] [CrossRef] [Google Scholar]
 Laskar, J. 1999, Phil. Trans. R. Soc. London Ser. A, 357, 1735 [NASA ADS] [CrossRef] [Google Scholar]
 Laskar, J., Robutel, P., Joutel, F., et al. 2004, A&A, 428, 261 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Laskar, J., Fienga, A., Gastineau, M., & Manche, H. 2011, A&A, in press [Google Scholar]
 Novaković, B., Knežević, Z., & Milani, A. 2011, Synthetic proper elements hamilton.dm.unipi.it/astdys/ [Google Scholar]
 Ostro, S. J., Magri, C., Benner, L. A. M., et al. 2010, Icarus, 207, 285 [NASA ADS] [CrossRef] [Google Scholar]
 Rice, S. 1945, Bell System Technical Journal, 24, 46 [Google Scholar]
 Standish, E. M. 1998, JPL IOM, 312.F,98 [Google Scholar]
 Vedder, J. D. 1996, Icarus, 123, 436 [NASA ADS] [CrossRef] [Google Scholar]
 Vedder, J. D. 1998, Icarus, 131, 283 [NASA ADS] [CrossRef] [Google Scholar]
 Yoshikawa, M., & Nakamura, T. 1994, Icarus, 108, 298 [NASA ADS] [CrossRef] [Google Scholar]
All Tables
All Figures
Fig. 1 Evolution of the Earth’s eccentricity uncertainty δe in S0_{a,b} (dotted line) and S5_{a,b} (solid line). In both cases, the initial difference at t = 0 is δe = 1.4 × 10^{12}. The thin solid line, with slope −1/10 represents the secular chaotic behavior (Laskar 1989, 1990). The eccentricity solution loses its validity when its uncertainty approaches its averaged value (δe ≈ 0.03). 

Open with DEXTER  
In the text 
Fig. 2 Relative variations in thetotal orbital energy of the Moon and planets for S0_{a} (red), S0_{b} (green), S5_{a} (blue), S5_{b} (purple). The large variations in energy of S5_{a,b} correspond to asteroid close encounters (Table 1). 

Open with DEXTER  
In the text 
Fig. 3 Decimal logarithm of the differences in longitude for the asteroids in the S0_{a,b} (dotted line) and S5_{a,b} (solid line) solutions. The initial conditions differ by about 1.4 × 10^{10} in longitude. 

Open with DEXTER  
In the text 
Fig. 4 Global stability of asteroid motion with a) and without b) the mutual interactions of Ceres, Pallas, Vesta, Iris, and Bamberga. The initial conditions are taken on a regular grid of 701 × 51 in (a,e). The color index, obtained through frequency analysis (Laskar 1993), is the measure of the diffusion in mean motion (and thus in semimajor axis through Kepler’s law) over two consecutive time intervals. Blue denotes very regular regions and red highly chaotic ones. 

Open with DEXTER  
In the text 
Fig. 5 Probability density of the close encounters for the pairs 1–2, 1–4, 1–7, 1–324. All approaches closer than 0.01 AU are collected in bins of 0.0001 AU. The number of hits in each bin is plotted versus the distance of closest approach Δ (in AU/1000). The solid line is the linear fit to this data, which provides, after normalization, the probability density function for the close approaches (Table 3). 

Open with DEXTER  
In the text 