About the dynamics of the evection resonance
J. Frouard^{1,2}  M. Fouchard^{1,2}  A. Vienne^{1,2}
1  Institut de Mécanique Céleste et de Calcul des
Ephémérides (IMCCE), Observatoire de Paris, UMR 8028, Avenue
DenfertRochereau, 75014 Paris, France
2 
Université de Lille 1, LALIMCCE, UMR 8028, 59000 Lille, France
Received 31 July 2009 / Accepted 10 November 2009
Abstract
Context. The evection resonance appears to be the outermost
region of stability for prograde satellite orbiting a planet, the
critical argument of the resonance indeed being found librating in
regions surrounded only by chaotic orbits. The dynamics of the
resonance itself is thus of great interest for the stability of
satellites, but its analysis by means of an analytical model is not
straightforward because of the high perturbations acting on the
dynamical region of interest.
Aims. It is thus important to show the results and the limits
inherent in analytical models. We use numerical methods to test the
validity of the models and analyze the dynamics of the resonance.
Methods. We use an analytical method based on a classical
averaged expansion of the disturbing function valid for all
eccentricities as well as numerical integrations of the motion and
surfaces of section.
Results. By comparing analytical and numerical methods, we show
that aspects of the true resonant dynamic can be represented by our
analytical model in a more accurate way than previous approximations,
and with the help of the surfaces of section we present the exact
location and dynamics of the resonance. We also show the additional
region of libration of the resonance that can be found much closer to
the planet due to its oblateness.
Key words: celestial mechanics  planets and satellites: general
1 Introduction
The problem of the stability of the orbits for real Solar System objects is a wellknown topic and has been thoroughly studied by using many different dynamical models. Concerning the satellites of planets, an important problem is to determine the orbital stability limit below which the motion is bound to the planet,and above which the escapes of satellites on heliocentric orbits eventually occur. In the framework of the restricted circular threebody problem, the wellknown Jacobi constant allows us to derive useful limits to bounded motion (see for example Murray & Dermott 1999), but for prograde orbits it can be shown that stable regions subsist outside these limits due to the socalled evection resonance. A detailed knowledge of the dynamics induced by this resonance is therefore interesting for objects orbiting far away from the planet, such as the irregular (outer) satellites of the giant planets.
In their description of the lunar theory, Brouwer and Clemence (1961) present the evection term in , being the longitude of pericentre of the satellite and the longitude of the perturbing body, which appears in the development of the solar disturbing function as the ``largest periodic correction in the Moon's longitude'' and the resonance have been of first importance to the dynamic of the moon in the past (Touma & Wisdom 1998). For distant orbits such those of the irregular satellites of the giant planets, this correction is crucial, as shown by recent analytic models (Yokoyama et al. 2003; Cuk & Burns 2004; Beaugé & Nesvorný 2007). The dynamics of the resonance itself was studied for the first time by Hénon in his numerical study of the restricted threebody problem (Hénon 1969, 1970) where he shows that in the Hill's approximation, this resonance appears as a bifurcation of a family of simple periodic orbits (named ``g'' in his papers) at a value of semimajor axis of where is the Hill's radius of the planet. However, this feature was not found in these studies for retrograde orbits. Hamilton & Krivov (1997) studied the orbital motions of distant satellites using a generalized Tisserand constant and found the basic features of evection: resonant prograde orbits are elongated toward the perturbing body and their model shows in polar coordinates where is the argument of the resonance, a characteristic eightshape centered on e=0. This specific dynamic appears at a value of semimajor axis of , and exhibits two elliptic points at and .
Retrograde orbits are shown to be elongated perpendicular to the perturber's direction, but the eightshape observed once again in this direction was thought to be an artefact of the method by the authors. The importance of the evection resonance on the stability of satellites was pointed out by Nesvorny et al. (2003), where it is shown that the effects of the perturber's tides on a satellite in the resonance can accumulate at each aphelion passage and induce a important drift of the orbits, potentially causing the escape of the body.
The dynamics of the resonance was studied with an analytic model using an expansion of the disturbing function for the first time by Yokoyama et al. (2008). The ``eightshape'' of the resonance for the planar problem was recovered by the authors, both for the prograde and retrograde case, and the apparition of the resonance was derived as and for the prograde and retrograde case respectively.
However, numerically, resonant orbits can be found closer to the planet than predicted by these analytical models and do not seem to follow exactly the eightshape found by previous authors. The aim of this paper is, first, to show and explain the limits associated with analytical models of the resonance, and secondly to resort to numerical methods, with the aim of localizing and studying the resonance itself. Since we are interested in the stability of the irregular satellites of Jupiter, the results of this paper will be applied in this context. We provide also an additional result concerning the dynamics of the resonance obtained when one takes into account the oblateness of the parent planet in the analytical model. Unlike the distant satellite orbits studied throughout this paper, we show that this modification only affects orbits very close to the planet.
In the following section, we construct and use an analytical model of the resonance and show the results that can be obtained. In Sect. 3, we resort to numerical studies (surfaces of section) to precisely localize the resonance and compare the results with those provided by the analytical model. In Sect. 4 we show how the oblateness of the parent planet affects the dynamics of the resonance. Finally, we present our conclusions and outline future work that should be performed in our last section.
2 Analytical model
To determine the limits of the use of an analytical model, we follow the method outlined by Yokoyama et al. (2008) and construct an analytical model of the resonance based on a development of the disturbing function in . In addition we use an alternative seminumerical method to verify the results.
2.1 Development of the disturbing function in Legendre polynomials
As in Yokoyama et al. (2008),
we use a development of the disturbing function using Legendre
polynomials. We consider in all the following calculations in this
paper the motion of a satellite orbiting around a planet and perturbed
by the Sun (indicated by subscript )
in a planetocentric system. The reference plane is the orbital plane
described by the Keplerian motion of the Sun around the planet. The
disturbing function
related to the Sun of mass
acting on the satellite's motion can be written using a Legendre polynomials development
using the Legendre polynomials
Here k is the Gaussian gravitational constant and S denotes the angle between the radius vectors of the satellite and the Sun in the planetocentric frame. The term in the Eq. (1) can be expanded using the cartesian coordinates of the bodies ( ), r and being the norm of the radius vectors
Using the elliptical elements defined by f the true anomaly, the argument of pericenter, the longitude of node, and I the inclination, one has
and similar expressions for , , and .
2.2 Comparison in the order of Legendre polynomials
Since the polynomials appear as powers of in the Eq. (1) and due to the (presumed) high distance in semimajor axis of the resonance from the planet, we must study the validity of the approximation's order before using any analytical result. To make a clear statement about the approximations made by using different orders of the polynomials with respect to the real motion, in Fig. 1 we represent stability maps for a Jovian satellite in the planar problem with the initial conditions , the initial semimajor axis and eccentricity being chosen in the ranges AU and respectively. We take as constants the semimajor axis of the Sun AU and the mass of Jupiter . For each orbit, the total integration time was divided in two consecutive samples and of the same length and the maximal eccentricity reached by the satellite for each sample was determined. The grey color code is given by , which is a stability criterion allowing the detection of chaotic and resonant orbits with respect to the regular ones (see Morbidelli 2002, for a discussion). White regions indicate initial orbits that have escaped the Hill's sphere of the planet and the grey scale indicates stable orbits in black and chaotic ones in grey levels (resonant orbits, despite their stability, are also highlighted in grey levels in the figures). Since we are interested in the stability of the outer satellites, these maps are of great importance and allow the easy location of the various resonances and chaotic zones. The integration was performed simulating the motion with a 2nd order Legendre polynomial solely (namely the quadrupole term), then adding the 3rd order Legendre polynomial (the octupole term) and finally the full equations of motion.
Figure 1: Stability maps for a prograde Jovian satellite for years using a 2nd order Legendre development ( top), 3rd order Legendre development ( middle), and the full equations of motion ( bottom). See text for comments. 

