© ESO, 2010

1. Introduction

Understanding the rotation of gravitationally interacting bodies with an ever increasing accuracy is an important challenge of celestial mechanics. Many authors have modeled with very high accuracy the rotation of the Earth for instance Woolard (1953), Kinoshita (1977). The other planets of the Solar system followed logically (see for instance Yoder 1997; Rambaux 2007). Today, thanks to new and very precise observational astrometric techniques, the dynamical study of celestial bodies extends to comets, asteroids and satellites (Meyer 2008; Sinclair 1977; Wisdom 1984, 1987).

Yet only a few studies have been made of the rotational motion of Phoebe, the ninth satellite of Saturn, since its discovery in 1899 by Pickering. In 1905, Ross etablished that its motion is retrograde and gave for the first time the orbital elements for the mean equinox and ecliptic of date. Jacobson (1998) determined new orbital elements thanks to Earth-based astrometric observations from 1904 to 1996. More recently Emelyanov(2007) made ephemerides of Phoebe from 1904 to 2027. Others studies have been made of the composition and the inertial parameters of Phoebe (Aleshkina et al. 2010). Phoebe is the only non-synchronous satellite of Saturn with rather well known physical parameters, such as the moments of inertia. Because Phoebe is a dissymetric body (large dynamical flattening and large triaxiality), it is interesting to study its rotation.

We propose here for the first time, to determine the combined motion of precession and nutation of Phoebe considered as a rigid body. After presenting the model (Sect. 2), we describe the orbital motion of Phoebe by fitting the curves of the temporal variations of the orbital elements a,e,M and L_s with polynomial functions. Thanks to a fast fourier transform (FFT), the large periodic variations around the mean elements are analyzed and show that the orbital motion of Phoebe is far from being Keplerian (Sect. 3). Then we determine the motion of precession and nutation of Phoebe by numerical integration using the value of the obliquity obtained in Sect. 4. Finally, after quantifying the validity of the developments done by Kinoshita (1977) for the Earth when applied to Phoebe, an analytical model is constructed.

All our results are compared to the results obtained by Kinoshita for the Earth. We show in Sect. 6 that the analytical model which is very accurate for the Earth is less acurate for Phoebe than the numerical approach to determine the precession and the nutation of Phoebe. This is due to the large periodic departure from the Keplerian motion of the Phoebe’s orbit. We note however that it gives a fairly good approximation for the precession, which is sensitive to the average value.

2. Theoretical equation of the precession-nutation motion

The present study of the rotation of Phoebe considered as a rigid body is based on the theoritical framework by Kinoshita (1977), who adopted an Hamiltonian formulation using Andoyer variables (1923) to evaluate the motion of the precession-nutation of the rigid Earth. Here a simplified model is described and only the variables and the equations used in this paper are detailed. The rotation of a rigid Phoebe involves three axes: the figure axis, which coincides with the axis of the largest moment of inertia, the angular momentum axis directed along , with G the amplitude of the angular momentum and a reference axis arbitrarly chosen as the axis perpendicular to the orbit of Phoebe at a given epoch t (see Fig. 1). The precession and nutation motions are the linear and the quasi-periodic parts of the motion of the figure axis or of the angular momentum axis with respect to the reference axis. Assuming that the angle J between the angular momentum axis and the figure axis is small, as it is for the Earth (for which its value is less than 1′′), only the motion of the angular momentum axis is considered here. The motion of this axis with respect to the reference axis is described by the angles h,I (Andoyer 1923), where I is the obliquity angle and h is characterizing the precession-nutation in longitude, which corresponds to the angle between the reference point γ_t and the node Q. The reference point γ_t is the intersection between the orbit and the equator of Phoebe at the epoch t , the so-called “departure point” (Capitaine 1986), and Q is the ascending node between the plane normal to the angular momentum and the orbital plan. The Andoyer variables g,l in Fig. 1 are the angle between the node Q and the node P and the angle between a meridian origin and the node P, where P is itself the ascending node between the plane normal to the angular momentum axis and the equatorial plane. The proper rotation of Phoebe is described by the angle l + g = Φ. The Hamiltonian related to the rotational motion of Phoebe is $\mathit{K}\mathrm{=}{\mathit{F}}_{\mathrm{o}}\mathrm{+}\mathit{E}\mathrm{+}{\mathit{E}}^{\mathrm{\prime }}\mathrm{+}\mathit{U}$(1)F_o is the Hamiltonian for the free rotational motion defined by ${\mathit{F}}_{\mathrm{0}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}\left(\frac{{\sin }^{\mathrm{2}}\mathit{l}}{\mathit{A}}\mathrm{+}\frac{{\cos }^{\mathrm{2}}\mathit{l}}{\mathit{B}}\right)\mathrm{(}{\mathit{G}}^{\mathrm{2}}\mathrm{-}{\mathit{L}}^{\mathrm{2}}\mathrm{)}\mathrm{+}\frac{\mathrm{1}}{\mathrm{2}}\frac{{\mathit{L}}^{\mathrm{2}}}{\mathit{C}}\mathit{,}$(2)where L is the component of the angular momentum axis along the figure axis and A,B,C are the principal moments of inertia of Phoebe. E + E′ is a component related to the motion of the orbit of Phoebe, which is caused by planetary perturbations and has been discussed in detail in Cottereau & Souchay (2009) when studying the rotation of Venus. U is the disturbing potential of Saturn, considered as a point mass and its disturbing potential is given by

