Laplace equilibria

Satellite dynamics on the Laplace surface

Scott Tremaine11affiliation: School of Natural Sciences, Institute for Advanced Study, Einstein Drive, Princeton, NJ 08540, USA; tremaine@ias.edu , Jihad Touma22affiliation: Department of Physics, American University of Beirut, PO Box 11–0236, Riad El-Solh, Beirut 1107 2020, Lebanon; jihad.touma@gmail.com , and Fathi Namouni33affiliation: Université de Nice and CNRS, Observatoire de la Côte d’Azur, BP 4229, 06304 Nice, France; namouni@obs-nice.fr
Abstract

The orbital dynamics of most planetary satellites is governed by the quadrupole moment from the equatorial bulge of the host planet and the tidal field from the Sun. On the Laplace surface, the long-term orbital evolution driven by the combined effects of these forces is zero, so that orbits have a fixed orientation and shape. The “classical” Laplace surface is defined for circular orbits, and coincides with the planet’s equator at small planetocentric distances and with its orbital plane at large distances. A dissipative circumplanetary disk should settle to this surface, and hence satellites formed from such a disk are likely to orbit in or near the classical Laplace surface. This paper studies the properties of Laplace surfaces. Our principal results are: (i) if the planetary obliquity exceeds there is a range of semimajor axes in which the classical Laplace surface is unstable; (ii) at some obliquities and planetocentric distances there is a distinct Laplace surface consisting of nested eccentric orbits, which bifurcates from the classical Laplace surface at the point where instability sets in; (ii) there is also a “polar” Laplace surface perpendicular to the line of nodes of the planetary equator on the planetary orbit; (iv) for circular orbits, the polar Laplace surface is stable at small planetocentric distances and unstable at large distances; (v) at the onset of instability this polar Laplace surface bifurcates into two polar Laplace surfaces composed of nested eccentric orbits.

planets and satellites: formation – planets and satellites: general

1 Introduction

In his study of Jupiter’s satellites, Laplace (1805) recognized that the combined effects of the solar tide and the planet’s oblateness induced a “proper” inclination in satellite orbits with respect to Jupiter’s equator. He remarked that this proper inclination increases with the distance to the planet, and defines an orbital plane for circular orbits that lies between the orbital plane of the planet’s motion round the Sun and its equator plane. This plane is called the Laplace plane.

More generally, the Laplace plane is usually defined as the plane normal to the axis about which the pole of a satellite’s orbit precesses111Unfortunately, the term is sometimes also applied to the invariable plane, the plane perpendicular to the total angular momentum of an -body system.. The “Laplace surface” is the locus of all orbits that do not precess (i.e., the secular motion of the node and apse vanishes).

In the most common situation, we consider circular satellite orbits around an oblate planet with non-zero obliquity, traveling around the Sun. The Laplace surface is then determined by the competition between the interior quadrupole potential from the equatorial bulge and the external quadrupole potential from the Sun. Close to the planet, the “classical” Laplace surface nearly coincides with the planetary equator, while at large distances it nearly coincides with the planetary orbital plane. The transition between these two orientations occurs near the “Laplace radius,” defined below in equation (24).

The Laplace surface is important because it traces the shape expected for a thin gas disk or dissipative particulate ring surrounding the planet. Thus it is not surprising that many planetary satellites orbit close to the Laplace surface; those that do probably formed from a circumplanetary gas disk while those that do not were presumably either captured from heliocentric orbits or experienced unusual events in their past history.

The purpose of this paper is to study the properties of the Laplace surface, including the stability of is generating orbits and its generalization to eccentric orbits. Although the Laplace surface has been known and studied for over two centuries, we believe that many of the results we present are new.

2 Secular equations of motion

2.1 The Hamiltonian

The Kepler Hamiltonian that describes the motion of a test particle orbiting an isolated point mass is

(1)

here is the position vector measured from the center of the planet, , , and is the semimajor axis of the test particle. The constants of motion are the Hamiltonian (or semimajor axis ), and the angular momentum and eccentricity vectors

