A&A 429, 731-738 (2005)
DOI: 10.1051/0004-6361:20041173

On the dynamics of some resonances of Phobos in the future

T. Yokoyama - M. R. Mana - C. do Nascimento - M. T. Santos - N. Callegari Jr

Universidade Estadual Paulista, IGCE-DEMAC Caixa Postal 178,CEP 13.500-970 Rio Claro (São Paulo), Brasil

Received 28 April 2004 / Accepted 8 September 2004

Abstract
A probable capture of Phobos into an interesting resonance was presented in our previous work. With a simple model, considering Mars in a Keplerian and circular orbit, it was shown that once captured in the resonance, the inclination of the satellite reaches very high values. Here, the integrations are extended to much longer times and escape situations are analyzed. These escapes are due to the interaction of new additional resonances, which appear as the inclination starts to increase reaching some specific values. Compared to classical capture in mean motion resonances, we see some interesting differences in this problem. We also include the effect of Mars' eccentricity in the process of the capture. The role played by this eccentricity becomes important, particularly when Phobos encounters a double resonance at $a \approx 2.619 R_{\rm M}$ . Planetary perturbations acting on Mars and variation of its equator are also included. In general, some possible scenarios of the future of Phobos are presented.

Key words: planets and satellites: general - solar system: general

1 Introduction

In Yokoyama (2002), hereafter PYO, we showed that in the future, Phobos will face some interesting resonance because its semi major axis a is continuously decreasing due to the action of tides. This means that the frequencies associated with the longitude of the node ( $\Omega$ ) and the argument of the pericenter (w) are also varying, so that in the neighbourhood of a specific value of the semi major axis, these frequencies can satisfy resonant relations like $k_1 \dot{\Omega}+k_2 \dot{w}+k_3 n_\odot \approx 0$ , where $n_\odot$ is the mean motion of the Sun, k_i are integers, and dots mean time derivatives. This kind of commensurability is rather different to usual mean motion resonances, since only one mean longitude (the Sun) is involved. Moreover, this mean motion is not affected by the tidal effects. The current semimajor axis (a), eccentricity (e) and inclination (I) of Phobos are: $a \approx 2.76 R_{\rm M}$ , $e \approx 0.015$ and $I \approx 1^{\circ}$ , where inclination is referred to Mars' equator and $R_{\rm M}$ is the equatorial radius of Mars. In our previous work we showed that Phobos can be captured in the $\dot{\Omega}+ 2n_\odot$ resonance in the neighbourhood of $a=2.149 R_{\rm M}$ . Even in the absence of tides, i.e., for the ordinary restricted three body problem, the effect of this resonance is clearly visible. Indeed, taking $a \approx 2.149 R_{\rm M}$ and integrating the Cartesian Eqs. (1) as given in PYO (Mars is in a Keplerian orbit), we get Fig. 1. The obliquity of the equator is artificially fixed to $\epsilon =60^{\circ }$ (Panels A, B). For this high value of $\epsilon$ , the inclination varies up to almost $5^{\circ}$ while for $\epsilon =28^{\circ }$ (Panels C, D), the maximum inclination is $\approx$ $2.73^{\circ}$ . If however the semi major axis a is slightly changed to a=2.145, the angle $\Phi =2 \lambda _{\odot }+ \Omega -3\Omega _{\odot }$ does not librate any more and the variation of the inclination is negligible.

$\begin{figure} \par\resizebox{\hsize}{!}{\includegraphics[clip]{1173fig1.eps}} % \end{figure}$

Figure 1: Left panels: variation of the inclination for resonant semi major axes. Right panels: libration of the resonant argument: $\Phi =2 \lambda _{\odot }+ \Omega -3\Omega _{\odot }$ . Initial conditions for the satellite: Top: $a=2.1492R_{\rm M}$ , e=0.015, $l=w=0^{\circ }$ , $\Omega =160^{\circ }$ , $\epsilon =60^{\circ }$ , Bottom: a=2.1498, e=0.015, $l=w=0^{\circ }$ , $\Omega =160^{\circ }$ , $\epsilon =28^{\circ }$ . Mars is assumed in a Keplerian orbit with initial conditions: a=1.5236 AU, e=0.0933, $w=0^{\circ }$ , $\Omega =0^{\circ }$ , $\lambda =0^{\circ }$ , $I=\epsilon$ .

Open with DEXTER