$\mathit{U}\mathrm{=}{\mathit{U}}_{\mathrm{1}}\mathrm{+}{\mathit{U}}_{\mathrm{2}}$$\begin{array}{ccc}{\mathit{U}}_{\mathrm{1}}& \mathrm{=}& \frac{{G}{\mathit{M}}^{\mathrm{\prime }}}{{\mathit{r}}^{\mathrm{3}}}[\left[\frac{\mathrm{2}\mathit{C}\mathrm{-}\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{2}}\right]{\mathit{P}}_{\mathrm{2}}\mathrm{\left(}\sin \mathit{\delta }\mathrm{\right)}\\ & \mathrm{+}& \left[\frac{\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{4}}\right]{\mathit{P}}_{\mathrm{2}}^{\mathrm{2}}\mathrm{\left(}\sin \mathit{\delta }\mathrm{\right)}\cos \mathrm{2}\mathit{\alpha }]\end{array}$(3)$\begin{array}{ccc}{\mathit{U}}_{\mathrm{2}}& \mathrm{=}& \sum _{\mathit{n}\mathrm{=}\mathrm{3}}^{\mathrm{\infty }}\frac{{G}{\mathit{M}}^{\mathrm{\prime }}{\mathit{M}}_{\mathit{P}}{\mathit{a}}^{\mathit{n}}}{{\mathit{r}}^{\mathit{n}\mathrm{+}\mathrm{1}}}\mathrm{\left[}{\mathit{J}}_{\mathit{n}}{\mathit{P}}_{\mathit{n}}\mathrm{\right(}\sin \mathit{\delta }\mathrm{)}\\ & \mathrm{-}& \sum _{\mathit{m}\mathrm{=}\mathrm{1}}^{\mathit{n}}{\mathit{P}}_{\mathit{n}}^{\mathit{m}}\mathrm{\left(}\sin \mathit{\delta }\mathrm{\right)}\mathrm{·}\mathrm{(}{\mathit{C}}_{\mathit{nm}}\cos \mathit{m\alpha }\mathrm{+}{\mathit{S}}_{\mathit{nm}}\sin \mathit{m\alpha }\mathrm{)}\mathrm{]}\mathit{,}\end{array}$(4)