Open with DEXTER 
By observing the differences, at the limit of stability, between the orbits in the three maps, we conclude from Fig. 1 that the 2nd order Legendre development is unable to correctly reproduce the dynamics at the limit of stability given by the full equations of motion. We note that in these maps, the evection resonance is represented by the stable zone detached from the large one, surrounded by escaping regions and ranging from a=0.145 AU to a=0.17 AU. On the other hand, the 3rd order development appears to be sufficient to approximate the real motion. We note that the 2nd order development can still be used as a starting basis for studying the longterm behaviour of real satellites (see Cuk & Burns 2004), since these objects are closer to the planet.
2.3 The planar restricted threebody circular problem
In the planar case, we automatically assume that thus . The Sun's orbit is chosen to be circular, thus and , implying that , being the mean longitude of the Sun. Since an analytical model based solely on the 2nd order Legendre polynomial is insufficient for our purpose, we extend the method of Yokoyama et al. (2008) to the third order.
2.3.1 Development and averaging
The reduction to the planar case, with the perturber on a circular orbit, of Eq. (2) infers that if we consider a prograde satellite orbit, implying that .
Considering first the 2nd order Legendre polynomial, the disturbing function is
(3) 
Since we are interested only in the critical argument of the evection resonance , we average R_{2} according to the mean anomaly M of the satellite : .
This averaging is performed using the exact closedform equations of the twobody problem : , and
By using these equations, we obtain an averaged expression of the
disturbing function that is valid for all values of the eccentricity of
the satellite, resulting in
(4) 
This is identical to Eq. (3.4) of Yokoyama et al. 2008.
If we consider a development in Legendre polynomials of up to order 3, the corresponding expansion is
(5) 
where it is necessary to expand and average the expressions and .
To this effect, we use the following averaging equations calculated with an algebraic manipulator :
The final expression for the averaged expression is
=  (6)  
(7) 
To study the dynamics of the resonance, we must place the problem in an integrable form, which can be achieved by choosing appropriate canonical variables within the Hamiltonian formalism of the problem. This is the aim of the next section.
Figure 2: Analytical averaged model. Dynamics of the evection resonance for a prograde satellite using a model up to the 2nd order polynomial for a=0.19 AU ( top left) and a=0.2 AU ( top Right), and up to the 3rd order polynomial for a=0.19 AU ( bottom left) and a=0.2 AU ( bottom right). 