Now, let us suppose a case when the semi major axis of Phobos decreases secularly due to tides and passes through the resonant value ( $a\approx 2.149$ ). According to Szeto (1983) the time rate variation of Phobos' semimajor axis is of the order of 10^-7 cm s^-1. The current obliquity of Mars is $\epsilon\approx 25.19^{\circ}$ , but this value is chaotically varying (Ward 1979; Touma & Wisdom 1993; Laskar & Robutel 1993) and over a long time scale it can sweep a large interval. In particular, in a frequency analysis over 45 Myr, Laskar & Robutel (1993) found a chaotic zone of $\epsilon$ which can range from $(0^{\circ}$ to $60^{\circ})$ . If at the moment when $a\approx 2.149$ , the current value of $\epsilon$ is larger than $\approx 19^{\circ}$ , then capture usually occurs and inclination grows secularly to high values, much greater than the present value. As mentioned in PYO, the reason for this large increase is due to the continuous outward displacement of the libration center (Fig. 2 in PYO) as Phobos' semi major axis decreases. In PYO we always considered Mars in a circular orbit and integrations were not long enough, so that no escape or possible rupture of the resonance were presented. This time we want to explore some escape situations. This problem was studied by many authors: Beaugé & Ferraz-Mello (1994), Gomes (1995a,b, 1997a,b, 1998), Liou & Zook (1997), Hamilton (1994), Jancart et al. (2003), etc. In particular, Gomes has presented a deep study of the phenomenon of resonance lock. He also deduced some interesting formulae that govern the evolution of the eccentricity, inclination and semi major axis not only in the restricted three body problem but in the mutual case as well. Gomes considers different kind of dissipative forces; however in all cases the resonance is mainly due to the combination of two mean longitudes. Therefore, for the restricted three body problem, in general, once a capture takes place, the semi major axis of the particle must stop its secular variation and remains librating around some resonant value. Depending on the dissipative force, the amplitude of the libration of the semi major can either decrease to a stationary value (stable lock) or increase beyond some tolerable limit (unstable lock), so that the escape is unavoidable and capture is only temporary (Gomes 1995b). A noteworthy point in the present problem is the resonant angle: instead of two mean motions, the resonance comes from Phobos' longitude of node and the mean motion of the Sun. Therefore, in opposition to mean motion resonances, in our problem, the satellite remains captured, even if its semi major axis is decreasing (Sect. 3). In PYO we list other similar resonances, e.g., the evection: $\dot{\varpi}-n_{\odot}$ . Although in our problem we are mostly interested in the capture of inclination, it is instructive to comment briefly on some points related to the dynamics of the evection resonance. This commensurability might have been important in the past of the Earth-Moon system (Touma & Wisdom 1998). Capture in this case was possible only for the increasing semi major axis of the Moon. The duration of this regime is interesting, since unlike our problem, the increase of the eccentricity generates some consequences for the evolution of the semi major axis: in the Touma & Wisdom problem, the eccentricity grows rapidly, and after some time, its value is high enough that the Moon's semi major axis stops increasing. This occurs because, when the orbit becomes very eccentric, at the pericenter, the angular motion of the Moon becomes faster than the angular rotation of the Earth. This is not the situation in our problem, that is, the capture in inclination does not dampen the decay of Phobos' semi major axis. In the Earth-Moon problem the dynamics of the evection is very complex, but as pointed out by Touma & Wisdom (1998), capture in this resonance is always temporary (for tidal dissipative forces). Compared to problems involving two mean motion resonance, or evection resonances, we show that the dynamics of our resonance presents some clear differences in the evolution during the capture.

This work is organized in the following way: Sect. 2 shows some situations of capture and escape. To explain these escapes and to put the problem in a workable form, we still keep Mars' equator as the reference plane. The dynamics of the problem when Mars' eccentricity is taken into account is described in Sect. 3. In Sect. 4 we include the planetary perturbations acting on Mars and also the precessional equations of Mars' equator. Section 5 is devoted to conclusions.

2 Escape from the 2n $_\odot+ \dot{\Omega}$ resonance

As mentioned before, many authors have studied the phenomenon of capture and escape in resonances involving mean motion commensurability. Very roughly speaking, due to the dissipative force, the semimajor axis of the perturbed body varies secularly so that at some moment a particular orbital resonance involving the mean motions can be satisfied. Then, provided that certain conditions are fulfilled (see for instance, Henrard 1982; Peale 1989; Malhotra 1990; Gomes 1995a, 1997a, etc.) a capture takes place. Usually, in this case, the ratio of semi major axes, or the semi major axis of the perturbed body (restricted three body problem), stops varying secularly. Also, under the capture, a specific angle connected to the current resonance starts to librate, whereas the eccentricity or inclination begins to vary secularly. Beaugé & Ferraz-Mello (1994) showed that in the case of drag forces, the eccentricity varies secularly during some time, but after that, it has the tendency to stabilize around a specific value which they called universal eccentricity. In spite of this stabilization, the capture may not be permanent, that is, if the temporal variation of the libration angle increases in amplitude, an escape from the resonance will appear and this rupture is essentially due to the unstable character of the resonance lock. Therefore, according to Gomes (1995b), a resonance lock is stable or unstable and this stability depends on the particular dissipative force involved. A detailed study yielding analytical criteria about the permanency in resonance trapping is given in his paper.

Similarly to universal eccentricity, the concept of universal inclination was introduced by Gomes (1997b) considering only one specific resonant librating angle.

However, it seems that the existence of such universal values is not completely clear, at least for arbitrary dissipative forces or when it is not possible to define one specific resonant angle in the libration regime. Actually, Gomes (1997b) says: although the concept of a universal inclination is presented here, I could not get any example where this feature showed in an unequivocal way, with the inclination coming to a value equal to the universal inclination.

In this section, we discuss the escapes in our problem. Even though an analytical proof is not given, we show through numerical experiments that all the escapes occur, without any previous stabilization around some possible universal value. We put the problem into the simplest workable form, so that we still keep a simplified model taking Mars equator as the reference plane. The Sun is in a Keplerian circular orbit, with its current inclination with respect to Mars' equator. Therefore, we consider Phobos disturbed by the Sun and the oblateness of Mars. The coordinate system is fixed at the center of the planet. Besides the mean motion of the Sun ( $n_\odot$ ), the principal frequencies are:

$\displaystyle \dot{w}$	$\textstyle \approx$	$\displaystyle {3 nJ_2R^2_{\rm M} \over 4a^2(1-e^2)^2}(5 \cos^2 I-1)$
$\displaystyle \dot{\Omega}$	$\textstyle \approx$	$\displaystyle {-3 nJ_2R^2_{\rm M} \cos I \over 2a^2(1-e^2)^2}$

where ( $a,e,I,n,w,\Omega$ ) are respectively: semi major axis, eccentricity, inclination, mean motion, argument of the pericenter and longitude of the node. J₂ is the oblateness coefficient, k², is the gravitational constant and $R_{\rm M}$ is the equatorial radius of the planet. The complete disturbing function of the problem is given in Eqs. (4) to (5) in PYO. We integrate the averaged equations of the motion and the tidal effect on the Phobos orbit is computed through the simplest model for the semi major axis:

$\displaystyle {{\rm d}a\over {\rm d}t}= -3 k_2 n a {m\over M} {\left (R_{\rm M} \over a \right)}^5 {1\over Q_{2200}}$

(1)

$\begin{figure} \par\resizebox{8.2cm}{!}{\includegraphics[clip]{1173fig2.eps}} \end{figure}$

Figure 2: Left panels: evolution of the resonant angle ( top), eccentricity ( center), inclination ( bottom) taking the Radau integrator with LL=12. Initial conditions: $a=2.169 R_{\rm M}$ , e=0.01, $I=1^{\circ }$ , $w=0^{\circ }$ , $\Omega =90^{\circ }$ , $\epsilon =25.19^{\circ }$ .Right panels: similar but with $a=2.165 R_{\rm M}$ , e=0.02, $I=1.35^{\circ }$ , $w=0^{\circ }$ , $\Omega =90^{\circ }$ , $\epsilon =25.19^{\circ }$ . As the semi major axis and inclination evolve, new resonances can appear. Whenever a new resonance is encountered the dynamics of the libration of $\Phi$ is clearly affected and either an escape or a different mode of libration occurs. The passage through the new resonances provokes some jumps in the eccentricity.

Open with DEXTER

where k₂ is the second order Love number of Mars, while Q₂₂₀₀ is the dissipation function (Szeto 1983; Yoder 1981).

Figure 2 shows two different capture/escapes in $2n_\odot+ \dot{\Omega}$ resonance. The panels on the left correspond to the variation of the resonant angle, eccentricity and inclination with respect to the semi major axis of a single orbit starting with $a=2.169 R_{\rm M}$ , e=0.01, $I=1^{\circ }$ , $w=0^{\circ }$ , $\Omega =90^{\circ }$ , $\epsilon =25.19^{\circ }$ . Similarly, on the right, the initial conditions are: $a=2.165 R_{\rm M}$ , e=0.02, $I=1.35^{\circ }$ , $w=0^{\circ }$ , $\Omega =90^{\circ }$ , $\epsilon =25.19^{\circ }$ . In these integrations, the precision parameter LL in the Radau (Everhart 1985) integrator is fixed to be LL=12. Panel A shows the evolution of the resonant angle $\Phi =2 \lambda _{\odot }+ \Omega -3\Omega _{\odot }$ . Once captured (at $a \approx 2.149 R_{\rm M}$ ), libration takes place and lasts until $a \approx 2.04 R_{\rm M}$ , when we have escape. In panel B, Phobos is captured in the same way as before, but it does not escape when $a \approx 2.04 R_{\rm M}$ . Instead, $\Phi$ changes its libration mode, and remains librating in this new mode until $a\approx 1.933R_{\rm M}$ . Again, as before, it acquires a new libration mode and it continues until $a\approx 1.82 R_{\rm M}$ when it escapes definitively. Connected to each of these mentioned semi major axes, the eccentricity undergoes some clear jumps, although the inclination seems to ignore such effects. Starting from different initial conditions, the qualitative results are similar, i.e., at some specific values of the semi major axis, $\Phi$ changes its libration mode. In some cases, this change causes the escape from resonance, so that inclination stops its secular variation. To understand these dynamics, we simply calculate the values of the semi major axis where different resonances $k_1\dot{\Omega}+k_2 \dot{w}+k_3 n_\odot$ can occur. Of course, in the beginning, for the current I of Phobos and $a \leq 2.2 R_{\rm M}$ , only $2n_\odot+ \dot{\Omega}$ is possible. Then capture takes place at $a \approx 2.149 R_{\rm M}$ , however, after that, since I is increasing, we note that other resonant combination can be fulfilled. Considering each of the arguments of the cosines given in PYO (Eq. (5)), the following additional resonances can occur together with $\dot{\Omega}+2n_\odot\approx 0$ :

$\begin{displaymath}\left. \begin{array}{ll} \dot{\Omega} +\dot{w}-n_\odot & \app... ...y}\right \} \mbox{ at $a=2.04 R_{\rm M}$ , $I=33.02^{\circ}$ } \end{displaymath}$