(2)

these are related to the eccentricity and semimajor axis of the orbit by and .

We now examine how the motion of the test particle is affected by additional forces from the equatorial bulge of the planet and the Sun. The quadrupole potential arising from an oblate planet is

(3)

where is the quadrupole gravitational harmonic, is the planetary radius, is a Legendre polynomial, and is the polar angle measured from the planet’s spin axis, which is oriented along the unit vector .

The quadrupole potential of the planet may be enhanced by inner satellites. If the planet hosts satellites with masses , , on circular orbits in the equatorial plane of the planet with semimajor axes , then at distances the gravitational potential due to the satellites can be accounted for by augmenting to , where

(4)

Values of and for the giant planets and Pluto are given in Table 1.

planet (AU) (km) obliquity
Jupiter 5.2029 71492 0.014696 0.045020 35.36 743.3
Saturn 9.5367 60330 0.016291 0.070561 48.40 1080.1
Uranus 19.189 26200 0.003343 0.018699 63.96 2675.1
Neptune 30.070 25225 0.00341 0.024069 93.20 4600.8
Pluto 39.482 1151 14.296 419.6 6935.8

Note. – The planet’s semimajor axis and radius are and . Obliquity is the angle between the planet’s spin and orbital axes, which in the notation of this paper is given by . The Laplace radius and Hill radius are defined by equations (24) and (6).

Table 1: Properties of the outer planets

Since the solar tide is assumed to be weak, we may estimate its effects by averaging over the solar orbital period. We assume that the planetary orbit has semimajor axis and eccentricity , and denote the normal to the orbit by the unit vector . The obliquity of the planet is then . Since we need keep only the quadrupole term in the averaged solar potential,

(5)

The solar tide also causes the spin axis of the planet to precess, through its torque on the equatorial bulge. We neglect this effect to keep the analysis simple, since the precession rate of the planetary spin due to solar tides is normally much smaller than the precession rate of the satellite orbit (see Goldreich 1965 and Boué & Laskar 2006 for treatments that include precession of the planetary spin).

The planet’s “radius of influence” or “Hill radius” is

(6)

The Hill radius is roughly the point at which ; beyond the Hill radius the gravitational force experienced by a satellite is dominated by the solar tide rather than the force from the planet and most orbits are not bound to the planet. The Hill radius also marks the location of the collinear Lagrange points of the planet (Murray & Dermott, 1999). Our use of the orbit-averaged solar potential requires that the satellite radius .

The total non-Keplerian potential due to the oblate planet, inner satellites, and Sun is .

Now average over the Keplerian orbit of the test particle, which has semimajor axis , eccentricity , and orientation specified by the unit vectors along the angular momentum vector, towards pericenter, and . We have

(7)
(8)

The averaged potential is , where

(9)

The averaged Hamiltonian is then

(10)

By virtue of the orbit averaging, this Hamiltonian is independent of the mean anomaly, so its conjugate momentum is a constant of motion, which in turn means that the semimajor axis may be treated as a constant.

Now set

(11)

The angular momentum and, as suggested by the notation, is the eccentricity vector (eq. 2). If we define a dimensionless potential , we have

(12)

an unimportant constant has been dropped.

To repeat, the Hamiltonian is based on the assumptions that (i) the precession rate of the planetary spin due to solar tides is negligible; (ii) the satellite is a massless test particle; (iii) the Sun is far enough from the planet that the solar tide can be approximated by a quadrupole; (iv) the satellite is far enough from the planet that the potential from the planet and the inner satellites can be approximated as a a monopole plus a quadrupole; (v) the perturbing forces due to are weak enough that the secular equations of motion can be used to describe its orbital evolution.

2.2 Equations of motion

Using equation (2) it is straightforward to show that the Poisson brackets of and are

(13)

where is the antisymmetric tensor.

The time evolution of any variable determined by the Hamiltonian is given by

(14)

In secular dynamics the semimajor axis is fixed and can be considered to be a function only of the shape of the orbit, as expressed by and . Then from the chain rule