where G is the gravitational constant, M′ is the mass of Saturn, r is the distance between its barycenter and the barycenter of Phoebe. α and δ are the longitude and latitude of Saturn (not to be confused with the usual equatorial coordinates), with respect to the mean equator of Phoebe and a meridian origin. The are the classical Legendre functions given by: ${\mathit{P}}_{\mathit{n}}^{\mathit{m}}\mathrm{\left(}\mathit{x}\mathrm{\right)}\mathrm{=}\frac{\mathrm{(}\mathrm{-}\mathrm{1}{\mathrm{)}}^{\mathit{m}}\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{x}}^{\mathrm{2}}{\mathrm{)}}^{\frac{\mathit{m}}{\mathrm{2}}}}{{\mathrm{2}}^{\mathit{n}}\mathit{n}\mathrm{!}}\frac{{\mathrm{d}}^{\mathit{n}\mathrm{+}\mathit{m}}\mathrm{(}{\mathit{x}}^{\mathrm{2}}\mathrm{-}\mathrm{1}{\mathrm{)}}^{\mathit{n}}}{{\mathrm{d}}^{\mathit{n}\mathrm{+}\mathit{m}}\mathit{x}}\mathrm{·}$(5)Note that the perturbations caused by, the other planets, satellites and the Sun, which are a priori of second order, will not be studied in this paper. The variation of the obliquity and the precession angle is given by $\begin{array}{ccc}\frac{\mathrm{d}\mathit{I}}{\mathrm{d}\mathit{t}}& \mathrm{=}& \frac{\mathrm{1}}{\mathit{G}}\left[\frac{\mathrm{1}}{\sin \mathit{I}}\frac{\mathit{\partial K}}{\mathit{\partial h}}\mathrm{-}\cot \mathit{I}\frac{\mathit{\partial K}}{\mathit{\partial g}}\right]\\ \frac{\mathrm{d}\mathit{h}}{\mathrm{d}\mathit{t}}& \mathrm{=}& \mathrm{-}\frac{\mathrm{1}}{\mathit{G}\sin \mathit{I}}\frac{\mathit{\partial K}}{\mathit{\partial I}}\mathrm{·}\end{array}$(6)To solve these equations, Kinoshita (1977) used Hori’s method. Because the order of the disturbing function is given by $\frac{{\mathit{U}}_{\mathrm{1}}}{{\mathit{F}}_{\mathrm{o}}}\mathrm{\approx }\frac{{\mathit{M}}^{\mathrm{\prime }}}{{\mathit{M}}^{\mathrm{\prime }}\mathrm{+}{\mathit{M}}_{\mathrm{p}}}\mathrm{\left(}\frac{\mathit{n}}{\mathit{\omega }}{\mathrm{)}}^{\mathrm{2}}\mathrm{\right(}\frac{\mathrm{2}\mathit{C}\mathrm{-}\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{2}\mathit{C}}\mathrm{)}\mathrm{\approx }\mathrm{3.46}\mathrm{\times }{\mathrm{10}}^{-8}$(7)the same method can be applied to Phoebe. Neglecting the very small contribution as the component E + E′ and applying the first order of Hori’s method, this yields $\begin{array}{ccc}\mathrm{\Delta }\mathit{I}\mathrm{=}\frac{\mathrm{1}}{\mathit{G}}\left[\frac{\mathrm{1}}{\sin \mathit{I}}\frac{\mathit{\partial }}{\mathit{\partial h}}\mathrm{\int }{\mathit{U}}_{\mathrm{1}}\mathrm{d}\mathit{t}\mathrm{-}\cot \mathit{I}\frac{\mathit{\partial }}{\mathit{\partial g}}\mathrm{\int }{\mathit{U}}_{\mathrm{1}}\mathrm{d}\mathit{t}\right]& & \end{array}$(8)$\begin{array}{ccc}\mathrm{\Delta }\mathit{h}\mathrm{=}\mathrm{-}\frac{\mathrm{1}}{\mathit{G}\sin \mathit{I}}\frac{\mathit{\partial }}{\mathit{\partial I}}\mathrm{\int }{\mathit{U}}_{\mathrm{1}}\mathrm{d}\mathit{t}\mathit{.}& & \end{array}$(9)U₁ given by (3) can be expressed as a function of the longitude λ and the latitude β of Saturn with respect to the orbit of Phoebe at the epoch t with the transformations described by Kinoshita (1977) and based on the Jacobi polynomials: $\begin{array}{ccc}{\mathit{U}}_{\mathrm{1}}& \mathrm{=}& \frac{{G}{\mathit{M}}^{\mathrm{\prime }}}{{\mathit{r}}^{\mathrm{3}}}[\frac{\mathrm{2}\mathit{C}\mathrm{-}\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{2}}(\mathrm{-}\frac{\mathrm{1}}{\mathrm{4}}\mathrm{(}\mathrm{3}{\cos }^{\mathrm{2}}\mathit{I}\mathrm{-}\mathrm{1}\mathrm{)}\\ & \mathrm{-}& \frac{\mathrm{3}}{\mathrm{4}}{\sin }^{\mathrm{2}}\mathit{I}\cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{)})\\ & \mathrm{+}& \frac{\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{4}}[\frac{\mathrm{3}}{\mathrm{2}}{\sin }^{\mathrm{2}}\mathit{I}\cos \mathrm{(}\mathrm{2}\mathit{l}\mathrm{+}\mathrm{2}\mathit{g}\mathrm{)}\\ & \mathrm{+}& \sum _{\mathit{ϵ}\mathrm{=}\mathrm{\pm }\mathrm{1}}\frac{\mathrm{3}}{\mathrm{4}}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{ϵ}\cos \mathit{I}{\mathrm{)}}^{\mathrm{2}}\mathrm{·}\cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{-}\mathit{ϵl}\mathrm{-}\mathit{ϵg}\mathrm{)}]]\mathit{.}\end{array}$(10)Finally the motions of precession-nutation in longitude and in obliquity are given starting from Eqs. (8)–(10) by $\begin{array}{ccc}\mathrm{\Delta }\mathit{h}& \mathrm{=}& \mathrm{-}\frac{\mathrm{1}}{\sin \mathit{I}}\frac{\mathit{\partial }}{\mathit{\partial I}}\mathrm{\int }\mathrm{(}{\mathit{W}}_{\mathrm{1}}\mathrm{+}{\mathit{W}}_{\mathrm{2}}\mathrm{)}\mathrm{d}\mathit{t}\end{array}$(11)$\begin{array}{ccc}\mathrm{\Delta }\mathit{I}& \mathrm{=}& [\frac{\mathrm{1}}{\sin \mathit{I}}\frac{\mathit{\partial }}{\mathit{\partial h}}\mathrm{\int }\mathrm{(}{\mathit{W}}_{\mathrm{1}}\mathrm{+}{\mathit{W}}_{\mathrm{2}}\mathrm{)}\mathrm{d}\mathit{t}\\ & \mathrm{-}& \cot \mathit{I}\frac{\mathit{\partial }}{\mathit{\partial g}}\mathrm{\int }{\mathit{W}}_{\mathrm{2}}\mathit{dt}]\mathit{,}\end{array}$(12)where $\begin{array}{ccc}{\mathit{W}}_{\mathrm{1}}& \mathrm{=}& \frac{{G}{\mathit{M}}^{\mathrm{\prime }}}{\mathit{G}{\mathit{r}}^{\mathrm{3}}}(\frac{\mathrm{2}\mathit{C}\mathrm{-}\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{2}}[\mathrm{-}\frac{\mathrm{1}}{\mathrm{4}}\mathrm{(}\mathrm{3}{\cos }^{\mathrm{2}}\mathit{I}\mathrm{-}\mathrm{1}\mathrm{)}\\ & \mathrm{-}& \frac{\mathrm{3}}{\mathrm{4}}{\sin }^{\mathrm{2}}\mathit{I}\cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{)}])\\ {\mathit{W}}_{\mathrm{2}}& \mathrm{=}& \frac{{G}{\mathit{M}}^{\mathrm{\prime }}}{\mathit{G}{\mathit{r}}^{\mathrm{3}}}(\frac{\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{4}}[\frac{\mathrm{3}}{\mathrm{2}}{\sin }^{\mathrm{2}}\mathit{I}\cos \mathrm{(}\mathrm{2}\mathit{l}\mathrm{+}\mathrm{2}\mathit{g}\mathrm{)}\\ & \mathrm{+}& \sum _{\mathit{ϵ}\mathrm{=}\mathrm{\pm }\mathrm{1}}\frac{\mathrm{3}}{\mathrm{4}}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{ϵ}\cos \mathit{I}{\mathrm{)}}^{\mathrm{2}}\\ & & \cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{-}\mathit{ϵl}\mathrm{-}\mathit{ϵg}\mathrm{)}])\mathit{.}\end{array}$Because W₁ does not depend on g, this variable does not appear in the second part of the right hand side of Eq. (12). Supposing that the components ω₁ and ω₂ of Phoebe’s rotation are negligible with respect to the component ω₃ along the figure axis, as it is the case for the Earth, this yields $\begin{array}{ccc}{\mathit{W}}_{\mathrm{1}}& \mathrm{=}& \frac{\mathrm{3}{G}{\mathit{M}}^{\mathrm{\prime }}}{\mathit{\omega }{\mathit{r}}^{\mathrm{3}}}(\frac{\mathrm{2}\mathit{C}\mathrm{-}\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{2}\mathit{C}}[\mathrm{-}\frac{\mathrm{1}}{\mathrm{12}}\mathrm{(}\mathrm{3}{\cos }^{\mathrm{2}}\mathit{I}\mathrm{-}\mathrm{1}\mathrm{)}\\ & \mathrm{-}& \frac{\mathrm{1}}{\mathrm{4}}{\sin }^{\mathrm{2}}\mathit{I}\cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{)}])\\ {\mathit{W}}_{\mathrm{2}}& \mathrm{=}& \frac{\mathrm{3}{G}{\mathit{M}}^{\mathrm{\prime }}}{\mathit{\omega }{\mathit{r}}^{\mathrm{3}}}(\frac{\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{4}\mathit{C}}[\frac{\mathrm{1}}{\mathrm{2}}{\sin }^{\mathrm{2}}\mathit{I}\cos \mathrm{(}\mathrm{2}\mathit{l}\mathrm{+}\mathrm{2}\mathit{g}\mathrm{)}\\ & \mathrm{+}& \sum _{\mathit{ϵ}\mathrm{=}\mathrm{\pm }\mathrm{1}}\frac{\mathrm{1}}{\mathrm{4}}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{ϵ}\cos \mathit{I}{\mathrm{)}}^{\mathrm{2}}\\ & & \cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{-}\mathit{ϵl}\mathrm{-}\mathit{ϵg}\mathrm{)}])\mathit{,}\end{array}$where and are the dynamical flattening and the triaxiality of Phoebe. Because the rotation of Phoebe is retrograde, the conventional notations of the precession and obliquity will be respectively $\mathrm{\Delta }\mathit{\psi }\mathrm{=}\mathrm{\Delta }\mathit{h}$(17)$\mathrm{\Delta }\mathit{ϵ}\mathrm{=}\mathrm{\Delta }\mathit{I,}$(18)which are the opposite of the conventional ones for the Earth (Kinoshita 1977; Souchay et al. 1999).

3. Phoebe’s ephemerides

The ephemerides used in this paper for Phoebe’s orbital motion around Saturn have been constructed by Emelyanov (2007), where the orbital parameters and the coordinates (X,Y,Z) of Phoebe are given in the planeto-equatorial reference frame considered as inertial from the first observation date in 1904 to 2027. The reference point is one of the intersections between Saturn and Earth equatorial planes at J2000.0. To apply the theoretical framework of Kinoshita (1977), the mean elements for the variables a,e,M and L_s are needed. In this section, these elements are calculated in the same way as in Simon et al. (1994) for the planets. The differences betwen the true motion of Phoebe and a Keplerian motion are also studied. Throughout the paper our study is restricted to the time span ranging from 2000.0 to the limit of the ephemeris 2027.0.