$\begin{displaymath}\left. \begin{array}{ll} \dot{\Omega} +\dot{w} & \approx 0 \... ... \right \} \mbox{ at $a=1.93 R_{\rm M}$ , $I=46.38^{\circ}$ } \end{displaymath}$

$\begin{displaymath}\left. \begin{array}{ll} \dot{w}-n_\odot &\approx 0 \nonumbe... ... \right \} \mbox{ at $a=1.82 R_{\rm M}$ , $I=56.06^{\circ}$ } \end{displaymath}$

$\begin{displaymath}\left. \begin{array}{ll} 2\dot{w}+2n_\odot+\dot{\Omega} & \a... ...\right \} \mbox{ at $ a=1.71 R_{\rm M}$ , $I=63.43^{\circ}$ } \end{displaymath}$

$\begin{displaymath}\left. \begin{array}{ll} \dot{w}-n_\odot-\dot{\Omega} &\appr... ...\right \} \mbox{ at $ a=1.60 R_{\rm M}$ , $I=69.02^{\circ}$ } \end{displaymath}$

$\begin{displaymath}\left. \begin{array}{ll} \dot{w}-\dot{\Omega}& \approx 0 \non... ... \right \} \mbox{ at $a=1.51 R_{\rm M}$ , $I=73.15^{\circ}$ } \end{displaymath}$

$\begin{displaymath}\left. \begin{array}{ll} \dot{w}+n_\odot-\dot{\Omega}&\approx... ...\right \} \mbox { at $a=1.43 R_{\rm M}$ , $I=76.23^{\circ}$ }. \end{displaymath}$

These values were obtained taking e=0.015, since usually the variation of the eccentricity is very small (see Fig. 2) compared to the variation of the inclination. At these solutions, two or more resonances occur simultaneously. Actually the problem becomes a non-autonomous system with two degree of freedom. Then, the adiabatic invariant theory should break down in the neighbourhood of the double resonance when they are of the same importance. Some numerical experiments show that in this situation, the motion can become very unstable. More precisely, in Sect. 4 we show a typical situation that illustrates the existence of two different separatrices in the case of two simultaneous resonances.

In the present situation, due to the interaction of these two resonances, a chaotic motion appears, but it is temporary since the semi major axis is varying and the double resonance occurs for specific semi major axes given in the above equations. This disturbs the original capture, so that the occurrence of an escape or some change in the libration mode, shown in Fig. 2, seems to be quite natural. However, so far, we are not able to predict those encounters that provoke escapes.

3 Theoretical variation of I and e

$\begin{figure} \par\resizebox{7.6cm}{!}{\includegraphics[clip]{1173fig3.eps}} \end{figure}$

Figure 3: A) circles are solutions of Eq. (2), while the continuous curve is obtained from direct integration of an orbit starting at $a=2.165 R_{\rm M}$ , e=0.01, $I=1^{\circ }$ , $w=0^{\circ }$ , $\Omega =90^{\circ }$ , LL=12, $\epsilon =25.19^{\circ }$ . B) circles are solution of Eq. (3), while the continuous curve results from numerical integration of an orbit starting at: $a=2.558 R_{\rm M}$ , e=0.01, $I=1^{\circ }$ , $\epsilon =25.19^{\circ }$ , $w=\Omega =0^{\circ }$ , LL=10. In both cases, analytical results match very well the numerical integration. ( C)-D)): even considering only one resonant term ( $\cos 2 \varpi-2 \lambda_\odot)$ in $R_\odot$ , both inclination and eccentricity undergo large variations. Initial conditions: a=1.98, e=0.01, $I=33.04^{\circ }$ , $\epsilon =20^{\circ }$ , $w=\Omega =0^{\circ }$ , LL=10. ( E)-F)): similar to ( C)-D)), but for $\dot{w}-n_\odot$ resonance. This time inclination decreases. Initial conditions: $a=1.69R_{\rm M}$ , e=0.01, $I=56.0646^{\circ }$ , $\epsilon =60^{\circ }$ , $w=\Omega =0^{\circ }$ , LL=10.

Open with DEXTER

Assuming only one specific resonant librating angle, Gomes (1997b) was able to find very interesting separate expressions for the evolution of the eccentricity and inclination. They are differential equations that give the variation of these orbital elements in the case of capture. Through these relations, the concepts of universal inclination and eccentricity can be derived.

Here we can also obtain very simple relations giving the evolution of I and e during the capture. On the other hand, the existence of the universal inclination in this problem seems to be not so clear. Jancart, Lemaître & Letocart (2003) working with general dissipation forces could not find a universal value for the inclination. Again, the problem analyzed by these authors involved mean motion resonances. In the present problem, the evolution of the capture may present some particular features, since the resonance relation and the dissipative force are different. In order to see that, we make some simplified assumptions: instead of considering the whole effect of the tidal dissipative force, we assume that the tidal effect is simply given through Eq. (1), that is, only semi major axis variation is taken. Therefore, if Phobos is to be trapped in a specific resonance, the resonant relations like $k_1 \dot{\Omega}+k_2 \dot{w}+k_3 n_\odot \approx 0$ must remain valid all during the capture. For instance, in the $2n_\odot+ \dot{\Omega}$ resonance, the following relation between eccentricity and inclinations should hold:

$\displaystyle a^{7/2} = {3 \sqrt{k^2M}J_2\cos I\over 4n_\odot (1-e^2)^2}\cdot$

(2)

Since this is an inclination type resonance, and provided that the eccentricity variation is small, it is easy to conclude that for a decreasing semi major axis, the inclination can increase with no bounds.

In the present problem, if the initial eccentricity is small, usually, its variation is also very small, so that Eq. (2) is an analytical relation that gives the inclination variation versus semi major axis. Unlike problems involving orbital resonances (Gomes 1997b), in the $2n_\odot+ \dot{\Omega}$ only one resonant angle is in the libration regime. Therefore the theoretical variation of the inclination given by Eq. (2) is always valid until an escape occurs.

In Fig. 3, panel A, the continuous curve comes from direct numerical integration. The circles are solutions of Eq. (2) but we interrupted at $I=46.38^{\circ}$ because the numerical integration shows an escape at this point. The agreement between the continuous and the circled curves is excellent.

The situation is similar for the evection resonance:

$\displaystyle a^{7/2} = {3 \sqrt{k^2M}J_2(5\cos^2 I-2 \cos I-1)\over 4n_\odot (1-e^2)^2}$

(3)

where eccentricity increases when the satellite is receding from the planet (the variation of the inclination is assumed to be negligible). However, in a real and complete tidal model, there will appear a kind of bound on the maximum eccentricity, since as the orbit becomes more eccentric, the angular motion of Phobos, at the pericenter, can be larger than Mars' rotation. This will stop the secular variation of the semi major axis and this is the typical scenario in the case of a more complete tidal model (Touma & Wisdom 1998). In other words, for increasing semi major axis and to study further consequences in a eccentricity resonance capture, the tidal model should be more complete than the simple one given in Eq. (1). In the present problem, Phobos' semi major axis is always decreasing, so that capture in evection resonance or eccentricity type resonances is not expected. Even so, it is instructive to test Eq. (3), since at least during some time, the semi major axis increases until the orbit becomes very eccentric. The theoretical variation of the eccentricity, predicted in Eq. (3), can be compared against the numerical integration. As before, in Fig. 3B the circled curve is obtained from Eq. (3), while the continuous curve corresponds to numerical integration. The inclination remains almost constant ( $\approx$ $1^{\circ}$ ) over all the integration time ( $8 \times 10^8$ years). Again the agreement is quite good.

Therefore, provided that eccentricity (or inclination) does not vary greatly, Eq. (2) (or Eq. (3)) gives the variation of inclination (or eccentricity) with respect to a.

Note, however, that if the initial inclination is not small, we have some cases where even in an eccentricity type resonance, the inclination suffers large variation during the evolution of the capture. Although Eq. (3) remains valid, both inclination and eccentricity can undergo large variation. In Fig. 3C-E we show the significant increase and decrease of the inclination in the $\dot{\varpi}-n_{\odot}$ and $\dot{w}-n_\odot$ resonances respectively. Panel D shows the variation of the eccentricity, and it is the corresponding pair of panel C. Analogously, panel F is the corresponding pair of panel E. In order to decouple and to get separate equations for eccentricity and inclination, we need an extra relation, involving a,e,I. In this sense, we consider, for instance, the action: $\oint p{\rm d}q=\rm const.$ , where p is the canonical momentum associated with the resonant librating angle q. Interesting interpretations and geometric idea of this quantity are given by several authors (Henrard 1982; Henrard & Lemaître 1986; Peale 1989; Gomes 1995a,b, etc.). Let $\nabla$ be the gradient operator while J₁ and J₂ are the two first integrals of the motion. If $\nabla J_1$ and $\nabla J_2$ are not colinear vectors in the (p, q) phase space, then we say that J₁ and J₂ are independent integrals. With this reasoning, although $\oint p{\rm d}q=\rm const.$ and Eq. (3) are not integrals of the motion, separated equations for I and e could be obtained (provided that the previous relations are independent). Usually, an explicit analytical expression for $\oint p{\rm d}q$ is not possible, and in that case we have to resort to approximate numerical methods (Henrard & Lemaître 1986).

4 Mars' eccentricity

In the two previous sections we considered Mars in a circular orbit. The current inclination of Phobos is $\approx$ $1^{\circ}$ . Under normal conditions this value should remain almost unchanged, so that capture of Phobos in the $2n_\odot+ \dot{\Omega}$ commensurability is almost certain (Yokoyama 2002), provided that Mars' obliquity is not less than $\approx$ $19^{\circ}$ . However before reaching this resonance, two previous commensurabilities are encountered: $\dot{\Omega}+ n_\odot$ (N), and $\dot{\varpi}-n_{\odot}$ (E). They appear almost simultaneously at $a \approx 2.619 R_{\rm M}$ so that we called them double resonance. The former, for a decreasing semi major axis, is favorable for capture while the latter (evection), provides capture only if a is increasing. Since they appear simultaneously, capture in (N) resonance, for a decreasing semi major axis, is not certain. Indeed, some numerical experiments in our previous work (PYO), showed that for the circular case, capture in the $\dot{\Omega}+ n_\odot$ resonance usually does not occur. However, a salient fact is noted: some slight jumps in the inclination of Phobos were observed during the passage through these resonances. In the circular case, these jumps are not significant, and therefore, capture in the next resonance, $2n_\odot+ \dot{\Omega}$ , should not be affected. But if the height of these jumps becomes large, the post passage inclination may exceed some threshold value so that the capture in the next $\dot{\Omega}+ 2n_\odot$ resonance may not be certain. So, in this section we show the effect of Mars' eccentricity. It can play an important role mainly in the presence of the double resonance. Besides the resonances given in PYO, the following additional combinations (of order $e _\odot$ ) will appear for $a \leq 2.76 R_{\rm M}$ :