(15)

where is the vector with a similar definition for . Replacing successively by and and using the relations (13), we find

(16)

Since and are constants of motion for the Kepler Hamiltonian , the contribution of to the right side of this equation must vanish, so we can replace by . Equations (16) date back to Milankovich (1939, eqs. 204 and 205; see also and references therein).

These equations admit three integrals of motion, , and . Physically meaningful solutions are restricted to the four-dimensional manifold on which

(17)

Replacing and by the dimensionless variables and defined in equation (11) we obtain222There is a gauge freedom in the definition of since the variables and are related by (17). It is shown in Appendix A that the secular equations of motion (18) are independent of gauge.

(18)

For the potential given by equations (12) we finally have

(19)

To avoid distraction by trivial cases, we shall always assume that (the planetary quadrupole is non-zero and positive, as expected for a planet that is oblate or has inner satellites), (solar perturbations are non-negligible), and is neither parallel, antiparallel, or perpendicular to (the planetary obliquity is not ). The plane defined by the planetary spin axis and the normal to the solar orbit will be called the principal plane.

We shall sometimes use cylindrical coordinates with the -axis oriented along , so the plane coincides with the principal plane. The positive -axis is chosen to coincide with . Then lies in the plane, so it may be specified by its azimuthal angle (the obliquity), which lies in the range .

In the limit equations (19) provide a vector description of Kozai oscillations, the secular oscillations in the eccentricity and inclination of (for example) a planet orbiting a member of a binary star (Kozai, 1962; Holman et al., 1997; Ford et al., 2000).

Equations (19) are invariant under the transformations

(20)

Invariance under implies that we can restrict the range of the obliquity from to .

Recall that equations (19) hold for arbitrary eccentricity, i.e., they do not represent an expansion that is valid only for .

We define the Laplace equilibria to be stationary solutions of equations (19), and the Laplace surface(s) to be the locus of all orbits that are Laplace equilibria.

3 Circular Laplace equilibria

The equations of motion (19) can be solved explicitly in the case of circular orbits (). We banish this discussion to Appendix B and focus in this section on the equilibrium solutions, with and constant, which we call the circular Laplace equilibria.

In circular Laplace equilibria the second of equations (19) is satisfied trivially. The first of equations (19) yields

(21)

Taking the scalar product of this equation with , and again with , we conclude that either (i) ; (ii) . In the first case is perpendicular to the principal plane, and we call this the “orthogonal” or “circular orthogonal” Laplace equilibrium; in the second case lies in the principal plane and we call this the “coplanar” or “circular coplanar” Laplace equilibrium.

In the coplanar Laplace equilibrium, lies in the plane, so it may be specified by its azimuthal angle . The equilibrium condition (21) becomes

(22)

This equation has four solutions for in a interval; if is a solution then and are also solutions.

Equation (22) may also be written as

(23)

where the Laplace radius is defined by

(24)

3.1 Stability

The stability of the circular Laplace equilibria can be determined by writing and expanding equations (19) to first order in and :

(25)

Note that the two equations are decoupled: the linearized evolution of the orientation of the orbital plane, specified by , is independent of the linearized evolution of the eccentricity and apse direction, specified by .

Figure 1: The regions in which stable, circular, coplanar Laplace equilibria exist are marked by heavy black stippling and the label “00”. In the more lightly stippled regions, circular coplanar equilibria exist but are unstable, according to equations (32) and (34). The first figure in each label is “0” or “1” according to whether the solution is stable or unstable to changes in the orbit plane orientation, described by . The second figure is “0” or “1” according to whether the solution is stable or unstable to changes in the eccentricity, described by . The vertical coordinate is the inclination of the orbit relative to the planetary equator, . The horizontal coordinate is the obliquity of the planet relative to the ecliptic, . The results are unchanged if or or .

