Librational response of a deformed 3-layer Titan perturbed by non-keplerian orbit and atmospheric couplings

Librational response of a deformed 3-layer Titan perturbed by non-keplerian orbit and atmospheric couplings

A. Richard N. Rambaux B. Charnay IMCCE, Observatoire de Paris, CNRS UMR 8028
77 avenue Denfert-Rochereau, 75014 Paris, France
Université Pierre et Marie Curie, UPMC - Paris 06 Laboratoire de Météorologie Dynamique, Paris

The analyses of Titan’s gravity field obtained by Cassini space mission suggest the presence of an internal ocean beneath its icy surface. The characterization of the geophysical parameters of the icy shell and the ocean is important to constrain the evolution models of Titan. The knowledge of the librations, that are periodic oscillations around a uniform rotational motion, can bring piece of information on the interior parameters.

The objective of this paper is to study the librational response in longitude from an analytical approach for Titan composed of a deep atmosphere, an elastic icy shell, an internal ocean, and an elastic rocky core perturbed by the gravitational interactions with Saturn. We start from the librational equations developed for a rigid satellite in synchronous spin-orbit resonance. We introduce explicitly the atmospheric torque acting on the surface computed from the Titan IPSL GCM (Institut Pierre Simon Laplace General Circulation Model) and the periodic deformations of elastic solid layers due to the tides. We investigate the librational response for various interior models in order to compare and to identify the influence of the geophysical parameters and the impact of the elasticity.

The main librations arise at two well-separated forcing frequency ranges: low forcing frequencies dominated by the Saturnian annual and semi-annual frequencies, and a high forcing frequency regime dominated by Titan’s orbital frequency around Saturn. At low forcing frequency, the librational response is dominated by the Saturnian gravitational torque and the atmospheric torque has a small effect. In addition, the libration amplitude in that case is almost equal to the magnitude of the perturbation. The modulation of the gravitational torque amplitude at the orbital frequency with periodic deformation induces long-period terms in the librational response which contain information on the internal structure. At high forcing frequency the libration depends on the inertia of the layers and the elasticity can strongly reduce its amplitude at orbital frequency. For example, the amplitude of diurnal libration for oceanic models goes from about meters if the icy shell is purely rigid to meters when the elasticity is included, i.e. a reduction of about . For models without ocean, diurnal libration goes from meters in a rigid case to meters for an elastic case, a very low reduction due to the weak deformation of an entirely solid satellite compared to the deformation of a thin icy shell. Oceanic models with elastic solid layers have the same order of libration amplitude than the oceanless models, which makes more challenging to differentiate them by the interpretation of librational motion.

Titan, libration, elasticity, dynamic, orbital perturbation, atmospheric torque
journal: Planetary and Space Science\biboptions



Titan exhibits a rich interior structure due to its large mean radius of 2575 kilometers. The recent measurements of the gravity field by Iess et al. (2010, 2012) reveal that Titan’s moment of inertia (MoI) is as low as . A large panel of internal structures, made of a low density core surrounded by icy layers, are consistent with this range of MoI and the mean density value of kg m deduced by Iess et al. (2010). In addition, the presence of a global internal ocean has been suggested from several techniques. First, Lorenz et al. (2008) deduced the putative internal ocean through the response of the rotational motion of the surface to the atmospheric coupling. A second piece of evidence has been obtained from the determination of the non-zero obliquity by Stiles et al. (2008) from the radar images of Cassini. They determined an obliquity equal to 0.3 degree that is larger than the degree obtained for a rigid Titan in a Cassini state (e.g. Bills and Nimmo (2008)). Then Bills and Nimmo (2008) and Baland et al. (2011) suggested that this deviation is related to the influence of an internal ocean on the obliquity. Additional arguments in favor of this internal ocean have been obtained by electrical (Béghin et al., 2012), topographical (Nimmo and Bills, 2010) and gravitational (Iess et al., 2012) analyses of Titan from Cassini-Huygens data.

Titan’s rotational motion has been measured by Stiles et al. (2008, 2010) that used the Cassini spacecraft’s radar images in order to follow landmarks at the surface. They obtained a near-synchronous rotation for Titan with a drift rate of degree per day, i.e. degree per year (Stiles et al., 2010). This approach has been recently revisited by Meriggiola and Iess (2012) that determined a synchronous spin-orbit motion within a residual of about 0.02 degree per year. The main advantage of their approach is the introduction of Titan’s librational motion in the reduction process. The librations describe the oscillations around the uniform rotational motion. Here we focus on the librations in longitude that correspond to the oscillations of the body principal axis projected onto the equatorial plane of the satellite. In the case of Titan, they have two distinct origins. The first one comes from the gravitational torque exerted by Saturn on the dynamical figure of Titan. The second one results from the coupling between the surface and the dense atmosphere of Titan, the atmospheric torque. Their amplitudes are modified by the interior parameters of Titan. However, Goldreich and Mitchell (2010) suggested that the elastic behavior of the icy shell is not negligible and may strongly reduce the librational motion. Such prediction has been confirmed by Van Hoolst et al. (2013) that computed the libration in longitude at the orbital frequency of Titan. In that case, Titan’s surface will deform instead of rotate since the ocean figure should always point toward the planet and exert an elastic torque on the icy shell.