$\displaystyle 2\dot{w}+\dot{\Omega}-3n_\odot$	=	$\displaystyle 0, \hspace {1cm} \mbox { $a \approx 2.619372 R_{\rm M} $ }$
$\displaystyle 2\dot{\Omega}+3n_\odot$	=	$\displaystyle 0, \hspace {1cm} \mbox { $a \approx 2.333045 R_{\rm M} $ }$
$\displaystyle 2\dot{w}+2\dot{\Omega}-3n_\odot$	=	$\displaystyle 0, \hspace {1cm} \mbox { $a \approx 2.332741 R_{\rm M} $ }$
$\displaystyle \dot{\Omega}+3n_\odot$	=	$\displaystyle 0, \hspace {1cm} \mbox {$a \approx 1.913880 R_{\rm M} $ } .$	(4)

The remaining higher order terms are negligible. Notice that $\dot{w}+\dot{\Omega}-3n_\odot$ is just a linear combination of evection (E) and (N) resonances. Therefore, at $a \approx 2.619 R_{\rm M}$ , we have only (N) and (E) resonances. Moreover, both appear not only here, but also in the circular case, so that it is clear that they are much more important than the remaining resonances listed in Eq. (4). Several numerical runs have shown that the consequences and variations caused by these commensurabilities are very weak compared to those caused by (N) and (E) resonances. Therefore, we concentrate our investigation only on this point. Figure 4 is the spectral analysis of the inclination variable considering the interval $(2.6188 R_{\rm M} ,2.6202 R_{\rm M} )$ of Phobos' semi major axis. A total of 120 initial conditions were taken in this interval and the solutions were analyzed via FFT. Due to the amount of computation time in obtaining this figure, we considered only two cosines, corresponding to evection (E) and (N) resonances, in the disturbing function. Moreover, Mars was considered in a circular orbit. For each semi major axis, the problem is integrated and the log of the fundamental frequencies and their harmonics are shown in Fig. 4. Those frequencies whose amplitudes are $\leq$ 0.5% of $\mathcal {M}$ are discarded, where $\mathcal {M}$ is the maximum value of the amplitude for each integration (see Callegari et al. 2004). The curves show the behaviour of the frequencies for different values of the semi major axis. In the case of regular motion, their paths are well defined and vary continuously with respect to a. In Fig. 4, we use E (or N) to identify the curve whose frequency is due to E (or N) resonance. At about $a \approx 2.6197 R_{\rm M}$ we see two well defined frequencies, almost symmetric with respect to this value of a. Both correspond to the N frequency although only the curve on the left is marked. The remaining curves not labelled with these two letters correspond to higher harmonics of the E frequency or linear combinations between E and N frequencies. In particular, higher harmonics of the N frequency have very small amplitudes, so that they do not appear in the figure. At $a \approx 2.6192 R_{\rm M}$ and $a \approx 2.6197 R_{\rm M}$ , it is clear that the two fundamental frequencies fall to zero. At these points, the corresponding periods of (E) resonance and (N) resonance go to infinity, which characterizes the separatrices of both resonances. In the neighbourhood and in the region where these separatrices overlap, we have some chaotic motion. Since these two separatrices are very close to each other, we expect the appearance of chaotic motion in more the general case.

$\begin{figure} \par\resizebox{8cm}{!}{\includegraphics[clip]{1173fig4.eps}} \end{figure}$

Figure 4: Evolution of the log of the frequencies versus semi major axis. At the point where the fundamental frequency goes to zero we have the corresponding semi major of the separatrix. Since these separatrices are close, they can overlap and we expect the onset of chaos. Each integration was extended to 1 996 050 years, generating 2¹⁶ outputs, sampled at every 30 years. Phobos initial conditions: e=0.051, $I=1.075^{\circ }$ , $\varpi =150.247^{\circ }$ , $\Omega =164.931^{\circ }$ . Mars initial conditions: e=0.0, $I=1.5236^{\circ }$ , $\varpi =0^{\circ }$ , $\Omega =0^{\circ }$ , $\lambda =0^{\circ }$ , $\epsilon =25.19^{\circ }$ .

Open with DEXTER

We examine the dynamics when the semi major of Phobos varies, passing through these two resonant values. This time, Mars' eccentricity is considered (up to fourth order terms) and also all cosines of the disturbing function given in Eq. (5) in PYO.

$\begin{figure} \par\resizebox{8.2cm}{!}{\includegraphics[clip]{1173fig5.eps}} \end{figure}$