It is easy to show using the equilibrium equation (21) that a trivial solution of the first of these equations is where is a constant. In other words, the eigenvalue equation in obtained by assuming always has one zero eigenvalue. This is an unphysical solution of the linearized equations since the constant of motion requires that .

3.1.1 The circular orthogonal Laplace equilibrium

In the orthogonal equilibrium equations (25) simplify to

(26)

These can be converted to eigenvalue equations by assuming . To analyze the first equation, set , . Taking the scalar product of with and we find

(27)

these can be combined to show that either (i) , ; this is the unphysical solution noted above; or (ii)

(28)

since in this case, is oscillatory, so the circular orthogonal equilibrium is linearly stable to variations in the angular momentum vector . The quantity in square brackets is just , where is the obliquity.

The second of equations (26) can be analyzed similarly. A trivial solution is , ; this solution is unphysical since the constant of motion requires . The other eigenvalues are given by

(29)

Thus the circular orthogonal equilibrium is stable to variations in eccentricity if and only if

(30)

3.1.2 The circular coplanar Laplace equilibrium

In the circular coplanar equilibrium, the first of the linearized equations (25) implies that where

(31)

Substituting the equilibrium condition (22) we have

(32)

The regions in which (instability) are marked by red and blue stippling and the labels “10” and “11” in Figure 1. All equilibria with are stable to perturbations of this kind (variations in but not ) and equilibria with are unstable.

The circular coplanar equilibria that are stable in this sense have as and as , so the corresponding Laplace surface coincides with the equator of the planet at small radii and with the planet’s orbital plane at large radii. We call this the “classical” Laplace surface since this was the surface discovered by Laplace and the one that has been the focus of most work on this subject.

The second of the linearized equations (25) implies that where

(33)

Substituting the equilibrium solution (22) we have

(34)

The regions in which (instability) are marked by green and blue stippling and the labels “01” and “11” in Figure 1.

The inclination to the planetary equator is plotted as a function of the strength of the planetary quadrupole in Figure 2, for obliquities . Solutions for can be obtained by the transformation .

Figure 2: Circular coplanar Laplace equilibria. The vertical axis is the angle between the pole of the planet and the orbital pole of the satellite (the inclination of the orbit relative to the planet’s equator). Solutions are shown for eight values of the planetary obliquity, . There are two equilibrium curves for each obliquity, one with and the other with . Solid lines denote stable equilibria, while dashed and dotted lines denote equilibria with one or two unstable roots respectively. All equilibria with are unstable. The classical equilibria are those with . The horizontal axis represents the relative strength of perturbations from the solar tide and the planetary quadrupole (eq. 11); the planetary quadrupole dominates on the left side of the figure and the solar quadrupole on the right.
Figure 3: Regions in which the classical Laplace equilibria are unstable are stippled. The top panel shows the unstable range of obliquities as a function of the ratio of the solar and planetary perturbation strengths (eq. 11) and the bottom panel shows the unstable obliquities as a function of semimajor axis (in units of the Laplace radius , eq. 24). All instabilities are in the eccentricity vector; the angular-momentum vector is stable.
Figure 4: Regions in which the classical Laplace equilibria are unstable are shaded.

Unstable equilibria are shown in Figure 2 by dotted or dashed lines. All equilibria with are unstable. In addition, some of the classical equilibria () are unstable to eccentricity growth. These are shown as stippled or shaded regions in Figures 3 and 4. Instability first appears at obliquity and is restricted to semimajor axes between about 0.9 and 1.25 times the Laplace radius (eq. 24).

4 Eccentric Laplace equilibria

We now look for stationary solutions to equations (19) in which the eccentricity is non-zero (recall that ). It can be shown (see Appendix C) that all such solutions either have both and in the principal plane defined by and (the “coplanar-coplanar” or “eccentric coplanar-coplanar” solution), or one of and in the principal plane and the other orthogonal to it (the “coplanar-orthogonal” or “orthogonal-coplanar” equilibrium if or , respectively, lies in the principal plane).

4.1 Eccentric coplanar-coplanar Laplace equilibrium

In this case and lie in the principal plane, and equations (19) with yield