The objective of this paper is to determine the librational reponse of an elastic Titan at various forcing frequencies resulting from its orbital motion. By including the atmospheric coupling, we want to decorrelate surface forcing from the internal geophysical properties related to the internal ocean and perturbation of Titan’s orbit. We investigate the wide spectrum of Titan’s librations on contrary to Van Hoolst et al. (2009, 2013) that focused on the diurnal frequency (Titan’s orbital frequency) and its harmonics, and on the Saturnian semi-annual frequency (twice Saturn’s orbital frequency) for the atmospheric coupling. The main interest is that in the orbital motion, there are librations at the Saturnian semi-annual frequency that comes from the interaction of Saturn with the Sun. Such librations have the same frequency than the main component of the atmospheric torque. In addition, the orbital frequency spectrum can be in, or close to, resonance with some proper frequencies as shown for the Galilean satellites (Rambaux et al., 2011). Finally, the periodic variation of the gravitational torque amplitude at orbital frequency provides long-period terms in the libration with amplitudes dependent on the interior model.

In this paper, the recent interior models of Titan are also used (e.g. Fortes (2012); Castillo-Rogez and Lunine (2010); McKinnon and Bland (2011)) in order to compute the values of the proper frequencies and to discuss the differences in the amplitude of librations. The elasticity is investigated by computing the radial deformations of surfaces due to the tides and responsible for the gravitational torques amplitude variations.

Finally, we use the 3D atmospheric model from Lebonnois et al. (2012) that predicts an atmospheric torque smaller than in the Tokano and Neubauer (2005) paper used in all previous studies (Lorenz et al., 2008; Karatekin et al., 2008; Van Hoolst et al., 2009, 2013; Goldreich and Mitchell, 2010).

In the first part of the paper, the internal structure models selected for the librational computation are described. The properties of these models depend on the history and energy sources of Titan (e.g. Tobie et al. (2012)). We selected a broad range of possible internal scenarios for Titan in order to characterize the impact of the geophysical parameters on the librations. The atmospheric torques of Charnay and Lebonnois (2012) (called here CH12) and of Tokano and Neubauer (2005) (called here TO05) are described and discussed in Sect. 2. The orbit of Titan is then analyzed in Sect. 3 by using the frequency analysis method providing the spectrum of the orbital motion. In Sect. 4, the elasticity is introduced in the librational equations through periodic torques amplitudes and the librations are analytically determined for the rigid and elastic cases. Finally, the behavior of the libration angle at different forcing frequency ranges is analyzed and the influence of the geophysical parameters of the different interior models is discussed in Sect. 5.

1 Interior models

Titan’s internal structure has been revealed by accurately tracking the trajectory of Cassini spacecraft approaching Titan during six flybys (Iess et al., 2010, 2012). They measured the gravity field and its variations allowing to infer information on the density profile of the satellite. Since the inversion between the gravity field and density profile is not unique, only a range of models can be determined. The add-on assumption that Titan is close to the hydrostatic equilibrium led Iess et al. (2010) to obtain an estimation of the moment of inertia (MoI) of Titan between and (where and are the mass and mean radius of Titan, respectively). Such small MoI implies an increase of density towards Titan’s center and requires a low-density core to match the mean density of kg m (Iess et al., 2010).

In parallel to these gravity measurements, two categories of thermal and chemical models have been developed to determine Titan’s internal structure (e.g. Castillo-Rogez and Lunine (2010, 2012); McKinnon and Bland (2011); Fortes (2012) and see the review of Tobie et al. (2012)). In the first category, as developed by Fortes (2012), the models described essentially the upper layers (dense and light ocean, comparison of solid models with pure water ice or methane clathrate layers), while in the second category the models focused on the inner core composition assuming the presence of a global ocean layer (e.g. Castillo-Rogez and Lunine (2012); McKinnon and Bland (2011)).

Castillo-Rogez and Lunine (2012) built interior models made of anhydrous silicate core surrounded by hydrated rock and water ice while McKinnon and Bland (2011) focused on hydrated silicates surrounded by mixture of rock and ice. These models of inner core allow the presence of a small iron core, and are compatible with an internal ocean as deduced by Iess et al. (2012).