Figure 5: Due to the action of the two simultaneous resonances in the neighbourhood of $a \approx 2.619 R_{\rm M}$ , the motion is very irregular (chaotic). Depending on the input values, the inclination and eccentricity are strongly excited during the passage through the double resonance. The final stabilized values are not predictable. On the left, the system leaves the double resonance with inclination $\approx$ $4^{\circ }$ and eccentricity $\approx$ 0.08. On the right, these values are smaller than those used in the input data. Initial conditions ( left panels): a=2.625, e=0.015, $I=1^{\circ }$ , $\epsilon =25.19^{\circ }$ , $w=\Omega =0^{\circ }$ , ( right panels): a=2.630, e=0.015, $I=1^{\circ }$ , $\epsilon =25.19^{\circ }$ , $w=\Omega =0^{\circ }$ .

Open with DEXTER

Figure 5 shows the situation when Phobos encounters (E) and (N) resonances. As expected, it is clear that the motion is quite irregular in the vicinity of $a \approx 2.619 R_{\rm M}$ . Because of this, in all integrations of this figure, we again fixed high precision in the Radau integrator (LL=12). In Fig. 5A the inclination goes through jumps during the passage through the critical value of the semi major axis, but after that it stabilizes at $I\approx 4^{\circ}$ . Also the eccentricity is excited to $e\approx 0.08$ (Fig. 5C). This significant increase in the inclination is important since unless some extra damping occurs, the next resonance ( $2n_\odot+ \dot{\Omega}$ ) will be faced with an inclination so large that capture in this resonance may not be possible. Indeed, certain capture is possible only when the initial action $\oint p{\rm d}q$ (see PYO, Eq. (9)), far from the resonance, is less than the area inside of the separatrix first formed, during the evolution of the semi major axis (Henrard 1982; Peale 1989). This means that we cannot guarantee capture for an arbitrary initial inclination. For instance, several numerical experiments show that for $I > 1.5^{\circ}$ capture depends on the phase, while for $I \leq 1^{\circ}$ we always can get capture no matter the initial values for w or $\Omega$ .

Figure 5B shows another interesting situation: as before, the motion is also very complicated due to the double resonance and the inclination is excited to high values. However, at the end, it decreases back to about $I\approx 0.08^{\circ}$ . The final eccentricity is also smaller than the initial value before entering the critical region (Fig. 5D). Since the motion is very irregular it seems that the precision parameter LL=12 in the Radau integrator is important. In this sense, a curious case appears if LL is relaxed to LL=9: we have a capture in the $\dot{\Omega}+ n_\odot$ and inclination increases to high values similar to those shown in Fig. 2. After a lot of numerical experiments, we prefer to say that this capture is improbable and that the percentage of initial conditions leading to capture is very small. Moreover, tests showed that the general dynamics is very sensitive to the starting point and to the precision parameter LL used in the integrator. In particular, the initial conditions we used to get the mentioned capture were: $a=2.622 R_{\rm M}$ , e=0.015, $I=1^{\circ }$ , $w=0^{\circ }$ , $\Omega =0^{\circ }$ , $\epsilon =25.19^{\circ }$ , with LL=9.

We believe that the inclusion of Mars' eccentricity enhances the effect of the double resonance. The final value of the inclination when the semi major axis crosses the region $a \approx 2.619 R_{\rm M}$ seems to be completely uncertain.

5 Planetary perturbations

Up to now, the orbit of Mars was considered Keplerian with a fixed equator. In a more realistic model, we should take into account the planetary perturbations which will disturb not only Mars' orbit but also its equator. However, an integration starting from present data until Phobos semi major axis decreases to $a \approx 2.619 R_{\rm M}$ or $a\approx 2.149$ , needs an enormous and probably a prohibitive computational effort (Hamilton 2002). Without a reliable tidal model, and corresponding reliable constants, we cannot be confident of the final result of this long time integration. Usually tidal models demand knowledge of important physical constants, but in general, the accuracy of these values is always questionable: the dissipation function Q, Love's number $k_{\rm s}$ , mass satellite (Phobos) and also rigidity and dissipation in Phobos. For the precession of Mars' equator, the problem is the same and the important constants involved are: the moment of inertia and spin rate. In spite of all these difficulties we think it is interesting to have a rough qualitative idea of the main effects when planetary perturbations are included. To avoid the problem of computational time, we artificially change the initial semi major of Phobos, and start our integration in the neighbourhood of $a \approx 2.149 R_{\rm M}$ and not $2.76 R_{\rm M}$ , which is the current value of Phobos' semi major axis. In order to consider the perturbation of the planets on Mars, we use Laskar's secular theory (Laskar 1988). We also take the variation of Mars equator through the equations (Woolard 1953):

$\displaystyle {\rm d}\tilde{\Omega} \over {\rm d}t$	=	$\displaystyle -K_0 \bigg \{ (1-3\cos^2 I_\odot)\cos \tilde{I}+ {\sin 2I_\odot \cos 2\tilde{I}\over {\sin \tilde{I}}}$
		$\displaystyle \cos(\tilde{\Omega}-\Omega_\odot)+ \sin^2 I_\odot \cos\tilde{I} \cos(2\tilde{\Omega}-2\Omega_\odot) \bigg \} K_1$