(35)

These can be used to solve for and , given and . The solutions are physical only if .

Figure 5: Regions in which eccentric coplanar-coplanar Laplace equilibria exist are stippled. Heavy blue and light red stippling mark regions in which the equilibria are respectively stable and unstable. Contours mark eccentricities of . Compare to Figure 1.

The regions in space in which eccentric coplanar-coplanar Laplace equilibria exist are shown in Figure 5 by stipples. These regions are bounded, in part, by the lines , , and . Heavy stippling marks regions in which the equilibria are stable.

The stable regions are small for so we will focus on the region . The eccentric coplanar-coplanar equilibria in this region are closely related to the circular coplanar equilibria discussed in §3. We found there that the circular coplanar Laplace equilibria were unstable for a range of when the obliquity . Figure 6 shows the range in which the circular coplanar equilibria are stable for as solid horizontal lines, with gaps marking the unstable range. Superimposed on these lines are the eccentric coplanar-coplanar equilibria, marked by heavy blue or light red curves depending on whether they are stable or unstable. The height of these curves is proportional to the eccentricity of the equilibrium. The figure shows that the eccentric coplanar-coplanar equilibria bifurcate from the circular coplanar equilibrium at the point where the circular coplanar equilibrium becomes unstable. The structure of the secular Hamiltonian at the bifurcation is that of the standard resonant Hamiltonian at orbital resonances ( in the notation of Borderies & Goldreich 1984).

For the eccentric coplanar-coplanar solution exists in a limited range of and is stable throughout this region. For the eccentric coplanar-coplanar solution is unstable for some part of the range of in which it exists. Let us imagine a satellite in the circular coplanar Laplace equilibrium at (say) . If slowly decreases (for example, because the satellite is slowly migrating toward the planet) then we expect that (i) for the satellite will always remain in a circular coplanar Laplace equilibrium; (ii) for the satellite will transfer onto the eccentric coplanar-coplanar equilibrium, its eccentricity will then grow as its spirals in, reach a maximum depending on , then shrink back to zero, at which point it will rejoin the sequence of circular coplanar Laplace equilibria and spiral into the planet on a circular orbit; (iii) for the satellite will transfer onto the eccentric coplanar-coplanar equilibrium, and its eccentricity will grow until the equilibrium orbit becomes unstable, at which point it presumably undergoes large oscillations in eccentricity and inclination. Numerical simulations of this evolution are described in §5; we shall find that (i) is correct, but that (ii) and (iii) need to be qualified in that orbits track the sequence of eccentric equilibria less well than they imply.

An interesting but difficult question, which we do not attempt to answer here, is the behavior of a dissipative gas or particulate disk in the region where the classical Laplace surface is unstable. There are at least three alternatives: (i) The dissipative forces suppress the secular instabilities in the circular coplanar Laplace equilibria, thereby allowing the disk to occupy the classical Laplace surface. (ii) The disk occupies the surface defined by the eccentric coplanar-coplanar equilibria. A possible problem is the large eccentricity gradients in this sequence of equilibria: a flat, cold disk composed of aligned eccentric orbits cannot exist if , because the orbits cross. This condition does not apply directly to the Laplace surface because it is not flat; nevertheless, the large eccentricity gradients may lead to strong shears that destabilize the disk. (iii) A gap opens in the disk where the classical Laplace surface is unstable.

Figure 6: The relation of the circular coplanar Laplace equilibria to the eccentric coplanar-coplanar equilibria. The gaps in the solid horizontal lines represent regions in which the circular coplanar equilibria are unstable. The curved lines represent the eccentric coplanar-coplanar Laplace equilibria with inclination . The -coordinate of these lines is where is the eccentricity of the solution; the curves are heavy blue or light red according to whether the solution is stable or unstable. The figure shows that the eccentric coplanar-coplanar equilibria bifurcate from the circular coplanar equilibria.

4.2 Eccentric coplanar-orthogonal Laplace equilibrium

In this case lies in the principal plane and is normal to it, so