Open with DEXTER 
2.3.2 Canonical transformations of variables
We place the problem in a Hamiltonian form by using the Delaunay variables
.
The Delaunay variables are defined in the planar case to be
The Hamiltonian corresponding to a disturbing function of up to 2nd order is
(8) 
where .
Our Hamiltonian in these variables is given by
The above Hamiltonian has two degrees of freedom . Indeed, since it is independant of l the corresponding conjugate momenta is constant. The additive term can thus be dropped from Eq. (9).
We then look for a canonical change of variables allowing us to reduce the Hamiltonian to one dimension, which is:
with the generating function , which implies that
In these variables, the momenta P_{2} is constant and the term can be suppressed from the Hamiltonian. The Hamiltonian in the resonant variables thus has one degree of freedom ( ) and can be written by changing its sign to
(10) 
(see Yokoyama et al. 2008).
Using the same method, we found the Hamiltonian corresponding to the expansion of the disturbing function up to third order
Tests have shown that inclusion of the subsequent orders (4,5,...) in the expansion of the averaged disturbing function leads to negligible differences in the study of the resonance, confirming the results for the comparison of Sect. (2.2) which concerns the nonaveraged case. We note that the Hamiltonian defined in Eq. (11) can possibly gives values of librating other than those presented in the following : however their corresponding momenta P_{1} are found to be negative, an impossibility for physical values of a and e (since ), so we restrict our study to . In additional, P_{1}(e) (with constant a) is bijective only in this case and we can show the dynamics of the model using the coordinates ( ) in spite of the coordinates ( ). In the following, is denoted by for simplicity.
In Fig. 2, we show the dynamical portraits of the resonance obtained with the 2nd and the 3rd order models for the prograde case in the polar coordinates ( , with . Despite the presence of two islands at and for the two approximations, we can discern a clear difference between their sizes and their formations; the 2nd order model predicts that the two islands appear at a=0.1878 AU, while in the 3rd order model, one of the islands ( ) appears at a=0 AU and the other () appears at a=0.1976 AU.
For retrograde orbits, we follow the convention of Saha & Tremaine (1993) concerning the definition of retrograde elements in the expansion. In Fig. 3, we show the dynamical portraits of the resonance in the retrograde case. Following Yokoyama et al. (2008), using a 2nd order Legendre polynomial the resonance appears at AU and the libration islands move to and . The 3rd order term slightly modifies the dynamics; the hyperbolic point now has a nonnull eccentricity (e=0.0121 for the semimajor axis in Fig. 3) with and the libration centers move to increasing values of with increasing semimajor axis ( in Fig. 3).
Figure 3: Analytical averaged model. Dynamics of the evection resonance for a retrograde satellite using a 2nd order polynomial for a=0.26 AU ( top), a 3rd order polynomial for a=0.26 AU ( bottom). 

Open with DEXTER 
The global behaviour of the resonance in terms of semimajor axis and eccentricity predicted by the analytical model up to order 3 is shown in Fig. 4. In the upper figure (prograde case), the location of the stable elliptic points (light curves) and their corresponding separatrixes (bold curves) for the two islands is plotted. The island begins at the origin of the semimajor axis at e=0 and only its upper separatrix is plotted; the lower separatrix remains at e=0. The second island start at a=0.1976 AU. Its center of libration and upper separatrix have a similar behavior to that of the previous island, and its lower separatrix, which is the hyperbolic point of the resonance, decreases asymptotically to e=0 with increasing a. For the retrograde case (bottom figure), the upper separatrix and the elliptic point are indicated by the upper and middle curves respectively, and the hyperbolic point by the curve of very low eccentricity.
Figure 4: Analytical averaged model. Localisation of the evection resonance for a prograde ( top) and retrograde ( bottom) satellite using the model up to the 3rd order polynomial. Top: elliptic points (red curves), separatrixes and hyperbolic point (blue curves). Bottom: upper separatrix (upper curve), elliptic point (middle curve), hyperbolic point (bottom curve). 

Open with DEXTER 
2.4 Numerical averaging method
To verify and expand the results obtained above using the analytical
model, we use a numerical averaging procedure of the disturbing
function. The disturbing function in cartesian coordinates
,
where
represent the position vectors of the satellite and the Sun respectively, is expanded in elliptical elements
,
where E
denotes the eccentric anomaly of the satellite. The dependance on the
longitude of the satellite is numerically averaged following the
procedure of Moons (1994) (where the change of variables comes from the Kepler equation)
(12) 
The remaining angles are fixed to values corresponding to the evection resonance : (prograde case), and (retrograde case). We have shown that for retrograde orbits, the libration centers move slightly and is no longer , but the difference in the localisation of the exact resonance is marginal, for example the error in e is smaller than 0.01 for a=0.5 AU. This method is applicable to all eccentricities and is not limited by the ratio , in contrast to the classical analytic Legendre expansion. The frequency can then be obtained using the corresponding ``secular form'' of the Lagrange equation
(13) 
which is calculated with a numerical partial derivative scheme and compared to the mean motion of Jupiter to obtain the elliptic and hyperbolic points of the resonance. The results represented in Fig. 5 show very good agreement with the stable and hyperbolic points given by the 3rd order analytical model, and the curves are merged with the corresponding ones of Fig. 4.
Figure 5: Numerically averaged model. Localisation of the elliptic and hyperbolic points of the evection resonance for prograde orbits ( top) and retrograde ones ( bottom). 

Open with DEXTER 
Figure 6: Time evolution of the evection resonant angle ( top) and eccentricity ( bottom) for a satellite with initial elements a=0.155 AU, e=0.5. The initial angles are set to zero. 

Open with DEXTER 
2.5 Discussion
We note several things about these results. The discrepancy between the dynamics caused by the 2nd order and 3rd order terms from the expansion of Eq. (1) is that despite being inferior in the ratio by an order of magnitude compared to the 2nd order term, the 3rd order polynomial therefore plays a major role due to the effect of strict e in its formulation, comparing to the e^{2} in the 2nd order one. We note that the ``eightshape'' resonant dynamics was also found by Hamilton & Krivov (1997), but the construction of their Generalized Tisserand constant was achieved using solely a 2nd order Legendre polynomial to introduce the mean longitude of the Sun in their model. For the model used in this paper, we have shown that the following terms in the expansion (of order >3 in the ratio ) become negligible in the interval of semimajor axis studied and do not introduce noticeable changes, as shown by Fig. 5.
The main consequence of inclusion of the 3rd order term is the loss of symmetry between the two islands in Fig. 2. The islands now have distinct evolutions, implying different sizes of the islands at given values of the semimajor axis. We note that the island ( ) can now be found at all values of semimajor axis.
However, despite the precision provided by the 3rd order term, the analytical model still fails to describe the real dynamics, since one can find numerically resonant orbits at semimajor axis and eccentricity values inferior to those predicted by the model, as shown for example in Fig. 6. In Fig. 1, we however show that a nonaveraged expansion up to the 3rd order term seems to be sufficient to correctly describe the true motion.
The discrepancies observed from the numerical experiments come from the fact that we are using a normal form given by Eq. (11), but without having explicitly written the corresponding generating function (see for example FerrazMello 2007) which allows to determine the analytical relations between osculating and averaged variables. The initial osculating elements are thus considered to be equal to the averaged elements used in the analytical model. If this simplification works well in many cases where the perturbation is small, here at the limit of stability it is nonnegligible (see a discussion in Yokoyama et al. 2003). To consider this problem, one could use a perturbation method that explicitly permits the calculation of the generating function, such as Lie series method. However, this applies to use an explicit development of the disturbing function in Fourier series to keep things simple, and thus introduces an expansion in eccentricity (as opposed to the averaging method used above that is valid for all eccentricities), which would limits the applications of the model for low to moderate excentricities. In the following, we use an alternative and more direct way by using surfaces of section to analyze the resonance. However, if the high perturbation prevents us from correctly use the analytical model, one can choose to apply this method to a region where it is ``safe'' as the one of Sect. 4. In this section, we investigate the evection resonance in a region much closer to the planet where the solar perturbation is small, taking into account the oblateness of the planet, and using the above analytical model.
In addition, as shown in the next section resonant retrograde orbits cannot be found numerically despite their existence being expected by the model. Here the problem originates from the fact that the retrograde resonance is localized in a large chaotic zone corresponding to escaping orbits.
3 Surfaces of section
To carefully study the planar problem, we take advantage of the
twodimensionality of the Hamiltonian of the problem expressed in the
rotating frame evolving with the mean motion of Jupiter, and we use the
Poincaré surface of section method (PSS). The corresponding Hamiltonian
in the rotating frame (see for example Valtonen & Karttunen 2006) is given by
(14) 
where x,y,P_{x},P_{y} are, respectively, the coordinates and momenta of the particle in the rotating frame centered on the barycenter of the system, is the mass of the planet, is the mass of the Sun, and we have the distances
The PSS is then constructed from the motion of the particle in the cartesian coordinates ( ) by choosing a surface S (for example y=0) and plotting the remaining coordinates () each time the particle crosses the surface S with . The dynamics of different orbits can then be represented for the same value of the Jacobi constant C, which is the only integral of motion of the system,
(15) 
We emphasize that the PSS corresponds to the complete perturbed motion of the particle as opposed to an analytical model constructed from an approximation. In Fig. 7, we show two stability maps integrated up to 1000 periods of Jupiter, corresponding to prograde orbits with the initial conditions (pericenter in conjunction with the Sun), and , (pericenter in opposition with the Sun). The upper map in Fig. 7 is similar to the last one of Fig. 1 but here we use the MEGNO as an indicator of chaos (Cincotta et al. 2003) as we have observed its qualities (Frouard et al. 2008). The MEGNO is a fast and reliable chaos indicator which uses the time evolution of a tangent vector computed from the variational equations of the motion. Apart from the detection of chaotic orbits, for ordered motions this indicator is able to clearly distinguish resonant orbits from nonresonant ones. In the maps of Fig. 7, the white regions indicate satellites that have escaped the Hill's sphere of Jupiter while having a hyperbolic orbit with respect to an elliptic Keplerian one. We identified the limit of stability for initial circular orbits at a semimajor axis of a=0.17 AU which is consistent with the value found by Alvarellos (1996). Dark blue indicates stable orbits, while purple to yellow show different levels of chaoticicity for unstable orbits.
From the studies already cited (and our own numerical experiments), we know that the evection resonance in the prograde case has two libration islands, situated at and . This explain the choice of the initial conditions for the stability maps, which are chosen to place the satellite in each island at the beginning of the integration. In the maps, we indicate four isoJacobian curves that correspond to the PSS showed in Fig. 8; curves from left to right in the maps have respectively the values C=3.041, 3.0395, 3.0392 and 3.0388. We can then explain the dynamic shown in the stability maps with the surfaces of section.
Figure 7: Stability maps obtained with the MEGNO for a prograde Jovian satellite and 4 isoJacobian curves. Initial conditions: pericenter in conjunction with the Sun ( top), pericenter in opposition to the Sun ( bottom). 

Open with DEXTER 
Figure 8: From top to bottom: surfaces of section for C=3.041, C = 3.0395, C = 3.0392, and C = 3.0388. 

Open with DEXTER 
For the higher value of the Jacobi constant (C = 3.041), we are able to recognize at the center of the figure the orbit corresponding to the family of simple periodic orbits ``g''^{} (Hénon 1969,1970) and the libration is still not allowed. We determined its appearance at C = 3.0402552. For C = 3.0395 (the second figure), we can see that several loworder resonances appear along with chaos due to the overlapping between them and the ``spreading'' of their separatrixes. The periodic orbit g (now the evection resonance ) has moved to lower x. We can see a good correspondence with the stability maps of Fig. 7, where the corresponding Jacobian curves is in the weak chaotic region. For C = 3.0392, the hyperbolic point of the resonance was created and corresponds to the ``critical orbit'' ``g_{1}'' in Hénon (1970), and produced the second island centered on x=1.033, which corresponds to the evection argument (we found the appearance of the island at C = 3.0392409). In Hénon (1970), the island corresponding to the original orbit g was divided at its center by the critical orbit g_{1}, creating two libration islands of similar size whose centers are the elliptic orbits g'. This is not the case here, where the bifurcation has created an island of a very small size comparing to the other. This difference can be caused by the Hill's approximation not being used in this paper. For the same value of C we can observe that after the bifurcation, the chaos related to the hyperbolic point is still confined, but is about to be mixed with the already present one, giving a nearly generalized chaos on the surface of section. Obviously, we observe this phenomenon on the stability maps, where the isoJacobian curve now passes through stable resonant and strongly chaotic zones. Quasiperiodic orbits survive only at high eccentricities. The map also shows, in accordance with the surface of section, that the isoJacobian curve passes through the detached stability zone, which corresponds to the newly created island. Finally, for C = 3.0388 (last surface of section), a value slightly higher than the value corresponding to the Lagrangian point L_{1} ( C_{L1} = 3.0387559), thus just before the possible escape of the particle, the chaos is generalized.
3.1 Discussion
We mention that the two islands are very sensitive to the initial conditions that we choose for the PSS; the motion will be resonant in one of the islands if the resonant angle at the beginning of the integration is very close to one of the exact values of the resonance (e.g., or ). In fact we can see in Figs. 7 and 8 that the islands are surrounded only by chaotic regions.
Can we draw conclusions about the analytical model with the help of the surfaces of section presented above? Figs. 2 and 8 cannot be rigorously compared because each of the four PSS correspond to a single value of the Jacobi constant, which implies that the orbits integrated in one PSS have different initial elements a and e with or (the particle always starts in the SunJupiter line, with a perpendicular velocity, thus the initial and are the same for all orbits). Despite this and the fact that different coordinates are used in Figs. 2 and 8, the similarity between the PSS and the dynamic given by the analytical model up to the 3rd order is striking. On the other hand, the surfaces of section obtained with the Hill's approximation (see Hénon 1970) are similar to the dynamics predicted by the analytical model up to the 2nd order. The main difference caused by the PSS is that the island is created at a specific value of C; the central periodic orbit for (see the upper graph in Fig. 8) corresponds to the periodic orbit g.
The evection for retrograde orbits cannot be found by means of surfaces of section. The family of simple periodic orbits ``f'' in Hénon (1970) indeed experiences no bifurcation. If orbits can numerically be found librating around the centers and , it is always temporarily, and they alternate between circulation and/or libration around the other elliptic point after a certain amount of time.
Figure 9 indicates the centers of libration given by the PSS for several initial osculating variables. Since we analyze the complete problem, the centers of the resonance are found to depend of the initial osculating variables and the use of a more sophisticated perturbation method than in Sect. 2 would have brought a similar phase dependency with respect to the initial variables. In the two figures, the left and right curves correspond to the librations and respectively. In Fig. 10, we present the location of the two islands of the resonance as given by the PSS in averaged elements, along with their counterparts found by the numerically averaged model (Sect. 2.4). The averaged elements were obtained by averaging the semimajor axis and the eccentricity of the orbits found by the PSS in Fig. 9 over 10 000 years. In this figure, the two left curves correspond to the orbits found with the PSS, superimposed on those of Fig. 5. We can appreciate the interval between the curves found by the two methods, which illustrates the lack of precision of the analytical model and emphasizes that the high solar perturbation prevent us from using it. In the next section, we show that the analytical model can still be used if the solar perturbation is sufficiently weak, that is very close to the planet.
Figure 9: Location of the evection resonance with the PSS in initial osculating elements. Top: . Left curve: libration with , right curve: libration with . Bottom: . Left curve: libration with , right curve: libration with . 

Open with DEXTER 
Figure 10: Location of the evection resonance for prograde satellites with the PSS (in averaged elements) (the two left curves) and by the numerically averaged model (the two right curves). The island is indicated by the dot line, the island by the solid one. 

Open with DEXTER 
4 Evection resonance taking into account the oblateness of the planet
To study the effect of the resonance very close to the planet, we start from the analytical model described in Sect. 2, and take into account the correction due to the second harmonic in the development of the potential of the planet. Indeed, close to the planet, the frequency of the pericenter of the satellite is insufficiently high to be comparable to the mean motion of the Sun, and does not allows the appearance of the evection resonance. The oblateness of the planet has the wellknown effect of increasing the frequency of pericenter of a satellite. One can note that this effect is only local, and we have checked that no apparent change is present for the outer region discussed in the previous sections. We have also verified how the description of the solar perturbation could be improved by using different orders of approximation as in Sect. 2.2 and conclude that due to the small ratio ( ), the 2nd order polynomial alone is sufficient. To obtain the additional term that depends on the oblateness of the planet, we proceed in the following with a method similar to that used at the beginning of this paper.
4.1 Development and averaging
The gravitational potential of a planet formulated in spherical coordinates (
), where the axis of revolution of the body coincides with its rotation axis, is written, following Duriez (2002), as
=  
(16) 
where m is the mass of the planet, is its equatorial radius, and the J_{n} are the zonal harmonics. The P_{n} are the Legendre polynomials defined in Sect. 2.1. This simplification limits the application of the model to planets with a sufficiently low value of obliquity to be considered negligible; this is the case for Jupiter ( ), Mercury ( ), and Venus (retrograde rotation: ). However, for the last two planets all values of eccentricity of the orbits affected by the resonance lead to pericenters that are inside the physical radius of the planet. For Jupiter, only the higheccentricity part of the resonance induces this effect. The physical coefficients of Jupiter are km and J_{2} = 0.01475. We note that the resonance for the cases of Mars and Earth was studied analytically by Breiter (2000). The perturbation caused exclusively by the 2nd order Legendre polynomial, and involving only the J_{2} coefficient is given by
(17) 
Following classical calculations and using the averaging relation in a similar way to that in Sect. 2.3.1, we obtain the expression for U_{J2} averaged over the mean anomaly of the satellite :
4.2 The planar case
When i=0, Eq. (18) expressed in resonant variables becomes
(19) 
The phase portrait of the resonance in polar coordinates for prograde orbits of semimajor axis a=0.00515 AU can be seen in Fig. 11 (retrograde resonant orbits were not found). The effect of the oblateness is to change the libration centers of the resonance by . The resonance also appears to be narrower.
Figure 11: Dynamical portrait of the evection resonance for a Jovian satellite with semimajor axis a = 0.00515 AU, taking into account the Jupiter's J_{2}. 

Open with DEXTER 
A numerical example of a resonant orbit integrated with the full equations of motion is shown in Fig. 12.
Figure 12: Time evolution of the resonant angle ( top) and eccentricity ( bottom) of an orbit integrated with the full equations of motion with initial elements a = 0.00515 AU, e=0.75, , and . 

Open with DEXTER 
Figure 13: Location of the evection resonance and its width for Jupiter. The locations in semimajor axis of the Jovian main massive satellites are indicated by vertical lines. 

Open with DEXTER 
In Fig. 13, the exact resonance and its width are shown. The location of the averaged semimajor axis of the four main Jovian massive satellites (Io, Europe, Ganymede and Callisto)^{} are also indicated by vertical lines. We note that the pericenter of the satellite's orbit is found to be inside the physical radius of Jupiter when the upper separatrix of the resonance exceed a = 0.01176 AU for e = 0.959. For the exact resonance, the critical values are a = 0.01628 AU with e = 0.971.
5 Conclusions
We have studied the evection resonance using an analytical model extending the computations of Yokoyama et al. (2008). The extended model has different dynamics compared to the previous one, mainly in terms of creation of the islands of libration. We have show with an independent method that additional extensions of the model do not alter the predicted dynamics. However, since this analytical model still appears to have limited precision compared to direct numerical integrations due to the chosen averaging method, we resort to surfaces of section to compare them with the model, precisely localize the resonance and study its features. The comparison shows that the extended analytical model, while remaining inaccurate in the localisation of the resonance, reproduces the main features of its dynamics. We determined the appearance of the resonance at the averaged values and for the librations and respectively. In addition, we showed that the resonance can be found much closer to the planet due to its oblateness and we indicated its location in terms of semimajor axis and eccentricity. However, important questions about the resonance remain unanswered, for example the dependencies of the resonance on the inclination of orbits and the eccentricity of the perturbing body are of particular interest for escape/capture problems and dynamical studies of outer satellites, and are currently explored. For the same reason, the persistence of the resonance when the motion is subjected to other perturbations, such as massive satellites, planets, or nongravitational forces, is also investigated and will be presented in a forthcoming paper. AcknowledgementsThe authors would like to thank the anonymous referee for important suggestions and comments.
References
 Alvarellos, J. L. A. 1996, Master's thesis, San José State Univ. Beaugé, C., & Nesvorný, D. 2007, AJ, 133, 2537 [NASA ADS] [CrossRef] [Google Scholar]
 Breiter, S. 2000, Celest. Mech. & Dyn. Astr., 77, 201 [NASA ADS] [CrossRef] [Google Scholar]
 Brouwer, D., & Clemence, G. M. 1961, Methods of Celestial Mechanics (New York: Academic Press) [Google Scholar]
 Cincotta, P. M., Giordano, C. M., & Simó, C. 2003, Physica D, 182, 151 [NASA ADS] [CrossRef] [Google Scholar]
 Cuk, M., & Burns, J. A. 2004, AJ, 128, 2518 [NASA ADS] [CrossRef] [Google Scholar]
 Duriez, L. 2002, Cours de mécanique céleste classique, http://lal.univlille1.fr/mecanique_celeste.html [Google Scholar]
 FerrazMello, S. 2007, Canonical Perturbation Theories  Degenerate Systems and Resonance, Astrophysics and Space Science Library (Springer) [Google Scholar]
 Frouard, J., Fouchard, M., & Vienne, A. 2008, SF2A2008, 121 [Google Scholar]
 Hamilton, D. P., & Krivov, A. V. 1997, Icarus, 128, 241 [NASA ADS] [CrossRef] [Google Scholar]
 Hénon, M. 1969, A&A, 1, 223 [NASA ADS] [Google Scholar]
 Hénon, M. 1970, A&A, 9, 24 [NASA ADS] [Google Scholar]
 Moons, M. 1994, Celest. Mech. & Dyn. Astr., 60, 173 [NASA ADS] [CrossRef] [Google Scholar]
 Morbidelli, A. 2002, Modern celestial mechanics : aspects of solar system dynamics (London: Taylor & Francis) [Google Scholar]
 Murray, C. D., & Dermott, S. F. 1999, Sol. Sys. Dynamics (Cambridge: Cambridge Univ. Press) [Google Scholar]
 Nesvorný, D., Alvarellos, J. L. A., Dones, L., & Levison, H. F. 2003, AJ, 126, 398 [NASA ADS] [CrossRef] [Google Scholar]
 Saha, P., & Tremaine, S. 1993, Icarus, 106, 549 [NASA ADS] [CrossRef] [Google Scholar]
 Touma, J., & Wisdom, J. 1998, AJ, 115, 1653 [NASA ADS] [CrossRef] [Google Scholar]
 Valtonen, M., & Karttunen, H. 2006, The Threebody problem (Cambridge: Cambridge University Press) [Google Scholar]
 Yokoyama, T., Santos, M. T., Cardin, G., & Winter, O. C. 2003, A&A, 401, 763 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Yokoyama, T., Vieira Neto, E., Winter, O. C., Sanchez, D. M., & de Oliveira Brasil, P. I. 2008, Mathematical Problems in Engineering, doi: 10.1155/2008/251978 [Google Scholar]
Footnotes
 ...''^{}
 The periodic orbit named g is the resonance with the critical argument in the inertial frame.
 ...Callisto)^{}
 taken from http://ssd.jpl.nasa.gov/?sat_elem