Six different models reported in Table 1 have been selected, coming from Fortes (2012) (the models called F1, F2 and the solid model F3), Castillo-Rogez and Lunine (2010) (models CA10 and FE10) and McKinnon and Bland (2011) (model MC11). The FE10 model is similar to the CA10 model with an additional small inner iron core as suggested by Castillo-Rogez and Lunine (2010). The icy shell thickness has been taken equal to kilometers for each model. Such value corresponds to the upper bound obtained by topographical model of (Nimmo and Bills, 2010). The lower bound has been obtained by Béghin et al. (2012) by using the Schumann’s resonance in Titan’s atmosphere. The icy shell thickness is an essential parameter for Titan’s models with a rigid shell (Van Hoolst et al., 2009) but the elasticity strongly diminish its influence (see Section 5). Fortes (2012) used oceans with a bottom mean radius of kilometers (models F1 and F2). To be able to compare the influence of the ocean density and the inner core structure on the libration, the same ocean bottom mean radius is used for the CA10, MC11 and FE10 models. For a given set of solid layers sizes and densities, the ocean density is then adjusted to conserve Titan’s MoI and mass.

F1 F2 F3
Layer (kg m) R (km) (kg m) R (km) (kg m) R (km)
Icy shell 930.9 2575.5 930.9 2575.5 932.8 2575.5
Ocean 1023.5 2475 1281.3 2475 - -
HP ice 1 - - - - 1193.3 2429
HP ice mantle 1272.7 2225 1350.9 2225 1268.7 2304
Silicate mantle 1338.9 2163 - - 1338.9 2185
Silicate core 2542.3 2116 2650.4 1984 2584.1 2055
(MR) 0.3414 0.3413 0.3413
(MR) 0.0357 0.0357 -
(MR) 0.2320 0.2133 -

CA10 MC11 FE10
Layer (kg m) R (km) (kg m) R (km) (kg m) R (km)
Icy shell 932.8 2575.5 932.8 2575.5 932.8 2575.5
Ocean 1317.6 2475 1315.5 2475 1224.1 2475
HP ice mantle 1350.9 2225 1350.9 2225 1350.9 2225
Rocky mantle 2520 2000 2200 2075 2520 2000
Silicate core 3400 900 3300 1300 3400 900
Iron core - - - - 6500 500
(MR) 0.3402 0.3378 0.3336
(MR) 0.0358 0.0358 0.0358
(MR) 0.2095 0.2073 0.2096
Table 1: Characteristics of the internal structure models selected for this paper. The first three models quoted Light Ocean (F1), Dense Ocean (F2) and Pure-water-ice (F3) models are from Fortes (2012). Model CA10 is from Castillo-Rogez and Lunine (2010) and MC11 from McKinnon and Bland (2011), both with dense ocean and high pressure ice layers. The last model (FE10) is a model similar to CA10 including an iron core. is the total moment of inertia of the body, while and are polar moments of inertia of the icy shell and solid interior respectively, expressed in terms of mass and radius of Titan ( kg and km following Fortes (2012)).

2 Atmospheric torque

Titan has a thick atmosphere extending up to 800 km with a surface pressure of nearly 1.5 bars. Exchanges of angular momentum happen between the atmosphere and the surface, producing an atmospheric torque influencing Titan’s rotation. The atmosphere dynamic is mainly driven by insolation, which creates a circulation zone called a Hadley cell between hotter and cooler regions. The circulation in Titan’s troposphere is essentially dominated by one Hadley cell extending from one pole to the other. At the equinox, the circulation reverses with the formation of two Hadley cells, and the ascending air zone (also called the Intertropical Convergence Zone) moves from one pole to the other. This reversal produces a seasonal angular momentum exchange with Saturn’s semi-annual period ( terrestrial days) and has been suggested as the main torque influencing Titan’s rotation (Tokano and Neubauer, 2005). The total atmospheric angular momentum (AAM) is given by (Tokano and Neubauer, 2005)


where km is Titan’s mean radius, ms is Titan’s surface gravity, is the surface pressure, is the zonal wind speed (i.e. the horizontal wind speed in the west-east direction and relative to the surface, is positive for eastward wind), is the latitude and is the pressure.

The atmospheric torque is then


The zonal wind speed of an air parcel in the ascending zone is approximately null due to a strong surface friction. Its angular momentum is conserved during its meridional transport from the ascending zone to the poles. Thus, the AAM is maximum at equinoxes (i.e. when the ascending zone is the farthest from Titan’s rotation axis), and minimum at solstices (i.e. when the ascending zone is the nearest to Titan’s rotation axis) (Tokano and Neubauer, 2005). This change corresponds to an exchange of angular momentum with the surface by friction of winds.