The first of equations (19) with yields

(36)

We introduce a vector defined by or ; lies in the principal plane and . Substituting for in the second of equations (19) with and yields

(37)

The vectors and are linearly independent, so the coefficient of each must be zero. The coefficient of vanishes if equation (36) is satisfied. The condition that the coefficient of vanishes is

(38)

Equations (36) and (38) can be combined to eliminate , , and :

(39)

This result determines the inclination between the satellite orbit and the planetary spin in terms of the obliquity , independent of , , or . Numerical solution of this equation for shows that there are two solutions for each value of , but the upper solution is unphysical, because these values of and yield in equation (36). Thus there is a unique inclination for each value of the obliquity , as shown in the top panel of Figure 7.

Figure 7: (Top) Inclination as a function of the obliquity for the eccentric coplanar-orthogonal Laplace equilibrium, as obtained by solving equation (39). (Bottom) Eccentricity of the coplanar-orthogonal equilibrium for obliquity (left to right). The dashed lines indicate that all of these equilibria have one unstable mode.

The bottom panel of Figure 7 shows the eccentricity of the coplanar-orthogonal Laplace equilibria. All of these equilibria are unstable.

4.3 Eccentric orthogonal-coplanar Laplace equilibrium

In this case lies in the principal plane and is normal to it, so

(40)

The solution of equations (19) with requires that

(41)

Solutions exist whenever .

Linear stability analysis shows that small perturbations in or grow as where

(42)

since the eccentric orthogonal-coplanar solutions are stable.

Comparison with the results of §3.1.1 shows that these (eccentric) equilibria bifurcate from the (circular) orthogonal equilibrium sequence at the semimajor axis where the latter becomes unstable.

Figure 8: Numerical integrations of the equations of motion (19). The -coordinate is where is the eccentricity of the solution. Each orbit is begun at semimajor axis and the semimajor axis then shrinks according to where is defined in equation (11). The satellite is initially in the circular Laplace equilibrium except for a seed eccentricity of . The -coordinate is . Four integrations are shown, for obliquities .

5 Orbital evolution

Consider the evolution of a satellite on a circular orbit that is slowly decaying, so the satellite is migrating toward the planet. Let us assume that the initial semimajor axis is much larger than the Laplace radius (eq. 24) and that the orbital plane coincides with the planetary orbit (i.e., the satellite is in the classical Laplace surface). As outlined at the end of §4.1, if the obliquity we expect that the satellite will remain in a circular coplanar Laplace equilibrium as it spirals in. For we expect that at semimajor axis (Fig. 3) the satellite will transfer onto the eccentric coplanar-coplanar equilibrium; its eccentricity will then grow as it migrates inward, reach a maximum depending on , then shrink back to zero, at which point () it will rejoin the classical Laplace surface and spiral into the planet on a circular orbit. Finally, for the satellite will transfer onto the eccentric coplanar-coplanar equilibrium, and its eccentricity will grow until its orbit becomes unstable and begins to execute large oscillations in eccentricity and inclination

Numerical integrations of the equations of motion (19) are shown for a migrating satellite in Figure 8. As expected, when the initial obliquity is , , or large eccentricity oscillations develop in the region where the circular coplanar Laplace equilibrium is unstable (compare Fig. 6). Smaller oscillations develop in the case ; these are unexpected since the eccentric coplanar-coplanar sequence is stable at this obliquity, so we would expect smooth growth and decline in the eccentricity as the semimajor axis declines, as in Figure 6. Further integrations show that the behavior of migrating satellites at this obliquity depends on the migration time and other parameters; for example, the orbit shown in Figure 9 follows the eccentric coplanar-coplanar sequence for a while and then rather suddenly jumps to an orbit with a chaotic appearance that reaches eccentricities as high as 0.45. These results suggest the presence of large chaotic regions and perhaps higher order resonances in the phase space, but we have not yet explored these features.