All Figures
Figure 1: Stability maps for a prograde Jovian satellite for years using a 2nd order Legendre development ( top), 3rd order Legendre development ( middle), and the full equations of motion ( bottom). See text for comments. 

Open with DEXTER  
In the text 
Figure 2: Analytical averaged model. Dynamics of the evection resonance for a prograde satellite using a model up to the 2nd order polynomial for a=0.19 AU ( top left) and a=0.2 AU ( top Right), and up to the 3rd order polynomial for a=0.19 AU ( bottom left) and a=0.2 AU ( bottom right). 

Open with DEXTER  
In the text 
Figure 3: Analytical averaged model. Dynamics of the evection resonance for a retrograde satellite using a 2nd order polynomial for a=0.26 AU ( top), a 3rd order polynomial for a=0.26 AU ( bottom). 

Open with DEXTER  
In the text 
Figure 4: Analytical averaged model. Localisation of the evection resonance for a prograde ( top) and retrograde ( bottom) satellite using the model up to the 3rd order polynomial. Top: elliptic points (red curves), separatrixes and hyperbolic point (blue curves). Bottom: upper separatrix (upper curve), elliptic point (middle curve), hyperbolic point (bottom curve). 

Open with DEXTER  
In the text 
Figure 5: Numerically averaged model. Localisation of the elliptic and hyperbolic points of the evection resonance for prograde orbits ( top) and retrograde ones ( bottom). 