The angular momentum exchange, and so the atmospheric torque, depends essentially on the extension of the Hadley cell in latitude and altitude (Mitchell, 2009). All previous studies on Titan’s non-synchronous rotation (Lorenz et al., 2008; Karatekin et al., 2008; Van Hoolst et al., 2009, 2013; Goldreich and Mitchell, 2010) used the atmospheric torque at Saturnian semi-annual frequency from the general circulation model (GCM) of Tokano and Neubauer (2005) ( kg m s; this torque is called here TO05). In this study, the atmospheric torque from the Titan IPSL GCM (Institut Pierre-Simon Laplace General Circulation Model) (Lebonnois et al., 2012) is also used. This model successfully reproduces the winds and the thermal structure observed by Cassini and Huygens (Lebonnois et al., 2012; Charnay and Lebonnois, 2012). Charnay and Lebonnois (2012) show that the Hadley cell is essentially trapped in the first two kilometers. This trapping, in addition to a lower latitudinal extension of the cell, reduces the angular momentum exchange by approximately a factor 10 compared to Tokano and Neubauer (2005).

The latitudinal extension of the Hadley cell depends a great deal on the thermal inertia of the ground (i.e. the ground thermal response to solar heating). Its value for Titan is unknown, it has been estimated to be around J m s K (Tokano, 2005) and used in Tokano and Neubauer (2005), but the value may be higher. Atmospheric torques have been calculated with a thermal inertia of and J m s K, the last one corresponding to a realistic maximal value. In that case, the Saturnian semi-annual frequency component of the Titan IPSL GCM torque (called here CH12) is of about kg m s and is reduced to kg m s for a model with a thermal inertia of J m s K. In this study, we focus on the first model providing the largest Saturnian semi-annual frequency torque amplitude.

The atmospheric torque CH12 with a thermal inertia of J m s K and the Saturnian semi-annual frequency torque TO05 are represented in Fig. 1. Many other frequencies than the Saturnian semi-annual one appear in the spectrum of the atmospheric torque CH12. The most important ones are frequencies corresponding to (linked to the eccentricity of the orbit of Saturn around the Sun), and (other harmonics from the seasonal cycle) Saturnian year with amplitudes below kg m s. The other frequencies correspond to atmospheric waves that are either free (baroclinic waves, kelvin waves,…) or forced by the thermal tides and the gravitational tides caused by Saturn (Tokano and Neubauer, 2002). The free atmospheric waves are very dependent on the model, but their amplitudes appear to be small. Concerning the thermal and gravitational tides, they have a very limited impact on the tropospheric winds and therefore on the exchange of angular momentum with the surface in the model from Lebonnois et al. (2012). Other effects may influence the atmospheric torque. Mitchell (2009) suggested the methane cycle could reduce the amplitude of the torque and produce a small drift while Tokano (2012) suggested the presence of mountains could increase the AAM a little without changing the torque very much. These effects remain poorly constrained but should only have a small impact.

Figure 1: Atmospheric torque amplitude along the Z-axis for the Titan IPSL GCM torque (dark plot) with a thermal inertia of J m s K. The amplitude of the Saturnian semi-annual frequency torque TO05 is plotted in red for comparison. Initial date is J2000.

3 Orbital forcing

The librational motion is the rotational response of Titan to the gravitational or atmospheric torques exerted on its dynamical figure. The gravitational torque depends on the relative distance between Titan and Saturn and on the true longitude of Titan (e.g. Rambaux et al. (2011)). The orbital elements are time-dependent and can be approximated by quasi-periodic Fourier series resulting from the interaction of Titan with Saturn, the Sun or Iapetus (Vienne and Duriez, 1995).

Here, we develop a linear theory of the librational motion of Titan that is approximated in small quantities, the amplitude of libration and the eccentricity. The orbital motion of Titan comes from the JPL ephemerides that contain the last observation from Cassini space mission (Giorgini et al., 1996). The ephemerides are given over the period of 07-Jan-1800 to 07-Jan-2200 corresponding to an interval of 400 years. The frequency analysis method developed by Laskar (1988, 2003) and implemented in the TRIP software (Gastineau and Laskar, 2012) is used in order to decompose the true longitude into Fourier series.

The decomposition of the true longitude is given in Table 2. The notation of Vienne and Duriez (1995) is used, where represents the mean longitude of Titan ( with the longitude of node, the mean anomaly and the argument of pericenter, see Figure 2), and the longitudes of pericenter () of Titan and Iapetus, respectively, and the longitudes of node of Titan and Iapetus, and the mean longitude of the Sun. The true longitude is given by with the true anomaly.

Figure 2: The geometry of the satellite orbit around the planet with respect to the reference plane.

The largest term of the Table 2 corresponds to the first term of the development of the difference between the true anomaly and the mean anomaly as function of the eccentricity , i.e. . The magnitude is equal to 11899.3237” corresponding to an eccentricity of 0.0288 and the frequency is equal to the orbital frequency. The following term corresponds to the second term of the development in the true anomaly, i.e. and oscillates at twice the orbital frequency. The other frequencies result from the motion of Saturn around the Sun (, , and here called Saturnian annual, semi-annual and ter-annual terms respectively, and the combination ) and interaction with Iapetus () (see Vienne and Duriez (1995)). The last term has a frequency close to the combination of () where is the mean longitude of Jupiter. The last two terms of the Table 2 with periods of and days respectively are small terms of the true longitude series with magnitudes below . These terms have small magnitudes but due to their frequency values, the librational response amplitudes are greater than the 10 meters limit which is used for the Table 4 (see Section 5).

