**Secular Evolution of Satellites by Tidal Effect**

Alexandre C. M. Correia *University of Aveiro Portugal*

#### **1. Introduction**

102 Space Science

doi.10.1029/2002GLO15685 pp. 27-1 to 27-4, 2002. Zeik/Gauntand editors *Astronomy* (IInd Edition), Cosmic Perspective, Zimmerman, Robert " Exo-Earths", *Astronomy,* August 2004, pp.42-47;

Oscillations",*Geophysical Research Letters*, Vol. 29, No 21, 2031,

Both the Earth's Moon and Pluto's moon, Charon, have an important fraction of the mass of their systems, and therefore they could be classified as double-planets rather than as satellites. The proto-planetary disk is unlikely to produce such systems, and their origin seems to be due to a catastrophic impact of the initial planet with a body of comparable dimensions (e.g. Canup, 2005; Canup & Asphaug, 2001). On the other hand, Neptune's moon, Triton, and the Martian moon, Phobos, are spiraling down into the planet, clearly indicating that the present orbits are not primordial, and may have undergone a long evolving process from a previous capture (e.g. Goldreich et al., 1989; Mignard, 1981).

The present orbits of all these satellites are almost circular, and their spins appear to be synchronous with the orbital mean motion, as well as being locked in Cassini states (e.g. Colombo, 1966; Peale, 1969). This also applies to the Galilean satellites of Jupiter, which are likely to have originated from Jupiter's accretion disk and additionally show orbital mean motion resonances (e.g. Yoder, 1979). All these features seem to be due to tidal evolution, which arises from differential and inelastic deformation of the planet by a perturbing body.

Previous long-term studies on the orbital evolution of satellites have assumed that their rotation is synchronously locked, and therefore limits the tidal evolution to the orbits (e.g. McCord, 1966). However, these two kinds of evolution cannot be dissociated because the total angular momentum must be conserved. Additionally, it has been assumed that the spin axis is locked in a Cassini state with very low obliquity. Although these assumptions are correct for the presently known situations, they were not necessarily true throughout the evolution.

In this article we model the orbital evolution of a satellite from its origin or capture until the preset day, including spin evolution for both planet and satellite, and we also regard its future evolution. We provide a simple averaged model adapted for fast computational simulations, as required for long-term studies, following Correia (2009). However, we present an improvement with respect to previous work, here we do not average the equations of motion over the argument of the periastron, as in Correia et al. (2011). Therefore, this model is more complete, and allows the eccentricity of the satellite to show secular variations due to the gravitational perturbations of the star on its orbit around the planet. We then apply this model to the Triton-Neptune system. The results do not differ much from those in Correia (2009) for the final stages of the orbital evolution, but can show some significant differences during the initial stages. In the last section we discuss the results obtained.

Since the total angular momentum is conserved, the contributions to the spin of the planet

Secular Evolution of Satellites by Tidal Effect 105

Because we are only interested in the secular evolution of the system, we further average the equations of motion over the mean anomalies of both orbits. The resulting equations for the

<sup>1</sup>) cos *<sup>I</sup>* **kˆ** <sup>1</sup> <sup>×</sup> **kˆ** <sup>0</sup> <sup>+</sup> <sup>5</sup>*γ*(**e**<sup>1</sup> · **kˆ** <sup>0</sup>) **<sup>e</sup>**<sup>1</sup> <sup>×</sup> **kˆ** <sup>0</sup> <sup>−</sup> ∑

<sup>1</sup>) cos *<sup>I</sup>* **kˆ** <sup>0</sup> <sup>×</sup> **kˆ** <sup>1</sup> <sup>+</sup> <sup>5</sup>*γ*(**e**<sup>1</sup> · **kˆ** <sup>0</sup>) **kˆ** <sup>0</sup> <sup>×</sup> **<sup>e</sup>**<sup>1</sup> <sup>−</sup> ∑

1 2

*<sup>α</sup><sup>i</sup>* <sup>=</sup> <sup>3</sup>*GMmiJ*<sup>2</sup>*<sup>i</sup>*

*<sup>β</sup><sup>i</sup>* <sup>=</sup> <sup>3</sup>*Gm*0*m*<sup>1</sup> *<sup>J</sup>*2*<sup>i</sup>*

*<sup>γ</sup>* <sup>=</sup> <sup>3</sup>*GMm*1*a*<sup>2</sup>

are the direction cosines of the spins and orbits: *ε<sup>i</sup>* is the obliquity to the orbital plane of the planet, *θ<sup>i</sup>* the obliquity to the orbital plane of the satellite, and *I* the inclination between orbital

*<sup>f</sup>*4(*e*1) **kˆ** <sup>1</sup> <sup>×</sup> **<sup>e</sup>**<sup>1</sup>

11

<sup>2</sup> *<sup>f</sup>*4(*e*1) cos *<sup>θ</sup><sup>i</sup>*

2*a*<sup>3</sup> <sup>0</sup>(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

2*a*<sup>3</sup> <sup>1</sup>(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

4*a*<sup>3</sup> <sup>0</sup>(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

cos *<sup>I</sup>* **kˆ** <sup>0</sup> <sup>×</sup> **<sup>e</sup>**<sup>1</sup> <sup>−</sup> <sup>2</sup> **kˆ** <sup>1</sup> <sup>×</sup> **<sup>e</sup>**<sup>1</sup> <sup>−</sup> <sup>5</sup>(**e**<sup>1</sup> · **kˆ** <sup>0</sup>) **kˆ** <sup>0</sup> <sup>×</sup> **kˆ** <sup>1</sup>

(<sup>1</sup> <sup>−</sup> 5 cos2 *<sup>θ</sup>i*) **kˆ** <sup>1</sup> <sup>×</sup> **<sup>e</sup>**<sup>1</sup>

*R*2 *i*

> *R*2 *i*

1

cos *<sup>ε</sup><sup>i</sup>* <sup>=</sup> **sˆ***<sup>i</sup>* · **kˆ** <sup>0</sup> , cos *<sup>θ</sup><sup>i</sup>* <sup>=</sup> **sˆ***<sup>i</sup>* · **kˆ** <sup>1</sup> , cos *<sup>I</sup>* <sup>=</sup> **kˆ** <sup>0</sup> · **kˆ** <sup>1</sup> , (15)

**L˙** <sup>0</sup> <sup>=</sup> 0 , **L˙** <sup>1</sup> <sup>=</sup> <sup>−</sup>**H˙** <sup>0</sup> <sup>−</sup> **H˙** <sup>1</sup> , (16)

*ωi n*1

− 9 *f*5(*e*1)

 **e**1 

, (17)

**H˙** *<sup>i</sup>* <sup>=</sup> <sup>−</sup>*α<sup>i</sup>* cos *<sup>ε</sup><sup>i</sup>* **kˆ** <sup>0</sup> <sup>×</sup> **sˆ***<sup>i</sup>* <sup>−</sup> *<sup>β</sup><sup>i</sup>* cos *<sup>θ</sup><sup>i</sup>* **kˆ** <sup>1</sup> <sup>×</sup> **sˆ***<sup>i</sup>* , (11)

**H˙** <sup>0</sup> + **H˙** <sup>1</sup> + **L˙** <sup>0</sup> + **L˙** <sup>1</sup> = 0 . (7)

*i*

*i*

*<sup>α</sup><sup>i</sup>* cos *<sup>ε</sup><sup>i</sup>* **sˆ***<sup>i</sup>* <sup>×</sup> **kˆ** <sup>0</sup> , (8)

*<sup>β</sup><sup>i</sup>* cos *<sup>θ</sup><sup>i</sup>* **sˆ***<sup>i</sup>* <sup>×</sup> **kˆ** <sup>1</sup> , (9)

, (10)

0)3/2 , (12)

1)3/2 , (13)

0)3/2 , (14)

and satellite can easily be computed from the orbital contributions:

conservative motion are (Boué & Laskar, 2006; Farago & Laskar, 2010):

cos *θ<sup>i</sup>* **sˆ***<sup>i</sup>* × **e**<sup>1</sup> +

For the dissipative tidal effects, we obtain (Correia et al., 2011):

 *Ri a*1 <sup>5</sup>

(**e**<sup>1</sup> · **sˆ***i*) **kˆ** <sup>1</sup> <sup>−</sup>