Open with DEXTER  
In the text 
Figure 6: Time evolution of the evection resonant angle ( top) and eccentricity ( bottom) for a satellite with initial elements a=0.155 AU, e=0.5. The initial angles are set to zero. 

Open with DEXTER  
In the text 
Figure 7: Stability maps obtained with the MEGNO for a prograde Jovian satellite and 4 isoJacobian curves. Initial conditions: pericenter in conjunction with the Sun ( top), pericenter in opposition to the Sun ( bottom). 

Open with DEXTER  
In the text 
Figure 8: From top to bottom: surfaces of section for C=3.041, C = 3.0395, C = 3.0392, and C = 3.0388. 

Open with DEXTER  
In the text 
Figure 9: Location of the evection resonance with the PSS in initial osculating elements. Top: . Left curve: libration with , right curve: libration with . Bottom: . Left curve: libration with , right curve: libration with . 

Open with DEXTER  
In the text 
Figure 10: Location of the evection resonance for prograde satellites with the PSS (in averaged elements) (the two left curves) and by the numerically averaged model (the two right curves). The island is indicated by the dot line, the island by the solid one. 

Open with DEXTER  
In the text 
Figure 11: Dynamical portrait of the evection resonance for a Jovian satellite with semimajor axis a = 0.00515 AU, taking into account the Jupiter's J_{2}. 

Open with DEXTER  
In the text 
Figure 12: Time evolution of the resonant angle ( top) and eccentricity ( bottom) of an orbit integrated with the full equations of motion with initial elements a = 0.00515 AU, e=0.75, , and . 

Open with DEXTER  
In the text 
Figure 13: Location of the evection resonance and its width for Jupiter. The locations in semimajor axis of the Jovian main massive satellites are indicated by vertical lines. 

Open with DEXTER  
In the text 
Copyright ESO 2010