Let us note that the frequency analysis of the orbit over years is not sufficient to identify the periods longer than several hundred years. Vienne and Duriez (1995) have identified longer periods in Titan’s orbital motion listed in the analytical ephemerides. Frequencies associated with ( years), ( years), ( years) and ( years) possess high magnitudes. Since they are the signature of time variations much larger than the 400 years of the ephemerides, these terms act like secular deviation on the true longitude. These secular terms are removed from the analysis in the determination of the linear deviation of the true longitude. However, for a longer time range, the long period terms from longitudes of nodes and pericenters have to be considered explicitly.

Freq. Period Magnitude Phase Identification
(rad/days) (days) (”) (degree)
0.394018 15.9464 11899.3237 163.3693
0.788036 7.9732 212.5868 -32.7941
0.394081 15.9439 56.6941 -68.1211
0.001169 5376.6331 43.7313 -66.0428
0.000584 10750.3648 37.5508 138.4821
0.392897 15.9919 31.5673 10.8789
0.001753 3583.9304 5.6147 250.1412
0.009810 640.4892 1.4983 -77.2905 -
Table 2: Frequency analysis of the true longitude of Titan from ephemerides. The difference between true longitude and mean anomaly is decomposed according to , with is the initial value of the satellite rotation angle measured from the line of ascending node, the magnitude, the frequency and the phase. The frequency of the last term is close to the combination of () where is the mean longitude of Jupiter. This table presents the main terms of the true longitude series that acts on the librational motion. Initial date is J2000.

4 Librational model

4.1 Rigid case

The librational equations describing the variations of Titan’s rotation are obtained from the angular momentum equation applied to each layer of the satellite. In the case of rigid solid layers, the hydrostatic equilibrium shape is composed of a static tidal bulge on which the planet exerts a gravitational forcing at each instant. In addition, the misalignment of the static tidal bulges of the shell and the inner core generates an internal gravitational torque coupling the upper layers to the inner core. Following the developments of Van Hoolst et al. (2008, 2009) and Rambaux et al. (2011), the equations governing the librational motions of the rigid triaxial shell and inner core perturbed by external and internal gravitational forces are written as


where subscripts and refer to the icy shell and to the inner core respectively, is the libration angle, the mean anomaly, is the true longitude, is the rotation angle of the satellite’s longest axis measured from the line of the ascending node and its initial value for the layer (where subscript refers to the shell or the inner core), is the mean motion. and are the amplitudes of the effective gravitational torques exerted by Saturn on the dynamical figure of each layer, while is the amplitude of the internal gravitational torque exerted by the shell on the inner core due to their misalignment. The gravitational torques amplitudes are given by and , where , and are the principal moments of inertia (defined as ) and , and are the ocean pressure effect on the solid layer expressed as increment of inertia. The internal gravitational torque amplitude also depends on the geophysical parameters of the body as described by Van Hoolst et al. (2009). In this study, the dynamical equations are linearized as function of the eccentricity and libration amplitude. The frequencies of the forcing torques acting on the two solid layers are given by the Fourier series of developed in the previous section.

The proper frequencies (or natural frequencies) of the satellite are described in terms of free frequencies of the harmonic oscillator as done by Dumberry (2011). First, the free frequencies corresponding to the system (3) free of internal gravitational torque () are


Then, the free mode of oscillation corresponding to the system (3) without Saturn’s gravitational torque () is expressed as


Finally, the proper frequencies of the system (3) forced by Saturn and internal gravitational torques can be written as


with . The values of and are given for each selected model in Table 3.

Using these expressions of the proper frequencies, the solutions of the system (3) are written as by setting , where and are respectively the amplitude of libration and the magnitude of perturbation, is the frequency of perturbation and the phase. The icy shell libration amplitudes for an oceanic model are then given by


The libration amplitudes for a rigid solid model are given in Appendix.

As seen in Sect. 2, the atmospheric torque exerted on Titan’s surface is acting like a forcing term with a main component at Saturnian semi-annual frequency and many components with amplitudes below kg m s. The atmospheric torque is introduced into the right-hand side of the first equation of the system (3) as , where is the atmospheric torque magnitude, the frequency and the phase difference with . The libration amplitude is then decomposed into a sine and a cosine term expressed as


where and are the amplitudes of the sine and cosine terms respectively of the libration under orbital and atmospheric influence, as previously is the perturbing frequency, its phase, and the phase difference with the orbital perturbation. The cosine term implies a difference of phase with the libration amplitude due to orbital perturbation. For frequencies of the atmospheric torque which are not present in the orbital motion, the librational motion is obtained by doing in Eq. (10).