3.1. The semi-major axis

The mean element for the semi-major axis for the Earth is a constant (Simon et al. 1994) and the relative variations around this mean value are on the order 10^-5. In comparaison, for Phoebe, approximating a by a constant does not yield a good accuracy. Figure 2 shows the relative variation of a around its mean value estimated at in a 9000 days’ time span. The variations clearly consist of periodic components with an amplitude on the order of 10^-3. The leading oscillations of the signal are determined with a fast fourier transform (FFT). Table 1 gives the leading amplitudes and periods of the sinusoids characterizing the signal in Fig. 2, where the central curve represents the residuals after substraction of these sinusoids.

The relative large periodic components caused by the additional perturbations on the orbit of Phoebe show that the motion of this satellite is not close to a Keplerian one, as it is the case for the planets of the Solar system. The effects of the attraction of the Sun, Jupiter, and the other satellites of Saturn on Phoebe’s orbit entail an important departure from the Keplerian motion. The individual study of the other orbital elements confirms this result.

3.2. The mean elements of the orbital parameters

Below the formula for the mean elements of Phoebes′ Keplerian motion e,M and L_S are given, where M,L_s are the mean anomaly and mean longitude with respect to the “departure point” γ_t. They were obtained by fitting the curves of the temporal variations of these elements given by Emelyanov (2007) with a polynomial expression at 6th order as was done for the planets of the solar system (Simon et al. 1994). Considering the limits of the ephemerides, truncating the polynomial functions at the 6th order looks a sufficient approximation. Thus $\begin{array}{ccc}\mathit{e}& \mathrm{\approx }& \mathrm{0.168480}\mathrm{-}\mathrm{0.0000594}\mathrm{·}\mathit{t}\mathrm{+}\mathrm{4.931508}\mathrm{\times }{\mathrm{10}}^{-8}\mathrm{·}{\mathit{t}}^{\mathrm{2}}\\ & \mathrm{-}& \mathrm{1.30162}\mathrm{\times }{\mathrm{10}}^{-11}\mathrm{·}{\mathit{t}}^{\mathrm{3}}\mathrm{+}\mathrm{1.1744}\mathrm{\times }{\mathrm{10}}^{-15}\mathrm{·}{\mathit{t}}^{\mathrm{4}}\\ & \mathrm{-}& \mathrm{5.2502}\mathrm{\times }{\mathrm{10}}^{-22}\mathrm{·}{\mathit{t}}^{\mathrm{5}}\mathrm{-}\mathrm{3.0861}\mathrm{\times }{\mathrm{10}}^{-24}\mathrm{·}{\mathit{t}}^{\mathrm{6}}\mathit{,}\end{array}$(19)$\begin{array}{ccc}\mathit{M}& \mathrm{\approx }& \mathrm{58.9233}\mathrm{+}\mathrm{0.656386}\mathrm{·}\mathit{t}\mathrm{-}\mathrm{0.000017}\mathrm{·}{\mathit{t}}^{\mathrm{2}}\\ & \mathrm{+}& \mathrm{9.78065}\mathrm{\times }{\mathrm{10}}^{-9}\mathrm{·}{\mathit{t}}^{\mathrm{3}}\mathrm{-}\mathrm{2.1541}\mathrm{\times }{\mathrm{10}}^{-12}\mathrm{·}{\mathit{t}}^{\mathrm{4}}\\ & \mathrm{+}& \mathrm{2.06592}\mathrm{\times }{\mathrm{10}}^{-16}\mathrm{·}{\mathit{t}}^{\mathrm{5}}\mathrm{-}\mathrm{7.2487}\mathrm{\times }{\mathrm{10}}^{-21}\mathrm{·}{\mathit{t}}^{\mathrm{6}}\end{array}$(20)$\begin{array}{ccc}{\mathit{L}}_{\mathit{s}}& \mathrm{\approx }& \mathrm{331.109}\mathrm{+}\mathrm{0.653712}\mathrm{·}\mathit{t}\mathrm{-}\mathrm{3.64487}\mathrm{\times }{\mathrm{10}}^{-7}\mathrm{·}{\mathit{t}}^{\mathrm{2}}\\ & \mathrm{+}& \mathrm{3.54186}\mathrm{\times }{\mathrm{10}}^{-10}\mathrm{·}{\mathit{t}}^{\mathrm{3}}\mathrm{-}\mathrm{9.42668}\mathrm{\times }{\mathrm{10}}^{-14}\mathrm{·}{\mathit{t}}^{\mathrm{4}}\\ & \mathrm{+}& \mathrm{1.03032}\mathrm{\times }{\mathrm{10}}^{-17}\mathrm{·}{\mathit{t}}^{\mathrm{5}}\mathrm{-}\mathrm{4.02263}\mathrm{\times }{\mathrm{10}}^{-22}\mathrm{·}{\mathit{t}}^{\mathrm{6}}\mathit{,}\end{array}$(21)where t is counted in Julian days, M and L_s are expressed in degrees. Figure 3 shows the variations of the eccentricity of Phoebe for a 9000 days’ time span. An important feature of Phoebes’ orbit is the high eccentricity as well as its very large relative variations (of the order of 10%). As we will see in Sect. 6, this will be very important in our study because of the large corresponding variations of that it causes for the orbit.

Figures 4–6 represent the residuals after substraction of the polynomial functions given by Eqs. (19)–(21)from the eccentricity, the mean anomaly and the mean longitude of Phoebe. As already observed for the semi-major axis, the residuals have periodic components, which reach an amplitude of 4° for the mean anomaly.

3.3. Frequencies analysis

Thanks to a fast fourier transform (FFT) analysis the leading ampitudes and the periods of the sinusoids characterizing the signals in Figs. 4–6 are determined and given in Table 2.

The central curves in Figs. 4–6 represent the residuals after substraction of the leading sinusoidal components of the mean elements. We highlight here that the motion of Phoebe is far from being a quasi-Keplerian motion as is the case for the planets of the solar system. Below the precession-nutation of Phoebe is determined both analytically and numerically. The differences between the results given by the two methods because of the non-Keplerian motion of Phoebe will be significantly increased when compared with the equivalent study for the Earth (Kinoshita 1977). All Phoebe’s parameters involved in the Eqs. (11) and (12), which enable one to calculate the precession and nutation are known except for the obliquity of Phoebe at J2000.0, which represents a fundamental initial condition. It is calculated in the following section.