$\displaystyle {\rm d}\tilde{I} \over {\rm d}t$	=	$\displaystyle K_0 \left \{ -\sin \tilde{I} \sin^2 I_\odot\sin(2\tilde{\Omega}-2\Omega_\odot) \right.$
		$\displaystyle \left. -\cos\tilde{I}\sin 2I_\odot\sin(\tilde{\Omega}-\Omega_\odot) \right \} K_1$	(5)

where

K₀	=	$\displaystyle -3k^2M_\odot M_M R^2_{\rm M}J_2\over {Cw_z}$
K₁	=	$\displaystyle {1\over 4 a_\odot^3(1-e^2_{\odot})^{3/2}} \cdot$

In the above equations, we take the Sun as the main disturber of Mars equator (Woolard 1953). $\tilde {I}$ and $\tilde {\Omega}$ are the inclination and longitude of the node of Mars' equator referred to an inertial reference plane (ecliptic of 2000). These equations are averaged with respect to the mean anomaly of the Sun. The notation of the variables and constants are the same as used in PYO, with two new constants: C=polar moment of inertia, w_z = spin rate. In this work we adopted $C/M_{\rm M} R_{\rm M}^2=0.366$ and the Mars rotation period was fixed to be 1.026 days.

The Laskar secular theory (Laskar 1988) is given through the series:

$\displaystyle \alpha = \sum_{i=1} ^N A_i \exp \sqrt(-1) (\nu_i t+\phi_i)$

(6)

where in general only 25 terms in the above series were considered (we are interested only in the qualitative behaviour of the solution). The frequencies ( $\nu_i$ ), phases ( $\phi_i$ ) and amplitudes (A_i) of each orbital element ( $\alpha$ ) are given in Laskar (1988), Table 9.

$\begin{figure} \par\resizebox{8.5cm}{!}{\includegraphics[clip]{1173fig6.eps}} \end{figure}$

Figure 6: A typical situation of capture in the $2n_\odot+ \dot{\Omega}$ resonance including planetary perturbation on Mars and variation of its equator. In panel A, Phobos' inclination with respect to the ecliptic of 2000 is shown, while in panel B, the same inclination is shown with respect to the moving Mars equator. Panel C, shows the behaviour of $\Omega-\tilde{\Omega}$ . Initially these two nodes are almost synchronous and after the capture, they start to diverge ( more details are given in the text). Panel D, exhibits the variation of the obliquity. Initial conditions a=2.152, e=0.015, $I=27^{\circ }$ , $\Omega =83.822^{\circ }$ , $\tilde{\Omega}=82.908^{\circ}$ , $\tilde{I}=26.718^{\circ}$

Open with DEXTER

The equations to be integrated are the same as before, but this time the orbit of Mars is given through series 6 and we also add Eq. (5). As usual, the initial conditions are chosen such that Phobos' inclination with respect to Mars' equator is $\leq$ $1^{\circ}$ and $a\approx 2.15$ . Figure 6 shows a case of a capture in the $2n_\odot+ \dot{\Omega}$ resonance, when planetary perturbations on Mars and motion of its equator are considered.

Panel A shows the variation of the inclination of Phobos with respect to the ecliptic plane. A fast oscillation around the initial inclination ( $27^{\circ}$ ) is noticed. The amplitude of this oscillation starts to increase as the capture occurs. In fact, this inclination may reach zero or values greater than $90^{\circ}$ (Gomes 1997b), so that regularized variables are needed in these cases. Panel B shows the corresponding inclination referred to the equatorial plane. This time, once capture occurs, instead of the fast oscillations of the previous case, we have the well-known slow and secular increasing behaviour of the inclination. Panel C shows the difference $\Omega-\tilde{\Omega}$ . Its behaviour is very similar to the inclination in panel A. Finally, the variation of the obliquity ( $\epsilon$ ) of Mars is shown in panel D. A few tests have shown that sometimes capture fails if the initial inclination with respect to the equator is $1^{\circ}$ , while we always found successful capture when we took it as $\approx$ $0.5^{\circ}$ . Therefore we see that even considering planetary perturbations and the variation of the equator, capture still occurs and the general path of the variation of the inclination with respect to the equator is quite similar to the cases without these two effects.

6 Conclusions

Considering the simplest tidal model for Phobos, we analyzed some situations of escape from the $2n_\odot+ \dot{\Omega}$ resonance. The escapes are due to new resonances that appear as the inclination reaches some specific values: $I= 33^{\circ}$ , $46.38^{\circ}$ , $56.06^{\circ}$ , $63.43^{\circ}$ , $69.02^{\circ}$ , $73.15^{\circ}$ , $76.23^{\circ}$ . However, we cannot predict the inclinations for which escapes will occur. After several numerical simulations, we conclude that the existence of a universal inclination is not clear in this problem, since no stabilization of the inclination was found. Some analytic equations that give the path followed by the inclination during the capture are derived for each resonance. We also show that for $a \approx 2.619 R_{\rm M}$ , the role played by the eccentricity of Mars is decisive, since Phobos' inclination can be excited to about $4^{\circ }$ , so that capture in $2n_\odot+\Omega$ is not possible any longer. Finally, including planetary perturbations on Mars and precession equations for the motion of its equator, we showed that capture of Phobos is still possible and probably, the general features of the dynamics should be similar to the simplified model.