4.2 Elastic case

The elasticity is introduced by modeling the radial deformations of the surfaces of the satellite layers and the resulting variations of the moments of inertia are computed. Each surface composing the satellite is parametrized by where is the mean radius of the equivalent sphere and the angles , are the colatitude and longitude, respectively. Then is the radial deformation induced by the centrifugal and tidal potential.

Hinderer et al. (1982) used a decomposition of into spherical harmonics. At the second order, we have


where and quantify the radial deformation normalized to and are associated Legendre polynomials. By using the definition of , the non-zero components of the inertia tensor can be written at first order in terms of , \linenomath


Love (1911) defined at the surface the Love number which characterizes the body response to second order perturbing potential. Here, we use the following definition for the radial deformation


where represents the surface gravity of Titan and is the perturbing potential developed at order . is the dimensionless radial function corresponding to at the surface. By using eqs. (12) and (17), a straight relation exists between the , and the response to the perturbing potential proportional to .

As shown by Giampieri (2004), the development at first order in eccentricity of the perturbing potential (centrifugal and tidal) can be decomposed into a secular and a periodic part \linenomath


where and are the mass and mean radius of the satellite, and with the mass of the planet and the semi-major axis of the satellite. The centrifugal potential has been supposed constant with a periodic variation of the potential only due to tides.

The secular part of the potential corresponds to static tides and is responsible for the static bulge of the satellite which is already included in our rigid case (see Sect. 4.1). The periodic part, which is of order of the eccentricity , governs time-dependent deformation of the surfaces.

The deformations of the satellite layers due to the tides induce a time-dependent gravitational potential. Mass distribution varies periodically and the inertia tensor of the body can be decomposed into a static part (corresponding to the static potential) and a periodic part associated to deformations around the static bulge. Using the definitions of and of the potential , the periodic components of the inertia tensor are given by


where the following definitions are used:




The total deformation factors are given by a combination of static deformations () and periodic deformations ().

The periodic components of the inertia tensor and are associated with and variations called radial tides and is associated with variations called zonal tides. The component corresponds to a periodic bulge shifted by in longitude from the principal axis of inertia (called librational tide)(Murray and Dermott, 1999; Greenberg et al., 1998).

The external gravitational torque exerted on a layer is given by where is the external gravitational potential of the deformed satellite layer given by (Jeffreys, 1976) \linenomath


where and are respectively the top and bottom radius of the layer . The z-component of this torque can be written at second order in eccentricity \linenomath


Here describes the torque amplitude exerted by the planet on the static bulge of the layer as in the rigid case (see Figure 3(a)). is the torque amplitude exerted on the radial tidal bulge of the layer with and is the contribution of the ocean pressure on the periodic figure based on equation (24).

The last two terms of equation (27) are the torque exerted by the planet on the librational bulge. It can be notice that the development in second order of for a Keplerian orbit cancels out , and in this case the radial tides terms vanishes. The signs of the remaining librational and static terms at first order in are in phase and opposed, which means that the torques are counter-acting as illustrated by the gravitational forces on Fig. 3(b).

The total internal gravitational torque exerted by outer layers on the interior is given by where is the potential exerted by outer layer on an inner point given by (Jeffreys, 1976) \linenomath


where refers to the inner core of mean radius . The z-component of the internal gravitational torque exerted on the inner core is then \linenomath


where and are defined such as \linenomath


with the deformation factor of the static bulge. In the internal gravitational torque, three main components are identified: the interaction between the static bulges of the shell and the inner core through like in the rigid case (see Figure 3(a)), the interaction between the periodic bulge of the shell and the static bulge of the core (subscript , Figure 3(c)), and the interaction between the static bulge of the shell and the periodic bulge of the core (subscript , Figure 3(d)). Here, the attraction between periodic bulges, which is of higher order in eccentricity, is neglected. Terms in are due to radial tides while terms in are due to librational tides. All the term but one have the same sign in the internal torque, which means that, for small and positive values of , they are acting in the same direction and tend to align the shell and the inner core figures. The last term is counter-acting however, it corresponds to the attraction of the static bulge of the core by the librational bulge of the shell (Fig. 3(c)).

Figure 3: Schematic position of the different forces acting on the bulges of the icy shell and the solid inner core when the satellite’s mean anomaly is small and positive. The angle with the planet position and the radial and librational bulges amplitude have been exaggerated. Red arrows are the gravitational forces exerted by the planet and blue arrows are the internal gravitational forces between solid layers. (a) The static bulges are attracted by the planet and attracting each other due to misalignment. (b) The planet exerts gravitational forces on the radial and librational bulges. The resulting torques tend to oppose themselves in the radial and librational components of a layer. Internal gravitational forces are acting between periodic bulges and the static bulge of the core (c) and of the shell (d). The gravitational interaction between the periodic bulges is not represented.