Figure 9: Numerical integration of the equations of motion (19) as a satellite migrates inward. The obliquity is and the satellite semimajor axis shrinks according to where is defined in equation (11). The satellite is initially in the circular Laplace equilibrium except for a seed eccentricity of . The -coordinate is . The satellite initially lies on the classical Laplace surface; at it transitions to the eccentric coplanar-coplanar sequence; and at it transitions to an orbit that exhibits large and irregular eccentricity oscillations.

6 Applications

Direct applications of these results to the solar system are rather few. Only Uranus and Pluto have obliquities that exceed the critical value at which the classical Laplace surface becomes unstable, and these do not have satellites close to the unstable range of semimajor axis.

The circular orthogonal Laplace equilibria have been suggested as possible sites for “polar” rings around Neptune (see Dobrovolskis et al. 1989 and references therein). In this case the solar quadrupole potential is much weaker than the potential from Neptune’s satellite Triton. However, the effects of Triton can be modeled approximately using the formalism we have derived, simply replacing the Sun by Triton. We find from equation (24) that the Neptune-Triton Laplace radius is . Equation (30) then implies that circular polar rings would be unstable at semimajor axes , and beyond this radius polar rings must be eccentric. This analysis is only approximate, since (i) the Laplace radius is 61% of Triton’s semimajor axis, , so the approximation of Triton’s potential by its quadrupole component is poor; (ii) Triton’s orbital axis and Neptune’s spin axis both precess around their mutual invariable plane, whereas we have assumed that they are fixed; this is probably only a minor correction, since the precession period of is much longer than the growth time of the instability, at ; (iii) our calculation neglects collective effects in the ring, such as interparticle collisions, which might suppress the instability.

Our original motivation for examining this problem was the orbit of Saturn’s satellite Iapetus, which has an eccentricity of 0.028 and an inclination of to the Laplace surface. If it formed from a circumplanetary disk, one might expect Iapetus to have zero eccentricity and inclination relative to this surface. Thus it appears that some process has pumped up Iapetus’s inclination while leaving its eccentricity near zero. Ward (1981) pointed out that the shape of the classical Laplace surface is affected by the mass in the circumplanetary disk, and suggested that the current orbit of Iapetus reflects its shape before the disk dispersed. However, this scenario requires the dispersal of the disk in yr—if the dispersal were slower, the inclination relative to the Laplace surface would be an adiabatic invariant and thus would remain near zero. The semimajor axis of Iapetus is not far from the Laplace radius (Table 1) so it is natural to wonder whether instabilities have excited the inclination. However, (i) at Saturn’s obliquity of , all circular orbits in the Laplace surface are stable; (ii) the instability we have found in the Laplace surface at high obliquity tends to excite eccentricity, not inclination.

The rich dynamics described in this paper is likely to play a larger role in scenarios in which the obliquities of the giant planets have changed substantially since the formation of the solar system (e.g., Tremaine, 1991).

A further application is to extrasolar planetary systems. It is well-known that Kozai oscillations can excite the eccentricity of a single planet orbiting one member of a binary star system (Kozai, 1962; Holman et al., 1997; Ford et al., 2000). Now consider a star that hosts two planets, one a hot Jupiter, and belongs to a binary system (see Takeda et al. 2008 for numerical simulations of such systems). The hot Jupiter augments the (otherwise negligible) quadrupole moment of the host star, just as inner satellites augment the quadrupole moment of a planet (eq. 4). The resulting Laplace radius (24) is

(43)

where and are the mass and semimajor axis of the hot Jupiter, and , , and are the mass, semi-major axis, and eccentricity of the binary companion. If, as we expect, the orbital axes of the hot Jupiter and the distant companion star are uncorrelated, 36% of binary stars will have obliquities that exceed the critical value at which instabilities set in at some semimajor axis; and if planetary migration is common we expect that many planets now in the habitable zone will have passed through the unstable region and thereby acquired substantial eccentricities (cf. Fig. 8).

7 Summary