4. Obliquity of Phoebe at J2000.0

With a similar method as in Cottereau & Souchay (2009), the obliquity is given by $\cos {\mathit{I}}_{\mathrm{0}}\mathrm{=}\begin{array}{}\mathrm{-}0\mathrm{\to }\\ {\mathit{P}}_{\mathrm{op}}\end{array}\mathrm{·}\begin{array}{}\mathrm{-}0\mathrm{\to }\\ {\mathit{P}}_{\mathrm{p}}\end{array}\mathit{,}$(22)where P_op and P_p are the unit vectors directed toward the orbital pole and the pole of rotation of Phoebe at t₀, which are: $\begin{array}{ccc}\begin{array}{}\mathrm{-}0\mathrm{\to }\\ {\mathit{P}}_{\mathrm{op}}\end{array}\mathrm{=}\left(\begin{array}{c}\\ & \\ & \\ \end{array}\right)\mathit{,}\begin{array}{}\mathrm{-}0\mathrm{\to }\\ {\mathit{P}}_{\mathrm{p}}\end{array}\mathrm{=}\left(\begin{array}{c}\\ & \\ & \\ \end{array}\right)& & \end{array}$(23)i₀, Ω₀ are the inclination and the longitude of the ascending node of Phoebe in the Saturn equatorial reference frame at J2000.0 given by the ephemerides of Emelyanov (2007) and , are the geo-equatorial coordinates of the north pole of Phoebe given by Seidelmann et al. (2007). The numerical values of these variables are given in the Appendix. Applying two rotations to , the two unit vectors are expressed in the Saturn equatorial reference frame and the obliquity at J2000.0 can be determined directly from Eq. (22).

After computation and taking into account the Phoebe retrograd rotation the numerical value of the obliquity at J2000.0 is ${\mathit{I}}_{\mathrm{0}}\mathrm{=}\mathrm{23}\mathrm{°̣}\mathrm{95.}$(24)Note that this obliquity is close to the Earth’s (). Knowing the obliquity at J2000.0, the precession-nutation of Phoebe can be determined with the help of Eqs. (11) and (12) in the following sections.

5. Numerical results of the precession-nutation of Phoebe

With the ephemeris of Emelyanov (2007), the precession and the nutation in longitude and in obliquity are determined by numerical integration with an extrapolation-algorithm based on the explicit midpoint rule. We remind here that the precession and the nutation in longitude and in obliquity are defined as the variation in space of the angular momentum axis with respect to the reference axis caused by Saturn’s gravitational perturbation. They correspond to the linear part of h and the quasi-periodic part of h and I. Before integration, the planeto-equatorial coordinates (X,Y,Z) taken from the ephemerides are converted into rectangular coordinates with respect to the orbital plane, which at their turn are converted into the spherical coordinates used in Eqs. (11) and (12). The values of the dynamical flattening and the triaxiality of Phoebe are taken from Aleshkina et al. (2010). Figures 7 and 8 show Phoebe’s precession-nutation in longitude and Phoebe’s nutation alone for a 9000 days’ time span (24.6 y). Figure 9 shows the nutation in obliquity for a 2000 days’ time span (5.47 y). Note that 9000 days correspond to the limit of the ephemerides. The 2000 days’ time span of the nutation in obliquity is chosen in order to clearly see the leading oscillations. The precession and the leading oscillations of the nutation in longitude and in obliquity are determined by a linear regression and a FFT. Our value obtained for the precession of Phoebe is $\mathrm{\Psi }\mathrm{=}\mathrm{5580.65}\mathit{t,}$(25)where t is counted in Julian centuries. Note that this value for the precession in longitude is very close to the precession for the Earth (cy). So the effect of the tidal torque of Saturn on the precession of Phoebe is roughly the same as the combined gravitational effect of the Moon and Sun on the precession of the Earth. The physical dissymmetry (high values of the dynamical flattening and of the triaxiality) of Phoebe which directly increase the amplitude of the precession, are compensated by its slow revolution and fast rotation (see Eq. (15)) .