The equations of librations are then linearized and written at second order in eccentricity:


where the second term of the left-hand side is due to zonal tides on the layer with . The internal and external gravitational torques at first order in eccentricity are dominated by the static and librational tides and the second order is dominated by the radial tides.

4.3 Analytical resolution of the elastic case

The system (32) is composed of two differential equations of second order in and . To solve this system, we transform it into a first order equation system in by setting , , , and . The system is then decomposed into a static and a periodic part such as where z is the vector of components , is the vector of forcing terms, and are the matrices of the static and periodic coefficients, respectively. The static coefficient matrix is diagonalized by defining a vector such that where corresponds to the eigenvector matrix of . The system is then reduced to


where is the diagonal matrix, is the transformation of and and are the transformation of decomposed into first order and higher order terms, respectively.

To solve this new system, the perturbation method used by e.g. Robutel et al. (2011) is followed and is decomposed in decreasing amplitude terms such that where subscripts and refer to first and second order solutions. We then have to solve \linenomath


The system is first solved and the solutions are substituted in the system to obtain the second order solutions. Then the solutions are transformed back to our libration angles also decomposed in and . The solution is not computed here since we are only interested in solutions of order 2 in eccentricity.

Since the deformation is also periodic with frequency , the solutions of the elastic case are given by


where the third subscript corresponds to the frequency. At the orbital frequency, the amplitude is and \linenomath


As suggested in the previous section, the torques variation and present in the solution at first order in eccentricity are only due to librational tides components.

The solution for an elastic solid body without ocean is given in Appendix.

5 Librational response and application to interior models

5.1 Discussion of the rigid solution

The libration amplitude obtained in the previous sections for the rigid and elastic cases are described and analyzed in this section. The libration responds to the perturbations according to three regimes depending on the values of the forcing frequencies with respect to the values of the proper frequencies: a high forcing frequency regime, a low forcing frequency regime, and a resonant regime where the forcing frequencies are close to the proper frequencies.

For all regimes, the solutions are simplified with the assumption that the icy shell thickness is small. is about 10 to 20 times larger than for the six models selected in Sect. 1. So the proper frequencies can be developed at first order in and such as


Using these developments, it can be noticed that for low forcing frequencies , the libration amplitude of the rigid case is simply , i.e. the librational response of the body or the icy shell will follow the magnitude of the perturbation. In that case, the librational response is dominated by the external gravitational torque and the internal structure has a negligible influence.

If the atmospheric torque excites the system then the amplitude of the librational response to low frequency forcing will be different from the magnitudes of the orbital perturbations. This difference would contain information on the atmospheric torque coupled to internal structure through the sine term


At low forcing frequency, the librational response due to the atmosphere is then inversely proportional to the external gravitational torque amplitude on the shell and dependent on the internal gravitational coupling through the constant which depends on the inner core composition.

At high forcing frequencies, the rigid solution (eq. (9)) can be simplified to


when is large in front of and . Here, the dynamics is dominated by if . By writing as function of the equatorial flattening of the layers (see Van Hoolst et al. (2009)), and by noticing that the flattening of the ocean is almost equal to the one of the shell from our interior models (the variations are only ), is written


where the icy shell thickness and equatorial flattening have been introduced. The free frequency and the libration amplitude are then dependent on the ratio between ocean and icy shell densities, and inversely proportional to the icy shell thickness.

Finally, if the forcing frequencies are close to the proper frequencies of the system, the librational response is dominated by the resonant behavior. The amplitude of the libration is then significantly increased due to the presence of a small divisor. If the libration amplitude becomes too large, the linearized equations used are not enough to describe the dynamics. For Titan, the proper frequencies listed in Table 3 compared to the forcing frequencies listed in Table 2 show that there is no close resonance between these frequencies.

5.2 Discussion of the elastic solution

Here, we assess the effect of the elasticity on the libration solution which is dominant only for high forcing frequencies. The elasticity induces periodic variations of the torques amplitudes that are identified in the diurnal libration amplitude through two terms denoted by and (eq. 36) corresponding respectively to the variation of the shell and the inner core inertia. The form of the elastic case diurnal solution at first order in eccentricity (eq. (36)) is similar to the rigid case solution (eq. (9)) and it can also be simplified as


with large in front of and , and is defined as

Thus, the libration amplitude at orbital frequency will be largely reduced compared to the rigid case if is large enough in front of , i.e. if the shell librational tide contribution is large enough in front of the external gravitational torque exerted on the static bulge. To evaluate this contribution, the reduction rates of the torques amplitudes are defined as for the shell and for the inner core. Their numerical values are given in the next section.