We have examined the orbital dynamics of planetary satellites under the combined influence of the quadrupole moment from the planet’s equatorial bulge and the tidal field from the Sun. The Laplace equilibria are orbits in which the secular evolution due to these forces is zero. They represent the orbits in which planetary satellites formed from a circumplanetary gas disk should be found.

Laplace equilibria exist for both circular and eccentric orbits. At a given semimajor axis, the orbit normals to the circular Laplace equilibria lie along three orthogonal directions, two of them in the principal plane defined by the planet’s spin axis and the normal to its orbit around the Sun (the “circular coplanar” equilibria). One of the two circular coplanar equilibria is always unstable. The warped surface swept out by the other circular coplanar equilibrium orbits as the semimajor axis varies is called the classical Laplace surface (eq. 23). The classical Laplace surface coincides with the equator of the planet at small distances and with the orbital plane of the planet at large distances. Orbits in the classical Laplace surface are stable if the planetary obliquity , while in the range there is a range of semimajor axes near the Laplace radius (eq. 24) in which the surface is unstable. The third Laplace equilibrium for circular orbits (the “orthogonal” or “polar” equilibrium) corresponds to orbits that cross over the planet’s pole, and these are stable if and only if the semimajor axis is less than (eq. 30).

At some obliquities and semimajor axes eccentric Laplace equilibria also exist. In the “coplanar-coplanar” equilibria both the eccentricity and angular-momentum vectors lie in the principal plane (Figures 5 and 6). The “orthogonal-coplanar” equilibria are stable, eccentric polar orbits that bifurcate from the circular orthogonal equilibria at the semimajor axis where these become unstable. Other eccentric equilibria exist but are unstable.

The use of the secular equations of motion (19) should be legitimate so long as the precession times for the eccentricity and angular-momentum vector are much longer than the orbital period of the satellite. A sufficient condition for this is that the satellite semimajor axis is much larger than the planetary radius and much smaller than the planet’s Hill radius (cf. Table 1). A possible limitation of the secular equations is that they neglect the evection resonance, where the apsidal precession rate equals the mean motion of the planet around the Sun. For the outer planets, the evection resonance typically occurs at . Touma & Wisdom (1998) have stressed the role of the evection resonance in the history of the lunar orbit.

Several unanswered questions remain: (i) What is the structure of the four-dimensional phase space of the equations of motion (19), and what fraction of the orbits are chaotic? These issues could be explored through surfaces of section, although the exploration would be laborious since there is a separate surface of section for each value of the Hamiltonian, semimajor axis, and planetary obliquity. (ii) What is the nature of the evolution of satellites in the classical Laplace surface that migrate into an unstable region (Fig. 8)? (iii) What are the shape and properties of dissipative disks in the range of semimajor axes where the classical Laplace surface is unstable? Does the dissipation stabilize the classical Laplace surface? Does the disk follow the eccentric coplanar-coplanar Laplace surface? Or is there no steady-state disk? Hints that eccentric structures can occur in dissipative disks include the eccentric structures seen in the planetary ring systems of Uranus and Saturn, debris disks around young stars, and galactic nuclei such as M31. (iv) What is the analog of the Laplace surface in accretion disks around black holes, where Lense–Thirring precession rather than a quadrupole potential is the dominant non-Keplerian perturbation from the central body?

We thank Yuri Levin for thoughtful comments. This research was supported in part by NASA grants NNX08AH83G and NNX08AH24G, and by NSF grant AST-0807432.

Appendix A Gauge dependence in the averaged potential

The orbit-averaged potential is a function of the six scalar components of and but is only physically meaningful on the four-dimensional manifold defined by the constraints (17). Thus may be replaced by , where the gauge function is arbitrary except that on the manifold (17). We now show that has no effect on the equations of motion (18), just as adding a constant to a Newtonian potential has no effect on the equations of motion .

We work in the coordinates with axes parallel to the unit vectors defined just before equation (7). Substituting for in the equations of motion (18) and observing that , , we have

(A1)

The gauge condition that on the manifold (17) is

(A2)

We write and in terms of the unit vectors