The nutation in longitude is given by $\begin{array}{ccc}\mathrm{\Delta \Psi }& \mathrm{=}& \begin{array}{}\sum \end{array}{\mathrm{\Phi }}_{\mathit{i}}\sin \mathrm{\left(}{\mathit{\nu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\mathrm{+}{\mathrm{\Phi }}_{\mathit{i}}^{\mathrm{\prime }}\cos \mathrm{\left(}{\mathit{\nu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\mathrm{+}{\mathrm{\Theta }}_{\mathit{i}}\mathit{t}\sin \mathrm{\left(}{\mathit{\nu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\\ & \mathrm{+}& {\mathrm{\Theta }}_{\mathit{i}}^{\mathrm{\prime }}\mathit{t}\cos \mathrm{\left(}{\mathit{\nu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\mathit{,}\end{array}$(26)where and ν_i are the leading amplitudes and the frequencies of the sinusoidal terms and can be found in Table 3. Similarly, the nutation in obliquity is given by $\begin{array}{ccc}\mathrm{\Delta }\mathit{I}& \mathrm{=}& \begin{array}{}\sum \end{array}{\mathrm{\Gamma }}_{\mathit{i}}\sin \mathrm{\left(}{\mathit{\mu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\mathrm{+}{\mathrm{\Gamma }}_{\mathit{i}}^{\mathrm{\prime }}\cos \mathrm{\left(}{\mathit{\mu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\mathrm{+}{\mathrm{\Pi }}_{\mathit{i}}\mathit{t}\sin \mathrm{\left(}{\mathit{\mu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\\ & \mathrm{+}& {\mathrm{\Pi }}_{\mathit{i}}^{\mathrm{\prime }}\mathit{t}\cos \mathrm{\left(}{\mathit{\mu }}_{\mathit{i}}\mathit{t}\mathrm{\right)}\mathit{,}\end{array}$(27)where and μ_i are given in Table 4. The corresponding arguments for a combination of Phoebe’s mean longitude and mean anomaly are determined empirically and given when they have been clearly identified. The residuals after substraction of Eqs. (26) and (27) to the results of the numerical integration are the flat curves in Figs. 8 and 9.

The nutation in obliquity is dominated by two frequencies associated with the arguments 2L_s and 2L_s + M with respective periods 275.13 d and 183.45 d. In addition to these two frequencies, the nutation in longitude also presents another leading component with argument 2L_s − M. The components that are not identified in Tables 3 and 4 those with period 497 d and 260 d, are already pointed out in Tables 1 and 2 (see Sect. 3). They probably come from the non-Keplerian motion of Phoebe. As for the precession, the nutations in longitude and in obliquity of Phoebe, which show peak to peak variations of 26′′ or 8′′ are on the same order as the nutation of the Earth (with a factor and ). Because the value of the triaxiality ishigh (Aleshkina et al. 2010), it is interesting to evaluate its effect on Phoebe’s nutation. This effect is determined by Eqs. (11) and (12) taking only into account the component depending on in Eqs. (15) and (16). Figures 10 and 11 show the nutation in longitude and in obliquity depending on the triaxiality for a 10 days’ times span beginning at 410 days. This short time span is chosen to clearly see the high frequency oscillation to Phoebe’s rotation with a period roughly half of the proper rotation of Phoebe (9 h) and the higher amplitude of the whole signal. Thus we remark that despite the large triaxiality of Phoebe, its effect on the nutation in longitude and in obliquity, at the level of a few 10^-3 arcsec amplitude, is very small. As explained already for the precession, the effect of the disymmetry of Phoebe, which is characterized by a very high value of its triaxiality, is compensated by its very fast rotation, which strongly decreases the amplitude of the signal after integration.

The goal now is to determine an analytical model to describe the precession-nutation of Phoebe and to compare it with our numerical results. Below the effect of the triaxiality on the nutation will be ignored and the in Eqs. (13) and (14) are replaced by .

6. Construction of an analytical model for the precession-nutation of Phoebe

Many authors have made an analytical model of Earth’s rotation in order to understand the causes of the observed motion and to be able to predict it (Kinoshita 1977; Bretagnon et al. 1997; Roosbeek & Dehant 1997; Souchay et al. 1999). In this approach, it is interesting to build an analytical model based on the model of Kinoshita for the Earth to determine the precession and the nutation of Phoebe.

To solve analytically Eqs. (11) and (12), Kinoshita (1977) developed and as a function of time, through the variables M and L_s taking the eccentricity as a small parameter: $\begin{array}{ccc}\frac{\mathrm{1}}{\mathrm{2}}{\left(\frac{\mathit{a}}{\mathit{r}}\right)}^{\mathrm{3}}& \mathrm{=}& \mathit{f}\mathrm{\left(}\mathit{M,e}\mathrm{\right)}\mathrm{+}\mathit{o}\mathrm{\left(}{\mathit{e}}^{\mathrm{3}}\mathrm{\right)}\\ {\left(\frac{\mathit{a}}{\mathit{r}}\right)}^{\mathrm{3}}\cos \mathrm{2}\mathit{\lambda }& \mathrm{=}& \mathit{g}\mathrm{\left(}\mathit{M,}{\mathit{L}}_{\mathit{s}}\mathit{,e}\mathrm{\right)}\mathrm{+}\mathit{o}\mathrm{\left(}{\mathit{e}}^{\mathrm{3}}\mathrm{\right)}\end{array}$(28)The detailed demonstration is presented in Kinoshita (1977). Because the eccentricity of the Earth is small, truncating the development at the 3rd order is a sufficient approximation. But the eccentricity of Phoebe is higher (by one order) than the eccentricity of the Earth (0.18 instead 0.0167). Moreover, as we have seen before, Phoebe’s orbital elements are undergoing intensive variations. Then the relative error made by replacing and by the developments such as Eq. (28) need to be quantified. To that purpose, and are developed as for the Earth, but with the real values of e,M,L_S taken in the ephemerides of Emelyanov (2007), instead of a linear value (for e,M ans L_s) for the Earth.

6.1. Test for the developments of and

In this section, the developments of and done by Kinoshita (1977) for the Earth and by Cottereau & Souchay (2009) for Venus are tested for Phoebe. We used the same methods as in Cottereau & Souchay (2009) to determine the development of . As the eccentricity of Phoebe is high, the development at the 3rd order as done by Cottereau & Souchay (2009) for Venus is insufficient because it yields relative errors on the order of 4 × 10^-3. This is because of the slow convergence of the development, which must be carried out up to 6th order in our case, if we want to reach a relative 10^-4 accuracy. Figure 12 shows the discrepancies between the numerical and the semi-analytical results using this developments both to the 3rd and 6th order. This shows that high order developments are needed when taking into account highly eccentric orbits as Phoebe’s one.

For the developments of , the same conclusion is obtained. Figure 13 shows the discrepancies between the numerical and the semi-analytical results using the developments of , to 3rd and 6th order. The difference descreases from 0$\mathrm{″̣}$3 to 0$\mathrm{″̣}$01 (peak to peak). The developments are given by $\begin{array}{ccc}\frac{\mathrm{1}}{\mathrm{2}}{\left(\frac{\mathit{a}}{\mathit{r}}\right)}^{\mathrm{3}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}\mathrm{+}\frac{\mathrm{3}}{\mathrm{4}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}\frac{\mathrm{15}}{\mathrm{16}}{\mathit{e}}^{\mathrm{4}}\mathrm{+}\frac{\mathrm{35}}{\mathrm{32}}{\mathit{e}}^{\mathrm{6}}& & \\ \mathrm{+}\left(\frac{\mathrm{3}}{\mathrm{2}}\mathit{e}\mathrm{+}\frac{\mathrm{27}}{\mathrm{16}}{\mathit{e}}^{\mathrm{3}}\mathrm{+}\frac{\mathrm{261}}{\mathrm{128}}{\mathit{e}}^{\mathrm{5}}\right)\cos \mathit{M}& & \\ \mathrm{+}\left(\frac{\mathrm{9}}{\mathrm{4}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}\frac{\mathrm{7}}{\mathrm{4}}{\mathit{e}}^{\mathrm{4}}\mathrm{+}\frac{\mathrm{141}}{\mathrm{64}}{\mathit{e}}^{\mathrm{6}}\right)\cos \mathrm{2}\mathit{M}& & \\ \mathrm{+}\left(\frac{\mathrm{53}}{\mathrm{16}}{\mathit{e}}^{\mathrm{3}}\mathrm{+}\frac{\mathrm{393}}{\mathrm{256}}{\mathit{e}}^{\mathrm{5}}\right)\cos \mathrm{3}\mathit{M}& & \\ \mathrm{+}\left(\frac{\mathrm{77}}{\mathrm{16}}{\mathit{e}}^{\mathrm{4}}\mathrm{+}\frac{\mathrm{129}}{\mathrm{160}}{\mathit{e}}^{\mathrm{6}}\right)\cos \mathrm{4}\mathit{M}\mathrm{+}\frac{\mathrm{1773}}{\mathrm{256}}{\mathit{e}}^{\mathrm{5}}\cos \mathrm{5}\mathit{M}& & \\ \mathrm{+}\frac{\mathrm{3167}}{\mathrm{320}}{\mathit{e}}^{\mathrm{6}}\cos \mathrm{6}\mathit{M}\mathrm{+}\mathit{O}\mathrm{\left(}{\mathit{e}}^{\mathrm{7}}\mathrm{\right)}& & \\ {\left(\frac{\mathit{a}}{\mathit{r}}\right)}^{\mathrm{3}}\cos \mathrm{2}\mathrm{(}\mathit{\lambda }\mathrm{-}\mathit{h}\mathrm{)}\mathrm{=}\left(\mathrm{1}\mathrm{-}\frac{\mathrm{5}}{\mathrm{2}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}\frac{\mathrm{13}}{\mathrm{16}}{\mathit{e}}^{\mathrm{4}}\mathrm{+}\frac{\mathrm{35}}{\mathrm{288}}{\mathit{e}}^{\mathrm{6}}\right)\cos \mathrm{2}{\mathit{L}}_{\mathit{S}}& & \\ \mathrm{+}\frac{\mathrm{4}}{\mathrm{45}}{\mathit{e}}^{\mathrm{6}}\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{-}\mathrm{6}\mathit{M}\mathrm{)}\mathrm{+}\frac{\mathrm{81}}{\mathrm{1280}}{\mathit{e}}^{\mathrm{5}}\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{-}\mathrm{5}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left({\mathit{e}}^{\mathrm{4}}\frac{\mathrm{1}}{\mathrm{24}}\mathrm{+}\frac{\mathrm{7}}{\mathrm{240}}{\mathit{e}}^{\mathrm{6}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{-}\mathrm{4}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left({\mathit{e}}^{\mathrm{3}}\frac{\mathrm{1}}{\mathrm{48}}\mathrm{+}\frac{\mathrm{11}}{\mathrm{768}}{\mathit{e}}^{\mathrm{5}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{-}\mathrm{3}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left(\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}{\mathit{e}}^{\mathrm{3}}\frac{\mathrm{1}}{\mathrm{16}}\mathrm{-}\frac{\mathrm{5}}{\mathrm{384}}{\mathit{e}}^{\mathrm{5}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{-}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left(\frac{\mathrm{7}}{\mathrm{2}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}{\mathit{e}}^{\mathrm{3}}\frac{\mathrm{123}}{\mathrm{16}}\mathrm{+}\frac{\mathrm{489}}{\mathrm{128}}{\mathit{e}}^{\mathrm{5}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{+}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left(\frac{\mathrm{17}}{\mathrm{2}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}{\mathit{e}}^{\mathrm{4}}\frac{\mathrm{115}}{\mathrm{6}}\mathrm{+}\frac{\mathrm{601}}{\mathrm{48}}{\mathit{e}}^{\mathrm{6}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{+}\mathrm{2}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left(\frac{\mathrm{845}}{\mathrm{48}}{\mathit{e}}^{\mathrm{3}}\mathrm{-}{\mathit{e}}^{\mathrm{5}}\frac{\mathrm{32525}}{\mathrm{768}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{+}\mathrm{3}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\left(\frac{\mathrm{533}}{\mathrm{16}}{\mathit{e}}^{\mathrm{4}}\mathrm{-}{\mathit{e}}^{\mathrm{6}}\frac{\mathrm{13827}}{\mathrm{160}}\right)\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{+}\mathrm{4}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\frac{\mathrm{228347}}{\mathrm{3840}}{\mathit{e}}^{\mathrm{5}}\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{+}\mathrm{5}\mathit{M}\mathrm{)}& & \\ \mathrm{+}\frac{\mathrm{73369}}{\mathrm{720}}{\mathit{e}}^{\mathrm{6}}\cos \mathrm{(}\mathrm{2}{\mathit{L}}_{\mathit{S}}\mathrm{+}\mathrm{6}\mathit{M}\mathrm{)}\mathrm{+}\mathrm{O}\mathrm{\left(}{\mathit{e}}^{\mathrm{7}}\mathrm{\right)}\mathit{.}& & \end{array}$

6.2. Analytical model of Phoebe

To express and by purely analytical functions, Kinoshita replaced e,M,L_S by their linear mean values. For Venus, the same method is used in Cottereau & Souchay (2009), who replaced e,M,L_S by their linear mean values given by Simon et al. (1994). Using the mean elements of Sect. (2) and the developments given by Eqs. (29) and (30), an analytical model is obtained. As explained in Kinoshita & Souchay (1990) and Cottereau et al. (2010) for the Earth and for Venus, the difference between the nutation computed by analytical tables and the numerical integration is small (on the order of 10^-5) and caused by the indirect planetary effects on the planet. For Phoebe, this difference is also caused by indirect perturbations, but cannot be ignored as will be shown below. Replacing the semi-major axis in Eqs. (11) and (12) by (Sect. 3), Δh and ΔI depend directly on the scaling factor K_S: ${\mathit{K}}_{\mathit{S}}\mathrm{=}\frac{\mathrm{3}{\mathit{n}}^{\mathrm{2}}}{\mathit{\omega }}\frac{\mathrm{2}\mathit{C}\mathrm{-}\mathit{A}\mathrm{-}\mathit{B}}{\mathrm{2}\mathit{C}}\mathrm{=}\mathrm{-}\mathrm{11674}\mathrm{″̣}\mathrm{4857} \mathrm{Julian}\mathit{.}\mathrm{cy}\mathit{,}$(31)where is the dynamical flattening of Phoebe, ω the angular velocity of Phoebe and n the mean motion of Phoebe given by ${\mathit{n}}^{\mathrm{2}}\mathit{a̅}\mathrm{=}\mathit{G}{\mathit{M}}^{\mathrm{\prime }}$(32)The numerical value is directly deduced from the values of the A,B,C and ω taken from Aleshkina et al. (2010) and Bauer et al. (2004). Note that this scaling factor is very close to Earth’s (11036$\mathrm{″̣}$6 by Julian cy). As seen in Sect. 2 the precession of Phoebe corresponds to the polynomial part of Eq. (11) and can be directly determined by $\begin{array}{ccc}\mathit{\psi }& \mathrm{=}& \mathrm{-}\mathrm{\left(}{\mathit{K}}_{\mathit{s}}\cos \mathit{I}\mathrm{\int }\left(\frac{\mathrm{1}}{\mathrm{2}}\mathrm{+}\frac{\mathrm{3}}{\mathrm{4}}{\mathit{e}}^{\mathrm{2}}\mathrm{+}\frac{\mathrm{15}}{\mathrm{16}}{\mathit{e}}^{\mathrm{4}}\mathrm{+}\frac{\mathrm{35}}{\mathrm{32}}{\mathit{e}}^{\mathrm{6}}\right)\mathrm{d}\mathit{t}\mathrm{\right)}\\ & \mathrm{=}& \mathrm{5541.18}\mathit{t,}\end{array}$(33)where e is replaced by the linear part of Eq. (19). The discrepancies between the precession as computed by analytical integration and that given by Eq. (25) are on the order of 1% and show that the analytical model describes accuratily the precession of Phoebe. Replacing e,M,L_S in the developments (29) and (30) by the linear part of Eqs. (19)–(21), the nutation of Phoebe is determined by analytical integration. Figures 14 and 15 show the differences between the nutation in longitude and in obliquity of the angular momentum axis from a purely analytical formulation and numerical integration on 9000 days’ time span. These differences are compared to the residuals obtained in Sect. 5 after substraction Eqs. (26) and (27) determined by FFT to the numerical signal (central curves in the Figs. 14 and 15).

Figures 14 and 15 show that the semi-analytical functions characterizing the nutation of Phoebe in longitude and in obliquity determined thanks to a FFT are more accurate than those obtained by a pure analytically integration (with linear values of e,M,L_s). The purely analytical model restricts the expressions of the functions that describe the orbital elements. For the Earth (Kinoshita & Souchay 1990) and for Venus (Cottereau et al. 2010), this model is very satisfying because their motions are close to Keplerian and consequently the orbital parameters can be replaced by linear functions. Because the motion of Phoebe is highly disturbedcompared to the planets of the solar system, the curves of the temporal variations of the orbital elements (Emelyanov 2007) cannot be fitted by simple functions as shown by the large variations around the mean elements in Sect. 3. Thus the analytical model does not describe the periodic components due to the indirect effect of the Sun or of the other satellites of Saturn. In contrast these effects are directly taken into account when using the functions given by Eqs. (26) and (27) because the FFT approach followed in Sect. 5 is better suited to describe quasi-periodic motion. However, the precession is less sensitive to these periodic variations, which are averaged during the calculation. This explains the similarity of the precession rate given by the two models. To conclude, the purely analytical model set by Kinoshita (1977) for the Earth gives a good first approximation of the nutation of Phoebe but a more sophisticated analytical model is needed if one wants a precision down to the 10^-5 level. We state here that this accuracy is iherent to our model and is independant of the uncertainties in the observations. This shows that we are limited only by the uncertainties on the input data.

7. Conclusion

The purpose of this paper was to calculate for the first time the combined motion of precession and nutation of Saturn’s satellite Phoebe.

One of the important steps in our study was to describe the orbital motion of Phoebe by fitting the curves of the temporal variations of the orbital elements a,e,M and L_s with polynomial functions. As shown in Fig. 2, the periodic variations of the semi-major axis around its mean value with a relative amplitude of the order on 10^-3 are clearly due to the perturbing effects of other celestial bodies like the Sun. The orbital motion is far from be a Keplerian one, as shown by the large polynomial expressions of e,M and L_s, and by the large sinusoidal amplitudes characterizing the residuals after substraction of these polynomials. Applying the theoritical framework already used by Kinoshita (1977) for the Earth, the precession and the nutation motion of Phoebe are determined both analytically and from numerical integration of the equations of motion. We found that the precession-nutation motion of Phoebe undergoing the gravitational perturbation of Saturn is quite similar to the Earth undergoing the gravitational effect of both the Moon and the Sun. Thus our value for the precession of Phoebe, that is to say 5580$\mathrm{″̣}$65 cy, is very close to the corresponding value for the Earth (5081"/cy) and the nutation in longitude and in obliquity of Phoebe with peak to peak variations of 26′′ and 8′′ are of the same order on amplitude as the nutation of the Earth (36′′ and 18′′ peak to peak). Moreover Phoebe’s obliquity () is roughly the same as the Earth’s (). Note that the physical dissymmetry (high value of the dynamical flattening and of the triaxiality) and the high eccentricity of Phoebe, which direclty increases the amplitude of the precession and the nutation, is compensated by its slow revolution and fast rotation.

In Cottereau & Souchay (2009) we have demonstrated that the effects of the large triaxiality of Venus on its nutation are on the same order of magnitude as the effects of the dynamical flattening. By contrast, the effects of the large triaxiality of Phoebe on its nutation are compensated by its fast rotation and decrease the of amplitude, which becomes negligible compared to the dynamical flattening.

We also investigated the possibility to construct analytical tables of Phoebe nutation, as was done for the Earth. Because Phoebe has a high eccentricity, the analytical developments at the 3th order given by Kinoshita (1977) for the Earth are not valid. We must carry out the developments up to 6th order to reach a relative 10^-4 accuracy.

Thus although the amplitude of the nutation motion is close to the Earth’s, we demonstrated that the analytical model used by Kinoshita for the Earth does not describe the nutation motion of a disturbed body like Phoebe with the same accuracy. This analytical model does not take into account the large perturbing effects of the celestial bodies on the orbit of the satellite.

To describe the nutation motion of Phoebe a FFT approach is better than a pure analytical integration done with linear expressions for e,M and L_S as shown Figs. 14 and 15. The FFT analysis is better able to describe the periodic variations characterizing the nutation signals of Phoebe. These periodic variations are averaged during the calculation of the precession, which lessens the discrepancies between the two models.

To conclude, we have shown that the analytical model set by Kinoshita (1977) gives a good first approximation of the precession-nutation of Phoebe, but further analytical developments are needed to reach the same accuracy as for the terrestrial planets.

This work can be a starting point for further studies such as the study of another very precise analytical model of the rotation of Phoebe by taking into account effects ignored in this paper, as for instance the direct effects of the Sun, of Titan and of Saturn’s dynamical flattening. Such a model is required to develop the long’s term ephemerides of Phoebe’s rotation, which should require long’s term orbital ephemerides, which are not yet available.

References

Appendix A Appendix A

All Tables

All Figures