The second order librational reponse to the forcing frequency is identical to the rigid case and the behavior is described in the previous section. Terms of frequencies also appear in the second order solution due to the modulation of the gravitational torque amplitude. Since can be small compared to , the low frequency forcing can contribute to the orbital frequency libration but these terms have small amplitudes below one meter. For the special case where is close to , the terms of frequency contribute to the libration with their amplitude being inversely proportional to the proper frequencies (see eq. 45). The solution can be written for the frequency , and


In that case, since the values of proper frequencies and are small (see Table 3), large values of magnitude can provide a non-negligible signature in the libration as detailed in the next section.

(rad day) (rad day) (rad day)
F1 -
F2 -
F3 - -
CA10 -
MC11 -
FE10 -
Table 3: Proper frequencies for the different interior models selected in this study and presented in Table 1. is defined in the Appendix.

5.3 Application

Now, the librational solutions are applied to the different interior models introduced in Sect. 1. The librations amplitudes projected onto Titan’s equator are given in Table 4. As discussed previously the librations can be classified into three regimes depending on the value of the forcing frequency with respect to the values of the proper frequencies.

For the low forcing frequencies, which are lower than the proper frequencies, the libration amplitude is almost equal to the magnitude of the orbital perturbation for any interior models. As seen in Table 4, the resulting libration amplitudes are large, around meters and meters for the Saturnian semi-annual and annual frequencies. The small differences in the amplitudes at Saturnian semi-annual and annual frequencies are due to a second order effect in the proper frequencies (eqs. (39) and (40)) and are related to the different inertia of the models, with a larger value of libration for the solid model induced by the small value of proper frequency .

Amplitudes (m)
Freq. Period F1 F2 F3 Identification
(rad/days) (days)
0.394018 15.946441 -62.061 -86.257 -50.661
0.001169 5376.6317 552.772 552.309 560.119
0.000584 10750.4115 470.287 470.183 471.838
0.001753 3584.9304 72.120 72.000 74.306
0.000063 99027.4111 47.477 50.171 31.500
0.001121 5606.2511 26.687 28.182 17.978
0.009810 640.4878 11.868 16.118 -24.275 -
0.001169 5376.6317 517.761 515.967 523.095 Tokano 05
0.001169 5376.6317 548.164 548.107 556.305 Charnay 12
0.394018 15.9464 -319.036 -384.527 -52.032 Rigid case
Amplitudes (m)
Freq. Period CA10 MC11 FE10 Identification
(rad/days) (days)
0.394018 15.946441 -72.751 -76.679 -75.420
0.001169 5376.6317 552.265 552.284 552.478
0.000584 10750.4115 470.171 470.174 470.217
0.001753 3584.9304 71.994 72.002 72.054
0.000063 99027.4111 51.963 51.996 50.714
0.001121 5606.2511 29.185 29.205 28.492
0.009810 640.4878 16.956 17.351 16.721 -
0.001169 5376.6317 515.887 515.995 516.814 Tokano 05
0.001169 5376.6317 548.084 548.055 548.028 Charnay 12
0.394018 15.9464 -391.131 -387.614 -358.397 Rigid case
Table 4: Analytical icy shell libration amplitudes of an elastic Titan (as equatorial deviation in meters) for different internal structure models with a kilometers icy shell thickness and different forcing frequencies. The second part of the table is the resulting librations under gravitational and atmospheric torques, with the given amplitude TO05 of Tokano and Neubauer (2005) and a recent value CH12 derived from Charnay and Lebonnois (2012) at Saturnian semi-annual frequency. The last part of the table is the libration amplitudes for models where solid layers are rigid. The libration amplitudes are truncated at meters.

In the case of rigid layers, the librational response at orbital frequency are expected to be strongly dependent on the satellite inertia (Van Hoolst et al., 2009). Differences between oceanic and solid models would then appear clearly through the variation of the amplitude by at least a factor six, as for example the values of 391 meters for CA10 and 52 meters for F3. In the rigid case, equation (43) shows that the amplitude would mainly depend on the ocean and icy shell densities ratio and on the icy shell thickness.

As described in 5.1, the total torque acting on the solid layers is modified by reduction rates and when elasticity is taken into account. To evaluate these rates, the radial function is computed by the numerical integration of the gravito-elastic equations set of Alterman et al. (1959) with the SatStress numerical code developed by Wahr et al. (2009). For the different interior models selected, the corresponding surface Love numbers are contained between and for oceanic models, and about for the solid model. In an elastic and solid Titan, only a small deformation is necessary for the shear stresses to balance the tidal force and thus a smaller value of . There are no shear stresses in an inviscid fluid layer allowing the Love number value to increase (Rappaport et al., 2008). For the oceanic models, values obtained in the shell are at least forty times larger than at the surface of the core. It means that the periodic radial deformation of the core is negligible compared to the shell as well as its radial and librational bulges. Thus the internal torque amplitude can be ignored in the forcing terms.

For the CA10 model, we find