**L˙** <sup>0</sup> <sup>=</sup> <sup>−</sup>*γ*(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*

**L˙** <sup>1</sup> <sup>=</sup> <sup>−</sup>*γ*(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*

**e˙**<sup>1</sup> <sup>=</sup> <sup>−</sup>*γ*(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*

− ∑ *i*

and

where

and

planes.

**e˙**<sup>1</sup> = ∑ *i*

> − ∑ *i*

15 2 *k*2*i n*1 *mj mi*

*Ki β*1*a*<sup>2</sup> 1 *<sup>f</sup>*4(*e*1) *<sup>ω</sup><sup>i</sup>* 2*n*<sup>1</sup>

2

2

2 1)


*βi* ||**L**1||

#### **2. The model**

We consider a hierarchical system composed of a star, a planet and a satellite, with masses *M* � *m*<sup>0</sup> � *m*1, respectively. Both planet and satellite are considered oblate ellipsoids with gravity field coefficients given by *J*<sup>20</sup> and *J*21 , rotating about the axis of maximal inertia along the directions **sˆ**<sup>0</sup> and **sˆ**1, with rotation rates *ω*<sup>0</sup> and *ω*1, respectively. The potential energy *U* of the system is then given by (e.g. Smart, 1953):

$$\mathcal{U} = -G\frac{Mm\_0}{r\_0} \left( 1 - \sum\_{i=0,1} f\_{2\_i} \frac{m\_i}{m\_0} \left( \frac{R\_i}{r\_0} \right)^2 P\_2(\mathbf{\hat{r}}\_0 \cdot \mathbf{\hat{s}}\_i) \right)$$

$$-G\frac{m\_0 m\_1}{r\_1} \left( 1 - \sum\_{i=0,1} f\_{2\_i} \left( \frac{R\_i}{r\_1} \right)^2 P\_2(\mathbf{\hat{r}}\_1 \cdot \mathbf{\hat{s}}\_i) \right)$$

$$-G\frac{Mm\_1}{r\_0} \left( \frac{r\_1}{r\_0} \right)^2 P\_2(\mathbf{\hat{r}}\_0 \cdot \mathbf{\hat{r}}\_1) \tag{1}$$

where terms in (*Ri*/*rj*)<sup>3</sup> have been neglected (*i*, *j* = 0, 1). *G* is the gravitational constant, *Ri* the radius of the planet or the satellite, **r***<sup>i</sup>* the distance between the planet and the star or the satellite, and *<sup>P</sup>*2(*x*)=(3*x*<sup>2</sup> <sup>−</sup> <sup>1</sup>)/2 the Legendre polynomial of degree two.

Neglecting tidal interactions with the star, the tidal potential is written (e.g. Kaula, 1964):

$$\mathcal{U}\_T = -\frac{G}{r\_1^3} \sum\_{i=0,1} k\_{2\_i} m\_{(1-i)}^2 \frac{R\_i^5}{r\_i^{73}} P\_2(\mathbf{\hat{r}}\_1 \cdot \mathbf{\hat{r}}\_i') \,, \tag{2}$$

where *k*2*<sup>i</sup>* is the potential Love number for the planet or the satellite, and **r**� *<sup>i</sup>* the position of the interacting body at a time delayed of Δ*ti*. For simplicity, we will adopt a model with constant Δ*ti*, which can be made linear (e.g. Mignard, 1979; Néron de Surgy & Laskar, 1997):

$$\mathbf{r}'\_{i} \simeq \mathbf{r}\_{1} + \Delta t\_{i} \left(\omega\_{i} \mathbf{s}\_{i} \times \mathbf{r}\_{1} - \dot{\mathbf{r}}\_{1}\right) \tag{3}$$

The complete evolution of the system can be tracked by the evolution of the rotational angular momentums, **<sup>H</sup>***<sup>i</sup>* � *Ciωi***sˆ***i*, the orbital angular momentums, **<sup>L</sup>***<sup>i</sup>* � *minia*<sup>2</sup> *<sup>i</sup>* (<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup> *<sup>i</sup>* )1/2**kˆ***i*, and the Laplace-Runge-Lenz vector, which points along the major axis in the direction of periapsis with magnitude *e*1:

$$\mathbf{e}\_1 = \frac{\mathbf{\dot{r}}\_1 \times \mathbf{L}\_1}{GMm\_1} - \frac{\mathbf{r}\_1}{r\_1} \,. \tag{4}$$

*ni* is the mean motion, *ai* the semi-major axis, *ei* the eccentricity, and *Ci* the principal moment of inertia. The contributions to the orbits are easily computed from the above potentials as

$$
\dot{\mathbf{L}}\_0 = -\mathbf{r}\_0 \times \mathbf{F}\_0 \, \prime \quad \dot{\mathbf{L}}\_1 = -\mathbf{r}\_1 \times \mathbf{F}\_1 \, \prime \,\tag{5}
$$

$$
\dot{\mathbf{e}}\_1 = \frac{1}{GMm\_1} \left( \mathbf{F}\_1 \times \frac{\mathbf{L}\_1}{m\_1} + \dot{\mathbf{r}}\_1 \times \dot{\mathbf{L}}\_1 \right) \tag{6}
$$

where **F***<sup>i</sup>* = −∇**r***<sup>i</sup> U*� , with *U*� = *U* + *UT* + *GMm*0/*r*<sup>0</sup> + *Gm*0*m*1/*r*1. Since the total angular momentum is conserved, the contributions to the spin of the planet and satellite can easily be computed from the orbital contributions:

$$
\dot{\mathbf{H}}\_0 + \dot{\mathbf{H}}\_1 + \dot{\mathbf{L}}\_0 + \dot{\mathbf{L}}\_1 = 0 \,. \tag{7}
$$

Because we are only interested in the secular evolution of the system, we further average the equations of motion over the mean anomalies of both orbits. The resulting equations for the conservative motion are (Boué & Laskar, 2006; Farago & Laskar, 2010):

$$\dot{\mathbf{L}}\_{0} = -\gamma (1 - e\_{1}^{2}) \cos I \hat{\mathbf{k}}\_{1} \times \hat{\mathbf{k}}\_{0} + \dots \hat{\mathbf{s}}\_{7} (\mathbf{e}\_{1} \cdot \hat{\mathbf{k}}\_{0}) \cdot \mathbf{e}\_{1} \times \hat{\mathbf{k}}\_{0} - \sum\_{i} a\_{i} \cos \varepsilon\_{i} \cdot \hat{\mathbf{s}}\_{i} \times \hat{\mathbf{k}}\_{0} \,\tag{8}$$

$$\dot{\mathbf{L}}\_{1} = -\gamma (1 - e\_{1}^{2}) \cos I \hat{\mathbf{k}}\_{0} \times \hat{\mathbf{k}}\_{1} + \dots \hat{\mathbf{s}}\_{7} (\mathbf{e}\_{1} \cdot \hat{\mathbf{k}}\_{0}) \hat{\mathbf{k}}\_{0} \times \mathbf{e}\_{1} - \sum\_{i} \beta\_{i} \cos \theta\_{i} \cdot \hat{\mathbf{s}}\_{i} \times \hat{\mathbf{k}}\_{1} \,\tag{9}$$

$$
\dot{\mathbf{e}}\_1 = -\frac{\gamma (1 - e\_1^2)}{||\mathbf{L}\_1||} \left[ \cos I \hat{\mathbf{k}}\_0 \times \mathbf{e}\_1 - 2 \hat{\mathbf{k}}\_1 \times \mathbf{e}\_1 - 5 (\mathbf{e}\_1 \cdot \hat{\mathbf{k}}\_0) \hat{\mathbf{k}}\_0 \times \hat{\mathbf{k}}\_1 \right]
$$

$$
$$

and

2 Will-be-set-by-IN-TECH

We consider a hierarchical system composed of a star, a planet and a satellite, with masses *M* � *m*<sup>0</sup> � *m*1, respectively. Both planet and satellite are considered oblate ellipsoids with gravity field coefficients given by *J*<sup>20</sup> and *J*21 , rotating about the axis of maximal inertia along the directions **sˆ**<sup>0</sup> and **sˆ**1, with rotation rates *ω*<sup>0</sup> and *ω*1, respectively. The potential energy *U* of

> *Ri r*0

> > <sup>2</sup>

*R*5 *i r*�3 *i*

*P*2(**rˆ**<sup>1</sup> · **rˆ** �

*<sup>i</sup>* � **r**<sup>1</sup> + Δ*ti* (*ωi***s***<sup>i</sup>* × **r**<sup>1</sup> − **r˙**1) . (3)

**L˙** <sup>0</sup> <sup>=</sup> <sup>−</sup>**r**<sup>0</sup> <sup>×</sup> **<sup>F</sup>**<sup>0</sup> , **L˙** <sup>1</sup> <sup>=</sup> <sup>−</sup>**r**<sup>1</sup> <sup>×</sup> **<sup>F</sup>**<sup>1</sup> , (5)

<sup>+</sup> **r˙**<sup>1</sup> <sup>×</sup> **L˙** <sup>1</sup>

<sup>2</sup>

*P*2(**rˆ**<sup>0</sup> · **sˆ***i*)

*P*2(**rˆ**<sup>0</sup> · **rˆ**1) , (1)

*<sup>i</sup>*) , (2)

*<sup>i</sup>* (<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

. (4)

, (6)

*<sup>i</sup>* the position of the

*<sup>i</sup>* )1/2**kˆ***i*, and the

*P*2(**rˆ**<sup>1</sup> · **sˆ***i*)

**2. The model**

with magnitude *e*1:

where **F***<sup>i</sup>* = −∇**r***<sup>i</sup>*

*U*�

the system is then given by (e.g. Smart, 1953):

*<sup>U</sup>* <sup>=</sup> <sup>−</sup>*<sup>G</sup> Mm*<sup>0</sup>

*r*0

<sup>−</sup>*<sup>G</sup> <sup>m</sup>*0*m*<sup>1</sup> *r*1

<sup>−</sup>*<sup>G</sup> Mm*<sup>1</sup> *r*0

*UT* <sup>=</sup> <sup>−</sup> *<sup>G</sup> r*3 1 ∑ *i*=0,1 *k*2*i m*2 (1−*i*)

> **r** �

**e˙**<sup>1</sup> <sup>=</sup> <sup>1</sup>

*GMm*<sup>1</sup>

 *r*<sup>1</sup> *r*0 <sup>2</sup>

satellite, and *<sup>P</sup>*2(*x*)=(3*x*<sup>2</sup> <sup>−</sup> <sup>1</sup>)/2 the Legendre polynomial of degree two.

where *k*2*<sup>i</sup>* is the potential Love number for the planet or the satellite, and **r**�

momentums, **<sup>H</sup>***<sup>i</sup>* � *Ciωi***sˆ***i*, the orbital angular momentums, **<sup>L</sup>***<sup>i</sup>* � *minia*<sup>2</sup>

<sup>1</sup> − ∑ *i*=0,1 *J*2*i mi m*<sup>0</sup>

<sup>1</sup> − ∑ *i*=0,1 *J*2*i Ri r*1

where terms in (*Ri*/*rj*)<sup>3</sup> have been neglected (*i*, *j* = 0, 1). *G* is the gravitational constant, *Ri* the radius of the planet or the satellite, **r***<sup>i</sup>* the distance between the planet and the star or the

interacting body at a time delayed of Δ*ti*. For simplicity, we will adopt a model with constant

The complete evolution of the system can be tracked by the evolution of the rotational angular

Laplace-Runge-Lenz vector, which points along the major axis in the direction of periapsis

*ni* is the mean motion, *ai* the semi-major axis, *ei* the eccentricity, and *Ci* the principal moment of inertia. The contributions to the orbits are easily computed from the above potentials as

> **L**1 *m*<sup>1</sup>

− **r**1 *r*1

**<sup>e</sup>**<sup>1</sup> <sup>=</sup> **r˙**<sup>1</sup> <sup>×</sup> **<sup>L</sup>**<sup>1</sup> *GMm*<sup>1</sup>

> **F**<sup>1</sup> ×

, with *U*� = *U* + *UT* + *GMm*0/*r*<sup>0</sup> + *Gm*0*m*1/*r*1.

Δ*ti*, which can be made linear (e.g. Mignard, 1979; Néron de Surgy & Laskar, 1997):

Neglecting tidal interactions with the star, the tidal potential is written (e.g. Kaula, 1964):

$$\dot{\mathbf{H}}\_{i} = -\boldsymbol{\omega}\_{i}\cos\varepsilon\_{i}\,\hat{\mathbf{k}}\_{0} \times \mathbf{s}\_{i} - \beta\_{i}\cos\theta\_{i}\,\hat{\mathbf{k}}\_{1} \times \mathbf{s}\_{i} \,\tag{11}$$

where

$$\alpha\_{i} = \frac{3GMm\_{i}I\_{2\_{i}}R\_{i}^{2}}{2a\_{0}^{3}(1-e\_{0}^{2})^{3/2}},$$

$$\beta\_i = \frac{3Gm\_0m\_1J\_{2\_i}R\_i^2}{2a\_1^3(1-e\_1^2)^{3/2}}\,\text{,}\tag{13}$$

$$\gamma = \frac{3GMm\_1a\_1^2}{4a\_0^3(1-e\_0^2)^{3/2}}\,\text{,}\tag{14}$$

and

$$\cos \varepsilon\_{\mathbf{i}} = \mathbf{\hat{s}}\_{\mathbf{i}} \cdot \mathbf{\hat{k}}\_{0} \,, \quad \cos \theta\_{\mathbf{i}} = \mathbf{\hat{s}}\_{\mathbf{i}} \cdot \mathbf{\hat{k}}\_{1} \,, \quad \cos I = \mathbf{\hat{k}}\_{0} \cdot \mathbf{\hat{k}}\_{1} \,, \tag{15}$$

are the direction cosines of the spins and orbits: *ε<sup>i</sup>* is the obliquity to the orbital plane of the planet, *θ<sup>i</sup>* the obliquity to the orbital plane of the satellite, and *I* the inclination between orbital planes.

For the dissipative tidal effects, we obtain (Correia et al., 2011):

$$
\dot{\mathbf{L}}\_0 = 0 \; \prime \quad \dot{\mathbf{L}}\_1 = -\dot{\mathbf{H}}\_0 - \dot{\mathbf{H}}\_1 \; \prime \tag{16}
$$

$$\begin{split} \dot{\mathbf{e}}\_{1} &= \sum\_{i} \frac{15}{2} k\_{2\_{i}} n\_{1} \left( \frac{m\_{j}}{m\_{i}} \right) \left( \frac{R\_{i}}{a\_{1}} \right)^{5} f\_{4}(e\_{1}) \,\hat{\mathbf{k}}\_{1} \times \mathbf{e}\_{1} \\ &- \sum\_{i} \frac{K\_{i}}{\beta\_{1} a\_{1}^{2}} \left[ f\_{4}(e\_{1}) \frac{\omega\_{i}}{2n\_{1}} (\mathbf{e}\_{1} \cdot \mathbf{s}\_{i}) \,\hat{\mathbf{k}}\_{1} - \left( \frac{11}{2} f\_{4}(e\_{1}) \cos \theta\_{i} \frac{\omega\_{i}}{n\_{1}} - 9f\_{5}(e\_{1}) \right) \mathbf{e}\_{1} \right], \end{split} \tag{17}$$

**3.1 Spin evolution**

attained for:

giving (Correia & Laskar, 2009):

*<sup>ω</sup>*˙ <sup>1</sup> <sup>=</sup> <sup>−</sup>*K*<sup>1</sup> *<sup>n</sup>*<sup>1</sup> *C*1

*ω*1 *n*1

*dt* <sup>=</sup> **H˙** *<sup>i</sup>* · (**kˆ** <sup>1</sup> <sup>−</sup> cos *<sup>θ</sup>i***sˆ***i*)

can be found from a single relationship (e.g. Ward & Hamilton, 2004):

 sin *I*<sup>0</sup> cos *I*<sup>0</sup> ± *λ*<sup>1</sup>

> sin *θ*<sup>1</sup>

The obliquity variations can be obtained from equation (15):

*d* cos *θ<sup>i</sup>*

solutions are approximately given by:

with ||**H**1|| � ||**L**1||, giving:

tan−<sup>1</sup>

˙ *<sup>θ</sup>*<sup>1</sup> � *<sup>K</sup>*1*n*<sup>1</sup> *C*1*ω*<sup>1</sup>

have four real roots approximately given by expressions (30).

*f*1(*e*1)

The variation in the satellite's rotation rate can be obtained from *<sup>ω</sup>*˙ *<sup>i</sup>* <sup>=</sup> **H˙** *<sup>i</sup>* · **sˆ***i*/*Ci* (Eq. 25),

Secular Evolution of Satellites by Tidal Effect 107

1 + cos2 *θ*<sup>1</sup> 2

For a given obliquity and eccentricity, the equilibrium rotation rate, obtained when *ω*˙ <sup>1</sup> = 0, is


For the conservative motion (Eqs. 9, 11), stable configurations for the spin can be found whenever the vectors (**sˆ**1, **kˆ** 1, **kˆ** 0) or (**sˆ**1, **kˆ** 1, **sˆ**0) are coplanar and precess at the same rate *g* (e.g. Colombo, 1966; Correia et al., 2011; Peale, 1969). The first situation occurs if *γ* � *β*<sup>0</sup> (outer satellite) and the second situation when *γ* � *β*<sup>0</sup> (inner satellite). The equilibrium obliquities

where *λ*<sup>1</sup> = *β*1/(*C*1*ω*1*g*) is a dimensionless parameter and *I*<sup>0</sup> is the inclination of the orbit of the satellite with respect to the Laplacian plane (*I*<sup>0</sup> � *I* and *g* � *γ* cos *I*/||**L**1|| for an outer satellite, and *I*<sup>0</sup> � *θ*<sup>0</sup> and *g* � *β*<sup>0</sup> cos *θ*0/||**L**1|| for an inner satellite) (e.g. Laplace, 1799; Mignard, 1981; Tremaine et al., 2009). The above equation has two or four real roots for *θ*1, which are known by *Cassini states*. In general, for satellites we have *I*<sup>0</sup> ∼ 0, and these

For a generic value of *I*0, when *λ*<sup>1</sup> � 1, which is often the case of an outer satellite, the first expression gives the only two real roots of equation (29), one for *θ*<sup>1</sup> � *I*<sup>0</sup> and another for *θ*<sup>1</sup> � *π* − *I*0. On the other hand, when *λ*<sup>1</sup> � 1, which is the case of inner satellites, we will

In turn, the dissipative obliquity variations are computed by substituting equation (25) in (28)

*f*1(*e*1) cos *θ*<sup>1</sup>

Because of the factor *n*1/*ω*<sup>1</sup> in the magnitude of the obliquity variations, for an initial fast rotating satellite, the time-scale for the obliquity evolution will be longer than the time-scale for the rotation rate evolution (Eq.26). As a consequence, it is to be expected that the rotation

, <sup>±</sup> cos−<sup>1</sup>

*ω*1 2*n*<sup>1</sup> −cos *<sup>I</sup>*<sup>0</sup> *λ*1

− *f*2(*e*1)

. (30)

. (31)

<sup>=</sup> *<sup>f</sup>*2(*e*1) *f*1(*e*1) *ω*1 *n*1

2 cos *θ*<sup>1</sup> 1 + cos2 *θ*<sup>1</sup>

− *f*2(*e*1) cos *θ*<sup>1</sup>

**L˙** <sup>1</sup> · (**sˆ***<sup>i</sup>* <sup>−</sup> cos *<sup>θ</sup>i***kˆ** <sup>1</sup>)

*λ*<sup>1</sup> cos *θ*<sup>1</sup> sin *θ*<sup>1</sup> + sin(*θ*<sup>1</sup> − *I*0) = 0 , (29)

, (27)


. (26)

and

$$\begin{aligned} \dot{\mathbf{H}}\_{i} &= K\_{i} n\_{1} \left[ f\_{4}(e\_{1}) \sqrt{1 - e\_{1}^{2}} \frac{\omega\_{i}}{2n\_{1}} (\mathfrak{k}\_{i} - \cos \theta\_{i} \,\hat{\mathbf{k}}\_{1}) \right. \\\\ &- f\_{1}(e\_{1}) \frac{\omega\_{i}}{n\_{1}} \mathfrak{s}\_{i} + f\_{2}(e\_{1}) \hat{\mathbf{k}}\_{1} + \frac{(\mathbf{e}\_{1} \cdot \mathfrak{s}\_{i})(\mathfrak{k} + e\_{1}^{2})}{4(1 - e\_{1}^{2})^{9/2}} \frac{\omega\_{i}}{n\_{1}} \mathfrak{e}\_{1} \right], \end{aligned} \tag{18}$$

where,

$$K\_{l} = \Delta t\_{l} \frac{3k\_{2\_{l}}Gm\_{(1-i)}^{2}R\_{l}^{5}}{a\_{1}^{6}} \,, \tag{19}$$

and

$$f\_1(e) = \frac{1 + 3e^2 + \frac{3}{8}e^4}{(1 - e^2)^{9/2}} \, , \tag{20}$$

$$f\_2(e) = \frac{1 + \frac{15}{2}e^2 + \frac{45}{8}e^4 + \frac{5}{16}e^6}{(1 - e^2)^6} \,, \tag{21}$$

$$f\_3(e) = \frac{1 + \frac{31}{2}e^2 + \frac{255}{8}e^4 + \frac{185}{16}e^6 + \frac{25}{64}e^8}{(1 - e^2)^{15/2}}\,\Big/\tag{22}$$

$$f\_4(e) = \frac{1 + \frac{3}{2}e^2 + \frac{1}{8}e^4}{(1 - e^2)^5} \,, \tag{23}$$

$$f\_5(e) = \frac{1 + \frac{15}{4}e^2 + \frac{15}{8}e^4 + \frac{5}{64}e^6}{(1 - e^2)^{13/2}} \,. \tag{24}$$

The first term in expression (17) corresponds to a permanent tidal deformation, while the second term corresponds to the dissipative contribution. The precession rate of **e**<sup>1</sup> about **kˆ** <sup>1</sup> is usually much faster than the evolution time-scale for the dissipative tidal effects. As a consequence, when the eccentricity is constant over a precession cycle, we can average expression (18) over the argument of the periapsis and get (Correia, 2009):

$$\mathbf{\dot{H}}\_{i} = -K\_{i}n\_{1}\left(f\_{1}(c\_{1})\frac{\mathbf{s}\_{i} + \cos\theta\_{i}}{2}\frac{\mathbf{\hat{k}}\_{1}}{n\_{1}} - f\_{2}(c\_{1})\mathbf{\hat{k}}\_{1}\right). \tag{25}$$

#### **3. Secular evolution**

In the previous section we presented the equations that rule the tidal evolution of a satellite in terms of angular momenta and orbital energy. However, the spin and orbital quantities are better represented by the rotation angles and elliptical elements. The direction cosines (Eq.15) are obtained from the angular momenta vectors, since **sˆ***<sup>i</sup>* <sup>=</sup> **<sup>H</sup>***i*/||**H***i*|| and **kˆ***<sup>i</sup>* <sup>=</sup> **<sup>L</sup>***i*/||**L***i*||, as well as the rotation rate *ω<sup>i</sup>* = **H***<sup>i</sup>* · **sˆ***i*/*Ci*. The eccentricity and the semi-major axis can be obtained from *<sup>e</sup>*<sup>1</sup> <sup>=</sup> ||**e**1||, and *<sup>a</sup>*<sup>1</sup> <sup>=</sup> ||**L**1||2/(*GMm*<sup>2</sup> <sup>1</sup>(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup> <sup>1</sup>)), respectively.

#### **3.1 Spin evolution**

4 Will-be-set-by-IN-TECH

(**sˆ***<sup>i</sup>* <sup>−</sup> cos *<sup>θ</sup><sup>i</sup>* **kˆ** <sup>1</sup>)

<sup>4</sup>(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

(1−*i*)*R*<sup>5</sup> *i*

8 *e*4

<sup>8</sup> *<sup>e</sup>*<sup>4</sup> <sup>+</sup> <sup>5</sup> <sup>16</sup> *<sup>e</sup>*<sup>6</sup>

<sup>8</sup> *<sup>e</sup>*<sup>4</sup> <sup>+</sup> <sup>5</sup> <sup>64</sup> *<sup>e</sup>*<sup>6</sup>

> *ωi n*1

<sup>1</sup>(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

<sup>−</sup> *<sup>f</sup>*2(*e*1)**kˆ** <sup>1</sup>

<sup>1</sup>)), respectively.

<sup>16</sup> *<sup>e</sup>*<sup>6</sup> <sup>+</sup> <sup>25</sup>

<sup>64</sup> *<sup>e</sup>*<sup>8</sup>

<sup>8</sup> *<sup>e</sup>*<sup>4</sup> <sup>+</sup> <sup>185</sup>

<sup>2</sup> *<sup>e</sup>*<sup>2</sup> <sup>+</sup> <sup>1</sup> 8 *e*4

<sup>4</sup> *<sup>e</sup>*<sup>2</sup> <sup>+</sup> <sup>15</sup>

The first term in expression (17) corresponds to a permanent tidal deformation, while the second term corresponds to the dissipative contribution. The precession rate of **e**<sup>1</sup> about **kˆ** <sup>1</sup> is usually much faster than the evolution time-scale for the dissipative tidal effects. As a consequence, when the eccentricity is constant over a precession cycle, we can average

> **sˆ***<sup>i</sup>* + cos *θ<sup>i</sup>* **kˆ** <sup>1</sup> 2

In the previous section we presented the equations that rule the tidal evolution of a satellite in terms of angular momenta and orbital energy. However, the spin and orbital quantities are better represented by the rotation angles and elliptical elements. The direction cosines (Eq.15) are obtained from the angular momenta vectors, since **sˆ***<sup>i</sup>* <sup>=</sup> **<sup>H</sup>***i*/||**H***i*|| and **kˆ***<sup>i</sup>* <sup>=</sup> **<sup>L</sup>***i*/||**L***i*||, as well as the rotation rate *ω<sup>i</sup>* = **H***<sup>i</sup>* · **sˆ***i*/*Ci*. The eccentricity and the semi-major axis can be

1)

*ωi n*1 **e**1 

(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*2)9/2 , (20)

(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*2)<sup>6</sup> , (21)

(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*2)15/2 , (22)

(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*2)<sup>5</sup> , (23)

(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*2)13/2 . (24)

. (25)

, (19)

, (18)

1)9/2

**sˆ***<sup>i</sup>* <sup>+</sup> *<sup>f</sup>*2(*e*1)**kˆ** <sup>1</sup> <sup>+</sup> (**e**<sup>1</sup> · **sˆ***i*)(<sup>6</sup> <sup>+</sup> *<sup>e</sup>*<sup>2</sup>

*a*6 1

and

where,

and

**H˙** *<sup>i</sup>* = *Ki n*<sup>1</sup>

 *f*4(*e*1) <sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup> 1 *ωi* 2*n*<sup>1</sup>

> *ωi n*1

> > *Ki* = Δ*ti*

*<sup>f</sup>*2(*e*) = <sup>1</sup> <sup>+</sup> <sup>15</sup>

*<sup>f</sup>*5(*e*) = <sup>1</sup> <sup>+</sup> <sup>15</sup>

expression (18) over the argument of the periapsis and get (Correia, 2009):

*f*1(*e*1)

**H˙** *<sup>i</sup>* <sup>=</sup> <sup>−</sup>*Ki <sup>n</sup>*<sup>1</sup>

obtained from *<sup>e</sup>*<sup>1</sup> <sup>=</sup> ||**e**1||, and *<sup>a</sup>*<sup>1</sup> <sup>=</sup> ||**L**1||2/(*GMm*<sup>2</sup>

**3. Secular evolution**

*<sup>f</sup>*3(*e*) = <sup>1</sup> <sup>+</sup> <sup>31</sup>

3*k*2*<sup>i</sup> Gm*<sup>2</sup>

*<sup>f</sup>*1(*e*) = <sup>1</sup> <sup>+</sup> <sup>3</sup>*e*<sup>2</sup> <sup>+</sup> <sup>3</sup>

<sup>2</sup> *<sup>e</sup>*<sup>2</sup> <sup>+</sup> <sup>255</sup>

*<sup>f</sup>*4(*e*) = <sup>1</sup> <sup>+</sup> <sup>3</sup>

<sup>2</sup> *<sup>e</sup>*<sup>2</sup> <sup>+</sup> <sup>45</sup>

−*f*1(*e*1)

The variation in the satellite's rotation rate can be obtained from *<sup>ω</sup>*˙ *<sup>i</sup>* <sup>=</sup> **H˙** *<sup>i</sup>* · **sˆ***i*/*Ci* (Eq. 25), giving (Correia & Laskar, 2009):

$$\dot{\omega}\_1 = -\frac{K\_1 n\_1}{\mathbb{C}\_1} \left( f\_1(e\_1) \frac{1 + \cos^2 \theta\_1}{2} \frac{\omega\_1}{n\_1} - f\_2(e\_1) \cos \theta\_1 \right) \,. \tag{26}$$

For a given obliquity and eccentricity, the equilibrium rotation rate, obtained when *ω*˙ <sup>1</sup> = 0, is attained for:

$$\frac{\omega\_1}{m\_1} = \frac{f\_2(e\_1)}{f\_1(e\_1)} \frac{2\cos\theta\_1}{1+\cos^2\theta\_1} \,\mathrm{\,}\tag{27}$$

The obliquity variations can be obtained from equation (15):

$$\frac{d\cos\theta\_{i}}{dt} = \frac{\dot{\mathbf{H}}\_{i} \cdot (\hat{\mathbf{k}}\_{1} - \cos\theta\_{i}\mathbf{\hat{s}}\_{i})}{||\mathbf{H}\_{i}||} + \frac{\dot{\mathbf{L}}\_{1} \cdot (\mathbf{\hat{s}}\_{i} - \cos\theta\_{i}\mathbf{\hat{k}}\_{1})}{||\mathbf{L}\_{1}||} \,. \tag{28}$$

For the conservative motion (Eqs. 9, 11), stable configurations for the spin can be found whenever the vectors (**sˆ**1, **kˆ** 1, **kˆ** 0) or (**sˆ**1, **kˆ** 1, **sˆ**0) are coplanar and precess at the same rate *g* (e.g. Colombo, 1966; Correia et al., 2011; Peale, 1969). The first situation occurs if *γ* � *β*<sup>0</sup> (outer satellite) and the second situation when *γ* � *β*<sup>0</sup> (inner satellite). The equilibrium obliquities can be found from a single relationship (e.g. Ward & Hamilton, 2004):

$$
\lambda\_1 \cos \theta\_1 \sin \theta\_1 + \sin(\theta\_1 - I\_0) = 0 \,, \tag{29}
$$

where *λ*<sup>1</sup> = *β*1/(*C*1*ω*1*g*) is a dimensionless parameter and *I*<sup>0</sup> is the inclination of the orbit of the satellite with respect to the Laplacian plane (*I*<sup>0</sup> � *I* and *g* � *γ* cos *I*/||**L**1|| for an outer satellite, and *I*<sup>0</sup> � *θ*<sup>0</sup> and *g* � *β*<sup>0</sup> cos *θ*0/||**L**1|| for an inner satellite) (e.g. Laplace, 1799; Mignard, 1981; Tremaine et al., 2009). The above equation has two or four real roots for *θ*1, which are known by *Cassini states*. In general, for satellites we have *I*<sup>0</sup> ∼ 0, and these solutions are approximately given by:

$$\tan^{-1}\left(\frac{\sin I\_0}{\cos I\_0 \pm \lambda\_1}\right) \quad \pm \cos^{-1}\left(-\frac{\cos I\_0}{\lambda\_1}\right) \tag{30}$$

For a generic value of *I*0, when *λ*<sup>1</sup> � 1, which is often the case of an outer satellite, the first expression gives the only two real roots of equation (29), one for *θ*<sup>1</sup> � *I*<sup>0</sup> and another for *θ*<sup>1</sup> � *π* − *I*0. On the other hand, when *λ*<sup>1</sup> � 1, which is the case of inner satellites, we will have four real roots approximately given by expressions (30).

In turn, the dissipative obliquity variations are computed by substituting equation (25) in (28) with ||**H**1|| � ||**L**1||, giving:

$$\dot{\theta}\_1 \simeq \frac{K\_1 n\_1}{\mathcal{C}\_1 \omega\_1} \sin \theta\_1 \left( f\_1(e\_1) \cos \theta\_1 \frac{\omega\_1}{2n\_1} - f\_2(e\_1) \right) \,. \tag{31}$$

Because of the factor *n*1/*ω*<sup>1</sup> in the magnitude of the obliquity variations, for an initial fast rotating satellite, the time-scale for the obliquity evolution will be longer than the time-scale for the rotation rate evolution (Eq.26). As a consequence, it is to be expected that the rotation

time [Myr]

1 10 100 1000 10000

Fig. 1. Secular tidal evolution of the Triton-Neptune system. We plot the semi-major axis ratios *<sup>a</sup>*1/*R*<sup>0</sup> and *<sup>a</sup> <sup>f</sup>* /*R*0, the periastron distance *<sup>p</sup>*1/*R*0, the inverse of the eccentricity *<sup>e</sup>*−<sup>1</sup>

its interior is still warm and geologically active considering its distance from the Sun (Schenk & Zahnle, 2007). Its composition also presents some similarities with Pluto (Tsurutani et al.,

These bizarre characteristics lead one to believe that Triton originally orbited the Sun, belonging to the family of Kuiper-belt objects. Most likely during the outward migration of Neptune, the orbits of the two bodies intercepted and capture occurred. This possibility is strongly supported by the fact that Triton's present orbit lies between a group of small inner prograde satellites and a number of exterior irregular satellites both prograde and retrograde. Nereid, with an orbital eccentricity around 0.75, is also believed to have been scattered from a

How exactly the capture occurred is still unknown, but some mechanisms have been proposed: gas drag (McKinnon & Leith, 1995; Pollack et al., 1979), a collision with a pre-existing regular satellite of Neptune (Goldreich et al., 1989), or three-body interactions (Agnor & Hamilton, 2006; Vokrouhlický et al., 2008). All these scenarios require a very close passage to Neptune, and leave the planet in eccentric orbits that must be damped by tides to the present one. Tides are thus the only consensual mechanism acting on Triton's orbit. The tidal distortion of Triton after a few close passages around Neptune, and the consequent dissipation of tidal energy, can account for a substantial reduction in the semi-major axis of its orbit, quickly bringing the planet from an orbit outside Neptune's Hill sphere (∼ 4700 *R*0)

ω1 /n<sup>1</sup>

<sup>1</sup> and

e1 −1

0.1

1990).

Triton's rotation rate ratio *ω*1/*n*1.

regular satellite orbit (McKinnon, 1984).

1

10

100

af/*R*<sup>0</sup>

*p*/*R*<sup>0</sup> <sup>1</sup>

a1/*R*<sup>0</sup>

Secular Evolution of Satellites by Tidal Effect 109

rate reaches its equilibrium value (Eq.27) earlier than the obliquity. Thus, replacing equation (27) in (31), we have:

$$\dot{\theta}\_1 \simeq -\frac{K\_1 n\_1}{\mathcal{C}\_1 \omega\_1} f\_2(e\_1) \frac{\sin \theta\_1}{1 + \cos^2 \theta\_1} \,. \tag{32}$$

We then conclude that the obliquity can only decrease by tidal effect, since ˙ *θ*<sup>1</sup> ≤ 0, and the final obliquity tends to be captured in low obliquity Cassini states.

#### **3.2 Orbital evolution**

The variations in the eccentricity are easily obtained from the Laplace-Runge-Lenz vector (Eq.17):

$$\dot{e}\_{1} = \frac{\dot{\mathbf{e}}\_{1} \cdot \mathbf{e}\_{1}}{e\_{1}} = \sum\_{i} \frac{9\mathbf{K}\_{i}}{m\_{1}a\_{1}^{2}} \left(\frac{11}{18}f\_{4}(e\_{1})\cos\theta\_{i}\frac{\omega\_{i}}{n\_{1}} - f\_{5}(e\_{1})\right)e\_{1} \tag{33}$$

while the semi-major axis variations are obtained from the eccentricity and the norm of the orbital angular momentum:

$$\frac{d\_1}{a\_1} = \frac{2\dot{e}\_1 e\_1}{(1 - e\_1^2)} + \frac{2\dot{\mathbf{L}}\_1 \cdot \mathbf{L}\_1}{||\mathbf{L}\_1||^2} = \sum\_i \frac{2\mathbf{K}\_i}{m\_1 a\_1^2} \left( f\_2(e\_1) \cos \theta\_i \frac{\omega\_i}{n\_1} - f\_3(e\_1) \right) \,. \tag{34}$$

For gaseous planets and rocky satellites we usually have *k*20Δ*t*<sup>0</sup> � *k*21Δ*t*1, and we can retain only terms in *K*1.

The ratio between orbital and spin evolution time-scales is roughly given by *C*1/(*m*1*a*<sup>2</sup> <sup>1</sup>) � 1, meaning that the spin achieves an equilibrium position (**H˙** <sup>1</sup> = 0) much faster than the orbit. Replacing the equilibrium rotation rate (Eq. 27) with *θ*<sup>1</sup> = 0 (for simplicity) in equations (33) and (34), gives:

*<sup>a</sup>*˙1 <sup>=</sup> <sup>−</sup> <sup>7</sup>*K*<sup>1</sup> *m*1*a*<sup>1</sup> *f*6(*e*1)*e* 2 <sup>1</sup> , (35)

$$\dot{e}\_1 = -\frac{7\mathcal{K}\_1}{2m\_1a\_1^2}f\_6(e\_1)(1-e\_1^2)e\_1\,\,\,\,\,\,\,\tag{36}$$

where

$$f\_{\delta}(e) = \frac{1 + \frac{45}{14}e^2 + 8e^4 + \frac{685}{224}e^6 + \frac{255}{448}e^8 + \frac{25}{1792}e^{10}}{(1 + 3e^2 + \frac{3}{8}e^4)(1 - e^2)^{15/2}}.\tag{37}$$

Thus, we always have *a*˙1 ≤ 0 and *e*˙1 ≤ 0, and the final eccentricity is zero. However, from this point onwards, the tidal effects on the planet cannot be neglected (Eq.34), and they govern the future evolution of the satellite's orbit. For *a <sup>f</sup>* < *as* or *θ*<sup>0</sup> ≥ *π*/2, where *a*3 *<sup>s</sup>* = *Gm*0/(*ω*<sup>0</sup> cos *θ*0)2, the semi-major axis continues to decrease until the satellite crashes into the planet, while in the remaining situations it will increase.

#### **4. Application to Triton-Neptune**

Neptune's main satellite, Triton, presents unique features in the Solar System. It is the only moon-sized body in a retrograde inclined orbit and the images taken by the Voyager 2 spacecraft in 1989 revealed an extremely young surface with few impact craters (e.g. Cruikshank, 1995). This satellite should have remained molten until about 1 Gyr ago and 6 Will-be-set-by-IN-TECH

rate reaches its equilibrium value (Eq.27) earlier than the obliquity. Thus, replacing equation

The variations in the eccentricity are easily obtained from the Laplace-Runge-Lenz vector

while the semi-major axis variations are obtained from the eccentricity and the norm of the

2*Ki m*1*a*<sup>2</sup> 1

For gaseous planets and rocky satellites we usually have *k*20Δ*t*<sup>0</sup> � *k*21Δ*t*1, and we can retain

meaning that the spin achieves an equilibrium position (**H˙** <sup>1</sup> = 0) much faster than the orbit. Replacing the equilibrium rotation rate (Eq. 27) with *θ*<sup>1</sup> = 0 (for simplicity) in equations (33)

<sup>18</sup> *<sup>f</sup>*4(*e*1) cos *<sup>θ</sup><sup>i</sup>*

*f*6(*e*1)*e* 2

2

<sup>448</sup> *<sup>e</sup>*<sup>8</sup> <sup>+</sup> <sup>25</sup>

<sup>1792</sup> *<sup>e</sup>*<sup>10</sup>

<sup>8</sup> *<sup>e</sup>*4)(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*2)15/2 . (37)

*f*6(*e*1)(1 − *e*

<sup>224</sup> *<sup>e</sup>*<sup>6</sup> <sup>+</sup> <sup>255</sup>

11

*i*

The ratio between orbital and spin evolution time-scales is roughly given by *C*1/(*m*1*a*<sup>2</sup>

*<sup>a</sup>*˙1 <sup>=</sup> <sup>−</sup> <sup>7</sup>*K*<sup>1</sup> *m*1*a*<sup>1</sup>

<sup>14</sup> *<sup>e</sup>*<sup>2</sup> <sup>+</sup> <sup>8</sup>*e*<sup>4</sup> <sup>+</sup> <sup>685</sup>

(1 + 3*e*<sup>2</sup> + <sup>3</sup>

Thus, we always have *a*˙1 ≤ 0 and *e*˙1 ≤ 0, and the final eccentricity is zero. However, from this point onwards, the tidal effects on the planet cannot be neglected (Eq.34), and they govern the future evolution of the satellite's orbit. For *a <sup>f</sup>* < *as* or *θ*<sup>0</sup> ≥ *π*/2, where

*<sup>s</sup>* = *Gm*0/(*ω*<sup>0</sup> cos *θ*0)2, the semi-major axis continues to decrease until the satellite crashes

Neptune's main satellite, Triton, presents unique features in the Solar System. It is the only moon-sized body in a retrograde inclined orbit and the images taken by the Voyager 2 spacecraft in 1989 revealed an extremely young surface with few impact craters (e.g. Cruikshank, 1995). This satellite should have remained molten until about 1 Gyr ago and

*<sup>e</sup>*˙1 <sup>=</sup> <sup>−</sup> <sup>7</sup>*K*<sup>1</sup> 2*m*1*a*<sup>2</sup> 1

*<sup>f</sup>*6(*e*) = <sup>1</sup> <sup>+</sup> <sup>45</sup>

into the planet, while in the remaining situations it will increase.

**4. Application to Triton-Neptune**

*<sup>f</sup>*2(*e*1) sin *<sup>θ</sup>*<sup>1</sup>

1 + cos2 *θ*<sup>1</sup>

*ωi n*1

*f*2(*e*1) cos *θ<sup>i</sup>*

− *f*5(*e*1)

*ωi n*1 − *f*3(*e*1)

<sup>1</sup> , (35)

<sup>1</sup>)*e*<sup>1</sup> , (36)

. (32)

*θ*<sup>1</sup> ≤ 0, and the

*e*<sup>1</sup> , (33)

. (34)

<sup>1</sup>) � 1,

˙

final obliquity tends to be captured in low obliquity Cassini states.

= ∑ *i*

<sup>2</sup> **L˙** <sup>1</sup> · **<sup>L</sup>**<sup>1</sup> ||**L**1||<sup>2</sup> <sup>=</sup> <sup>∑</sup>

*<sup>e</sup>*˙1 <sup>=</sup> **e˙**<sup>1</sup> · **<sup>e</sup>**<sup>1</sup> *e*1

<sup>=</sup> <sup>2</sup>*e*˙1*e*<sup>1</sup> (<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup> 1) + *<sup>θ</sup>*<sup>1</sup> � − *<sup>K</sup>*1*n*<sup>1</sup> *C*1*ω*<sup>1</sup>

We then conclude that the obliquity can only decrease by tidal effect, since ˙

9*Ki m*1*a*<sup>2</sup> 1

(27) in (31), we have:

**3.2 Orbital evolution**

orbital angular momentum:

*a*˙1 *a*1

only terms in *K*1.

and (34), gives:

where

*a*3

(Eq.17):

Fig. 1. Secular tidal evolution of the Triton-Neptune system. We plot the semi-major axis ratios *<sup>a</sup>*1/*R*<sup>0</sup> and *<sup>a</sup> <sup>f</sup>* /*R*0, the periastron distance *<sup>p</sup>*1/*R*0, the inverse of the eccentricity *<sup>e</sup>*−<sup>1</sup> <sup>1</sup> and Triton's rotation rate ratio *ω*1/*n*1.

its interior is still warm and geologically active considering its distance from the Sun (Schenk & Zahnle, 2007). Its composition also presents some similarities with Pluto (Tsurutani et al., 1990).

These bizarre characteristics lead one to believe that Triton originally orbited the Sun, belonging to the family of Kuiper-belt objects. Most likely during the outward migration of Neptune, the orbits of the two bodies intercepted and capture occurred. This possibility is strongly supported by the fact that Triton's present orbit lies between a group of small inner prograde satellites and a number of exterior irregular satellites both prograde and retrograde. Nereid, with an orbital eccentricity around 0.75, is also believed to have been scattered from a regular satellite orbit (McKinnon, 1984).

How exactly the capture occurred is still unknown, but some mechanisms have been proposed: gas drag (McKinnon & Leith, 1995; Pollack et al., 1979), a collision with a pre-existing regular satellite of Neptune (Goldreich et al., 1989), or three-body interactions (Agnor & Hamilton, 2006; Vokrouhlický et al., 2008). All these scenarios require a very close passage to Neptune, and leave the planet in eccentric orbits that must be damped by tides to the present one. Tides are thus the only consensual mechanism acting on Triton's orbit. The tidal distortion of Triton after a few close passages around Neptune, and the consequent dissipation of tidal energy, can account for a substantial reduction in the semi-major axis of its orbit, quickly bringing the planet from an orbit outside Neptune's Hill sphere (∼ 4700 *R*0)

**6. References**

Agnor, C. B. & Hamilton, D. P. (2006). Neptune's capture of its moon Triton in a binary-planet

Secular Evolution of Satellites by Tidal Effect 111

Anderson, J. D., Schubert, G., Jacobson, R. A., Lau, E. L., Moore, W. B. & Sjogren, W. L. (1998).

Canup, R. M. & Asphaug, E. (2001). Origin of the Moon in a giant impact near the end of the

Chyba, C. F., Jankowski, D. G. & Nicholson, P. D. (1989). Tidal evolution in the Neptune-Triton

Correia, A. C. M. (2006). The core-mantle friction effect on the secular spin evolution of

Correia, A. C. M. (2009). Secular Evolution of a Satellite by Tidal Effect: Application to Triton,

Correia, A. C. M. & Laskar, J. (2009). Mercury's capture into the 3/2 spin-orbit resonance

Correia, A. C. M., Laskar, J., Farago, F. & Boué, G. (2011). Tidal evolution of hierarchical and inclined systems, *Celestial Mechanics and Dynamical Astronomy* 111: 105–130.

Farago, F. & Laskar, J. (2010). High-inclination orbits in the secular quadrupolar three-body

Goldreich, P., Murray, N., Longaretti, P. Y. & Banfield, D. (1989). Neptune's story, *Science*

Jacobson, R. A. (2009). The Orbits of the Neptunian Satellites and the Orientation of the Pole

Kaula, W. M. (1964). Tidal dissipation by solid friction and the resulting orbital evolution, *Rev.*

McCord, T. B. (1966). Dynamical evolution of the Neptunian system, *Astron. J.* 71: 585–590.

McKinnon, W. B. & Leith, A. C. (1995). Gas drag and the orbital evolution of a captured Triton.,

Mignard, F. (1979). The evolution of the lunar orbit revisited. I, *Moon and Planets* 20: 301–315. Mignard, F. (1981). Evolution of the Martian satellites, *Mon. Not. R. Astron. Soc.* 194: 365–379. Néron de Surgy, O. & Laskar, J. (1997). On the long term evolution of the spin of the Earth.,

Pollack, J. B., Burns, J. A. & Tauber, M. E. (1979). Gas drag in primordial circumplanetary

Schenk, P. M. & Zahnle, K. (2007). On the negligible surface age of Triton, *Icarus* 192: 135–149.

Stern, S. A. & McKinnon, W. B. (2000). Triton's Surface Age and Impactor Population Revisited

in Light of Kuiper Belt Fluxes: Evidence for Small Kuiper Belt Objects and Recent

envelopes - A mechanism for satellite capture, *Icarus* 37: 587–611.

Smart, W. M. (1953). *Celestial Mechanics.*, London, New York, Longmans, Green.

McKinnon, W. B. (1984). On the origin of Triton and Pluto, *Nature* 311: 355–358.

Boué, G. & Laskar, J. (2006). Precession of a planet with a satellite, *Icarus* 185: 312–330. Canup, R. M. (2005). A Giant Impact Origin of Pluto-Charon, *Science* 307: 546–550.

Colombo, G. (1966). Cassini's second and third laws, *Astron. J.* 71: 891–896.

including the effect of core-mantle friction, *Icarus* 201: 1–11.

Cruikshank, D. P. (1995). *Neptune and Triton*, University of Arizona Press.

problem, *Mon. Not. R. Astron. Soc.* 401: 1189–1198.

Laplace, P. S. (1799). *Traité de Mécanique céleste*, Paris: Gauthier-Villars.

Peale, S. J. (1969). Generalized Cassini's Laws, *Astron. J.* 74: 483–489.

Geological Activity, *Astron. J.* 119: 945–952.

of Neptune, *Astron. J.* 137: 4322–4329.

terrestrial planets, *Earth Planet. Sci. Lett.* 252: 398–412.

Europa's Differentiated Internal Structure: Inferences from Four Galileo Encounters,

gravitational encounter, *Nature* 441: 192–194.

Earth's formation, *Nature* 412: 708–712.

system, *Astron. Astrophys.* 219: L23–L26.

*Science* 281: 2019–2022.

*Astrophys. J.* 704: L1–L4.

245: 500–504.

*Geophys.* 2: 661–685.

*Icarus* 118: 392–413.

*Astron. Astrophys.* 318: 975–989.

to a bounded orbit. Therefore, it cannot be ruled out that Triton was simply captured by tidal interactions with Neptune after a close encounter in an almost parabolic orbit (McCord, 1966).

Here, we simulate the tidal evolution of the Triton-Neptune system using the complete model described in Sect.2. Triton is started in a very elliptical orbit with *e*<sup>1</sup> = 0.99 and a semi-major axis of *a*<sup>1</sup> = 40 *R*0, corresponding to a final equilibrium *a <sup>f</sup>* � 14.4 *R*0, close to the present position of 14.33 *R*0. These specific values give a closest approach at a periapse of 8 *R*0.

For the radius of the bodies we use *R*<sup>0</sup> = 24 764 km and *R*<sup>1</sup> = 1 353 km (Thomas, 2000), while for the masses, the *J*<sup>2</sup> of Neptune and the remaining orbital and spin parameters we take the present values as determined by Jacobson (2009). For Triton we adopt *<sup>J</sup>*<sup>2</sup> <sup>=</sup> 4.38 <sup>×</sup> <sup>10</sup>−4, the value measured for Europa (Anderson et al., 1998), and *C*<sup>22</sup> = 0, since our model does not take into account spin-orbit resonances. This choice is justified because Triton's observed topography never varies beyond a kilometer (Thomas, 2000). In addition, Triton should have undergone frequent collisions either with other satellites of Neptune, or with external Kuiper-belt objects, and any capture in a spin-orbit resonance different from the synchronous one, may not last for a long time (Stern & McKinnon, 2000). For tidal dissipation we adopt the same parameters as previous studies, that is, *k*<sup>20</sup> = 0.407 and *Q*<sup>0</sup> = 9000 (Zhang & Hamilton, 2008), and *k*<sup>21</sup> = 0.1 and *Q*<sup>1</sup> = 100 (Chyba et al., 1989; Goldreich et al., 1989), where *Q*−<sup>1</sup> *<sup>i</sup>* = *ωi*Δ*ti*.

As for the orbit, the initial spin of Triton is unknown. We tested several possibilities, but tides acting on the spin always drive it in the same way: the rotation rate quickly evolves into the equilibrium value given by equation (27), while the obliquity is trapped in a Cassini state. In our standard simulation (Fig.1) we start Triton with a rotation period of 24 h. The semi-major axis and the eccentricity always decrease, as predicted by equations (35) and (36), and the quantity *<sup>a</sup> <sup>f</sup>* <sup>=</sup> *<sup>a</sup>*1(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup> <sup>1</sup>) is preserved during the first stages of the evolution, with the reduction observed being caused by tides on Neptune. The eccentricity is very high during the first evolutionary stages, but it decreases rapidly as the satellite approaches its present orbit. Finally, the rotation of the satellite decreases as the satellite orbit shrinks into Neptune and ultimately stabilizes in the synchronous resonance, the presently observed configuration.

#### **5. Conclusion**

The numerical results presented here are very similar to those shown in Correia (2009) for *a*<sup>1</sup> ≤ 40 *R*0. For these values of the semi-major axis, Triton can still be considered as an "inner satellite", that is, the inclination with respect to the equatorial plane of Neptune, *θ*0, is approximately constant. However, for higher values of the semi-major axis, we observe exchanges between the inclination and the eccentricity of Triton. Since the eccentricity is no longer constant, it is not possible to average over the argument of the periastron as in Correia (2009). Therefore, the results with the non-averaged model presented here will show some differences. In particular, the perturbation on Triton's orbit will cause the eccentricity to vary around the mean value, allowing the periapse to attain lower values. As a consequence, tidal effects will be stronger for close encounters with Neptune, and the migration of Triton may occur in much faster time-scales. In a future study we will analyze in detail the evolution of Triton's orbit for *a*<sup>1</sup> > 40 *R*0.

Our study should also apply to the Moon, Charon and the satellites of Mars, although in this case we need to take into account the quadropole moment of inertia *C*<sup>22</sup> �= 0 (Correia, 2006).

#### **6. References**

8 Will-be-set-by-IN-TECH

to a bounded orbit. Therefore, it cannot be ruled out that Triton was simply captured by tidal interactions with Neptune after a close encounter in an almost parabolic orbit (McCord, 1966). Here, we simulate the tidal evolution of the Triton-Neptune system using the complete model described in Sect.2. Triton is started in a very elliptical orbit with *e*<sup>1</sup> = 0.99 and a semi-major axis of *a*<sup>1</sup> = 40 *R*0, corresponding to a final equilibrium *a <sup>f</sup>* � 14.4 *R*0, close to the present position of 14.33 *R*0. These specific values give a closest approach at a periapse of 8 *R*0.

For the radius of the bodies we use *R*<sup>0</sup> = 24 764 km and *R*<sup>1</sup> = 1 353 km (Thomas, 2000), while for the masses, the *J*<sup>2</sup> of Neptune and the remaining orbital and spin parameters we take the present values as determined by Jacobson (2009). For Triton we adopt *<sup>J</sup>*<sup>2</sup> <sup>=</sup> 4.38 <sup>×</sup> <sup>10</sup>−4, the value measured for Europa (Anderson et al., 1998), and *C*<sup>22</sup> = 0, since our model does not take into account spin-orbit resonances. This choice is justified because Triton's observed topography never varies beyond a kilometer (Thomas, 2000). In addition, Triton should have undergone frequent collisions either with other satellites of Neptune, or with external Kuiper-belt objects, and any capture in a spin-orbit resonance different from the synchronous one, may not last for a long time (Stern & McKinnon, 2000). For tidal dissipation we adopt the same parameters as previous studies, that is, *k*<sup>20</sup> = 0.407 and *Q*<sup>0</sup> = 9000 (Zhang & Hamilton, 2008), and *k*<sup>21</sup> = 0.1 and *Q*<sup>1</sup> = 100 (Chyba et al., 1989; Goldreich et al., 1989),

As for the orbit, the initial spin of Triton is unknown. We tested several possibilities, but tides acting on the spin always drive it in the same way: the rotation rate quickly evolves into the equilibrium value given by equation (27), while the obliquity is trapped in a Cassini state. In our standard simulation (Fig.1) we start Triton with a rotation period of 24 h. The semi-major axis and the eccentricity always decrease, as predicted by equations (35) and (36),

reduction observed being caused by tides on Neptune. The eccentricity is very high during the first evolutionary stages, but it decreases rapidly as the satellite approaches its present orbit. Finally, the rotation of the satellite decreases as the satellite orbit shrinks into Neptune and ultimately stabilizes in the synchronous resonance, the presently observed configuration.

The numerical results presented here are very similar to those shown in Correia (2009) for *a*<sup>1</sup> ≤ 40 *R*0. For these values of the semi-major axis, Triton can still be considered as an "inner satellite", that is, the inclination with respect to the equatorial plane of Neptune, *θ*0, is approximately constant. However, for higher values of the semi-major axis, we observe exchanges between the inclination and the eccentricity of Triton. Since the eccentricity is no longer constant, it is not possible to average over the argument of the periastron as in Correia (2009). Therefore, the results with the non-averaged model presented here will show some differences. In particular, the perturbation on Triton's orbit will cause the eccentricity to vary around the mean value, allowing the periapse to attain lower values. As a consequence, tidal effects will be stronger for close encounters with Neptune, and the migration of Triton may occur in much faster time-scales. In a future study we will analyze in detail the evolution of

Our study should also apply to the Moon, Charon and the satellites of Mars, although in this case we need to take into account the quadropole moment of inertia *C*<sup>22</sup> �= 0 (Correia, 2006).

<sup>1</sup>) is preserved during the first stages of the evolution, with the

where *Q*−<sup>1</sup>

**5. Conclusion**

*<sup>i</sup>* = *ωi*Δ*ti*.

and the quantity *<sup>a</sup> <sup>f</sup>* <sup>=</sup> *<sup>a</sup>*1(<sup>1</sup> <sup>−</sup> *<sup>e</sup>*<sup>2</sup>

Triton's orbit for *a*<sup>1</sup> > 40 *R*0.


**Part 4** 

**Cosmology and CMB Physics** 

Thomas, P. C. (2000). NOTE: The Shape of Triton from Limb Profiles, *Icarus* 148: 587–588.

