Forced libration of tidally synchronized planets and moons

Tidal dissipation of kinetic energy, when it is strong enough, tends to synchronize the rotation of planets and moons with the mean orbital motion, or drive it into long-term stable spin-orbit resonances. As the orbital motion undergoes periodic acceleration due to a finite orbital eccentricity, the spin rate oscillates around the equilibrium mean value too, giving rise to the forced, or eccentricity-driven, librations. Both the shape and amplitude of forced librations of synchronous viscoelastic planets and moons are defined by a combination of two different types of perturbative torque, the tidal torque and the triaxial torque. Consequently, forced librations can be tidally dominated (e.g., Io and possibly Titan) or deformation-dominated (e.g., the Moon) depending on a set of orbital, rheological, and other physical parameters. With small eccentricities, for the former kind, the largest term in the libration angle can be minus cosine of the mean anomaly, whereas for the latter kind, it is minus sine of the mean anomaly. The shape and the amplitude of tidal forced librations determine the rate of orbital evolution of synchronous planets and moons, i.e., the rate of dissipative damping of semimajor axis and eccentricity. The known super-Earth exoplanets can exhibit both kinds of libration, or a mixture thereof, depending on, for example, the effective Maxwell time of their rigid mantles. Our approach can be extended to estimate the amplitudes of other libration harmonics, as well as the forced libration in non-synchronous spin-orbit resonances.


INTRODUCTION
Most of the moons and satellites in the Solar system have a small orbital eccentricity and are locked in the 1:1 spinorbit resonance, when the period of sidereal rotation is exactly equal to the orbital period. The commonly accepted explanation for both phenomena is sufficiently strong tidal interactions with their host planets, which whittle down the angular momentum and kinetic energy over the course of billions of years and drive the system into the global potential well. Mercury represents a notable exception with its relatively high eccentricity and a stable 3:2 spin-orbit resonance. Tidal dissipation of orbital energy continues in resonantly rotating bodies as long as the eccentricity is nonzero. An analytical solution for energy dissipation in a closed form was obtained for the simplest model case of a uniform (unstructured) spherical body (Peale & Cassen 1978;, and references therein). This analysis implicitly assumes that the spin rate of the body is nearly constant over one orbit, leading to a zero contribution ⋆ E-mail: vvm@usno.navy.mil from the leading (resonant) tidal mode. In reality though, the spin rate of a resonant body is subject to periodic oscillations, called librations in longitude. In most cases, librations are caused by the strictly periodic torque acting on the perturbed body (i.e., the body with tides) resulting from the misalignment of the "long" axis of the smallest moment of inertia and the instantaneous direction to the center of the perturbing body. The misalignment comes from the periodic acceleration of the perturber on the eccentric orbit. Bodies with complete rotational symmetry about the principal axis of inertia, such as completely liquid planets and stars, do not have a permanent figure and therefore, are not subject to this type of libration. Since the main type of libration is related to a certain force (in this case, the force of gravitation from the perturber acting on the triaxial figure), it is called forced libration. Wisdom (2004) suggested that forced librations of resonantly rotating bodies boost the dissipation of orbital energy and heating of the perturbed body. This possibility was also discussed by  in the framework of a general viscoelastic tidal model. The boost may happen because the secular component of the tidal torque depends c 2015 RAS on the instantaneous spin rate, rather than on the instantaneous libration angle, and may result in a negative work, being in opposite-phase with the velocity. Physically, the tidal bulge is in a small-angle periodical motion across the surface of the body with a certain phase lag, giving rise to additional internal friction. This surmise lacks a detailed, quantitative analysis. Qualitatively, this extra dissipation is expected to be small for nearly circular orbits and small libration amplitudes (such as the lunar forced libration), and the traditional approach to tidal heating estimation (Peale & Cassen 1978) is perhaps justified.
Numerous exoplanets have been found close enough to their host stars for the tidal forces to be important in the dynamical heating and orbital evolution. Close-in super-Earths, in particular, often seem to have unrelaxed eccentricities. They should be locked in the 1:1 or higher order spin-orbit resonances. The vividly discussed issue of possible habitability of such planets can not be decoupled from spin-orbit dynamics in general and libration in particular. Dobrovolskis (2007) considered the important role of optical librations for climate and habitability of exoplanets, but ignored the physical forced librations. Optical libration is a merely geometrical effect related to the changing aspect angle of the satellite surface relative to the line through the centers even in the absence of spin rate oscillations. Even on synchronized planets, a small amount of physical rocking motion may generate rather extensive zones of temperate irradiation suitable for biological life. For now, the parameters of exoplanet libration remain a matter of theoretical speculation, as there are no observational means to detect or measure it.
Small satellites of strongly asymmetric or deformed shapes represent another interesting area for libration theory applications, this time with a distinct possibility of verification by observational evidence. Thomas et al. (2016) measured a large libration longitude angle of 0.12 degrees for Enceladus, which they explained by the presence of a massive subsurface ocean. A large, nearly decoupled molten core in Mercury is betrayed by the observed amplitude of the forced libration (Margot et al. 2007), which is roughly twice as large as what is expected from a uniformly rigid planet. A non-vanishing friction between the solid mantle and the liquid core, or between the crust and a subsurface ocean in case of smaller satellites, can change not only the magnitude but also the spectrum of observed short-period librations. Additional complications arise for small satellites with a large degree of elongation when the natural oscillation frequency becomes commensurate with the frequency of forced excitation. Both the capture of such bodies into the spin-orbit resonance (Aleshkina 2009) and the post-capture rotation (Correia & Robutel 2013) may be chaotic. Goldreich & Mitchell (2010) investigated a more complex case of attenuated forced librations in a liquid body covered by a crust. The dynamic elongation of the liquid component (e.g., a subsurface ocean) interacts with the elastically deformable crust, the latter resisting the motion of the tidal bulge across the surface, resulting in additional dissipation of kinetic energy and reduction in the amplitude of libration. This model was further worked out and applied to realistic icy satellite scenarios by Van Hoolst et al. There the authors pointed out that tidal interactions and possible resonances between them can change the spectrum of observed librations. This elaborate model has been compared with lunar ranging observations.
In this paper, we consider the simplest case of a settled binary system with tides in the 1:1 (synchronous) resonance with a small eccentricity. We consider a uniform body neglecting subtle effects that may be caused by a multi-layered structure. The novelty of our approach is in an accurate analysis of the interplay between the gravitational and tidal perturbative torques, the latter having been often ignored in the literature on libration. We will show that even within this simple set-up, a new type of forced libration emerges, drastically different from the relatively well studied model of Lunar and Mercury's libration.

EQUATIONS OF PERIODIC LIBRATION
Ultimately, the perturbing torque generating forced librations is related to the orbital acceleration of the perturber, which is strictly periodic with the orbital period as long as long-term orbital evolution or precession of the equator is ignored. The terms of forced libration in this simplified setup are harmonics of the main orbital frequency. We can therefore represent the libration angle in the general form as a Fourier series of the mean anomaly M: with the coefficients αs and βp to be determined from the equation of forces. The constant term p = 0 is omitted in this expansion because its value is dependent on the choice of the origin for the angular coordinate θ. It is practical to reckon θ from the direction to the pericenter of the perturber as seen from the perturbed body to the instantaneous direction of the longer principal axis associated with the smallest moment of inertia A (let us call it "A-axis"). Still, the action of tidal torque on a triaxial body tends to accelerate the prograde rotation of a synchronized body (in other words, the time-average tidal torque is positive). The equilibrium is achieved through a negative time-average triaxial torque in the counter direction, which requires the A-axis to be leading the mean direction of the perturber (Makarov 2013). For the Moon, this constant tilt amounts to 67.753 arcsec (Rambaux & Williams 2011). Being in itself an important observable parameter, this tilt does not contribute to the energy dissipation or orbital evolution. Differentiating Eq. 1 in time and usingṀ = n obtainṡ The equation of forces is with TTRI and TTIDE being the triaxial and tidal torques, respectively, M2 the mass of the perturbed body, R its ra- dius, and ξ the moment of inertia coefficient usually ranging between 0.3 and 0.4.

COMPONENTS OF TIDAL TORQUE
We are using equations for the polar tidal torque (i.e., the component directed along the axis of rotation) from (Makarov 2012). Those expressions are based on the theory of tidal torque developed by Efroimsky (2012a,b) for a homogeneous near-spherical body of an essentially arbitrary rheology. As in ibid., we chose a realistic rheological law combining viscoelastic (Maxwell) response and inelastic (Andrade) creep caused by defect unpinning. Both components of the model are strongly frequency-dependent in the vicinity of spin-orbit resonances, and they can not be decoupled in a linear approximation. The inherent nonlinearity of the functional dependence on frequency in this model leaves little room for analytical applications, compared with the commonly used Constant Time Lag (CTL) model. However, we have to exploit more complex and accurate models in specific applications rather than the CTL model, which can not correctly explain the spectrum of lunar libration, for example, or the rate of eccentricity change in the Earth-Moon system. Among more advanced tidal models, Ferraz-Mello (2013) proposed a theoretical framework emphasizing the role of the Andrade dislocation and creep. Here, on the contrary, we will exploit a purely Maxwell approximation, i.e., a uniform viscoelastic body with rigidity and self-gravitation. Such approximations may shed light on the tidal behavior of planets and satellites of relatively low effective viscosity, for example, with post-solidus mantles including partial melt.
As was first pointed out by Efroimsky (2012a), the tidal torque includes both secular and oscillating components of time (or mean anomaly M): where ξ is the coefficient of inertia, C/(M2R 2 ), G lpq are Kaula's functions of eccentricity related to Hansen's coefficients through G 20(j−2) = X −3, 2 j , the tidal mode ω220q = (2 + q)n − 2θ, and the quality functions Kc and Ks are described in (Makarov 2015). The secular terms can be seen emerging from the overlap of tidal modes, i.e., q = j. In Section 6, we shall explain why the oscillating terms q = j can not be always ignored in the consideration of forced libration.

TRIAXIALITY-DRIVEN LIBRATION
When the mass of the perturbed body is much smaller than the mass of the (primary) perturber, M2 ≪ M1, the instantaneous acceleration of the deformed, triaxial body can be written as (e.g., Noyelles et al. 2014) Note that the summation in Eq. 5 is over j = −∞, . . . , +∞.
The three principal moments of inertia A, B, and C are unequal, with C > B > A by definition. The angle of rotation θ in this equation is not a linear function of time (or M), but varies periodically for a synchronous body in accordance with Eq. 1. Rigorous substitution of the Fourier series into Eq. 5 results in intractable infinite trigonometric series with the Bessel function coefficients. Fortunately, for the most important case of small librations, when θ − M is small, we can write, approximately, where Here we also limit our analysis to the case of slightly deformed bodies with (B − A)/C ≪ 1. Also in the following, we will only consider the main term proportional to sin M. This is justified as long as the orbital eccentricity is small, because the amplitude of j = 2 terms is approximately proportional to G 20 ( For a derivation of the complete Fourier series of triaxialitydriven libration see (Comstock & Bills 2003).
In the absence of other torques, an estimate for the main component of forced libration is readily obtained from Eq. 8, This formula was derived by Eckhard (1993). This simple calculation provides a coefficient of −15.5 arcsec with our best knowledge physical parameters for the Moon listed in Table 1. We find close agreement with the observed values from (Rambaux & Williams 2011), as well as with our own more accurate computation including the tidal torque, see Table 2. The tides on the Moon raised by the Earth result in the emergence of cosine harmonics of libration, most notably, of a positive (leading) tilt constant in time, but they little change the amplitudes of the sine harmonics. In this respect, the influence of the tidal torque on the Moon is small. Why?
The secular tidal torque in the viscoelastic Maxwell model with self-gravitation considered in this paper is a strongly nonlinear function of tidal mode in the vicinity of the 1:1 spin-orbit resonance (Fig. 1). As explained in detail in (Makarov & Berghea 2014), this function can be broken into two components, an antisymmetric kink shape and a nearly constant bias term, defined by the contribution of non-resonant tidal modes. The bias, or the offset, is barely noticeable in Fig. 1, but it does shift the function up by 2.2 · 10 −5 yr −2 at zero frequency. The secular prograde acceleration coming from the tides has to be counterbalanced by a secular triaxial torque acting in the opposite direction. The balance is automatically achieved and maintained, as the accelerating rotation increases the lead of the long axis relative to the average direction to the perturber until the perfect equilibrium of torques is reached. The bias thus generates the average tilt, i.e., the cos 0M harmonic of libration. The amplitude of this tilt is especially sensitive to the Maxwell time parameter of the model, which helps to constrain it to ≈ 11 yr for the Moon. The kink is characterized by a couple of opposed peaks and a roughly linear segment between them. Within this segment, and only there, this model becomes very similar to the much-explored Constant Time Lag model of tides.
Note that the main tidal mode and frequency are semidiurnal for the synchronous rotation, i.e., |ω2200| = ν2200 = |2(n −θ)|. The dimensionless parameter Λ2, often called effective rigidity and denominated A2 in the literature, emerges from the definition of the quadrupole static Love number (MacDonald 1964). The peak tidal quality at this wave number is The loci of the peak torque depend on both the Maxwell time and effective rigidity Λ2, whereas the magnitude depends only on the latter. The effective rigidity varies in a relatively narrow range between ∼ 1 for super-Earth exoplanets, which are possibly the most massive planets of terrestrial composition, and several tens for the smaller satellites of the Solar system (for the Moon, it is approximately equal to 65). The range of τM is many orders of magnitude, because the viscosity of silicate minerals is a steep function of temperature (e.g., Běhounková et al. 2011). Therefore, the spread of the tidal peaks near a resonance is quite sensitive to the effective viscosity of the dissipating layers. Cold, inviscid planets such as Earth and the Moon have Maxwell times of tens to hundreds of years. Semiliquid or icy satellites such as Io and Titan are likely to have this parameter in the order of days. Consequently, the peaks can be tightly packed at the point of resonance, as we see for the Moon in Fig. 1, or they can spread out as far as the adjacent resonances, or even farther. The forced libration remain triaxiality-driven when the peaks are well within the range of regular velocity oscillation. Combining Eqs. 8 and 10, the resulting condition for synchronous rotation is When this condition is fulfilled, the sin M variation of velocity takes the tidal torque well over the peak values to the domain where it becomes rather insignificant. This happens because the perturbed body, rocked by the periodically accelerating perturber, passes the area of tidal peaks very quickly at M ≈ 0 and π (i.e., at the pericenter and apocenter), spending more time in quadratures, where the spin rate deviates the most from the resonance and the tidal torque is weak. The critical value of τM for the Moon beyond which the forced libration is triaxiality-driven, is 1.2 yr.

LIBRATION OF HIGHLY ELONGATED SYNCHRONOUS BODIES
Equation 8 is accurate enough for rotating nearly symmetric bodies with small eccentricities, but smaller satellites often have large triaxiality parameters (B − A)/C. A more accurate treatment may be required, taking into account the possible interferences with the natural mode of free libration. The classical equation of the spin-orbit problem with a trixial torque reads (Wisdom et al. 1984;Murray & Dermott 1999) To compute the evolution of the libration angle in the synchronous resonance, we make the change of variables γ = θ − pM , choose p = 1 to select the synchronous resonance, and geẗ where ω 2 0 = 3 B−A C Gmp a 3 . If we choose to neglect γ on the r.h.s. of the above equation (because of its smallness close to the exact resonance), we can expand a 3 r 3 sin 2(f − M ) in eccentricity and obtain the Fourier decomposition of γ(t) in multiples of nt (see the complete Fourier decomposition in Comstock & Bills (2003)). On the other hand, we can obtain a more complete formulation in case ω0 is not negligible compared to n, by expanding Eq.14 in Taylor series of γ and keeping only the terms of zeroth and first order:    Expanding this equation in terms of the eccentricity functions F + q , F − q of (Comstock & Bills 2003) 1 , we havë where we kept only the zeroth order term in the expansion of a 3 r 3 cos 2(f −M ) = 1+3e cos M +e 2 (− 5 2 + 17 2 cos 2M )+O(e 3 ). This equation can be solved to give the Fourier decomposition of the libration angle γ(t) = A0 cos(ω0t + φ0) where the C −3,2 q , S −3,2 q are the Cayley coefficients (see Comstock & Bills (2003)) and A0, φ0 depend on initial conditions. We can merge terms with similar frequencies to obtain The first term of this equation represents the regular free libration at the frequency ω0 with arbitrary initial parameters. Free libration is rather quickly damped at resonances and is not considered in this paper. The second term is a refinement of the main mode of forced libration (8), with a correction coefficient n 2 /(n 2 − ω 2 0 ). This correction implies that the amplitude of forced libration becomes infinitely large when n = ω0, or (B − A)/C ≈ 1/3 for the perturber's mass mp much greater than the mass of the perturbed body. This only means that the above equation does not work at the points of singularity ω0 = q n for positive integer q (Wisdom et al. 1984). Some satellites have triaxiality parameters close to 1/3, for example, Epimetheus, with a (B − A)/C ≈ 0.30 (Tiscareno et al. 2009) or 0.28 (Caudal 2013). Should we expect a very large sin M forced libration for Epimetheus in accordance with Eq. 18? More detailed dynamical analysis indicates that instead of singularity and a sudden change of libration phase, another island of equilibrium emerges in the phase space (Melnikov & Shevchenko 2010) with the average rotational velocity still locked in the 1:1 resonance.

TIDALLY DRIVEN LIBRATION
We have considered the case of relatively cold and highviscosity bodies with Maxwell times significantly larger than the critical value defined in Eq. 12. This model can be used for the Moon, Earth and Mercury. At the other end of the gamut, semiliquid, semi-molten planets or icy satellites are 1 Mind the small misprint in Table 2 of Comstock & Bills (2003) where F − 1 = G 2,0,−1 = − e 2 + e 3 16 − 5e 5 384 characterized by low viscosity but finite sheer modulus, resulting in a short Maxwell time well below the critical value.
We will see in the following that this circumstance may lead to a very different spectrum of forced libration. So, if the "secular" component of tidal acceleration can be well approximated as where see (Makarov 2015). This approximate equation just reflects the fact that the segment of tidal torque between the peaks resembles a straight line with negative slope (Fig. 1). Since the main term of libration in angle emerging from the triaxial torque is proportional to − sin M, the main term in spin rate is proportional to − cos M, and the acceleration θs, by Eq. 20, is proportional to + cos M. Therefore, this component of tidal torque gives rise to libration angle terms proportional to − cos M. Being in phase with the triaxialitydriven spin rate oscillation, it defines the rate of energy dissipation and the orbital evolution of synchronized bodies. The oscillating part of the tidal torque includes both sin jM and cos jM harmonics (Efroimsky 2012b;Makarov 2015). When the eccentricity of a synchronously rotating body is small, e < 0.1, and the τM is longer than the orbital period, the former lot can be neglected, because the kvalitet function Kc(2, ω20q) is usually very small at nonresonant tidal modes (q = 0). This can not be said about the sine-harmonics of the oscillating torque. The corresponding acceleration is, specifically, The kvalitet function is, specifically, where we dropped the indices 20q for all frequencies. At zero (resonant) tidal frequency, Ks(0) = 3/2, which is at least twice as large than the peak quality for the secular part and sine-harmonics (Eq. 11). It follows that the sine terms can be quite significant.
Introduce ∆K = Ks(2, 0) − Ks(2, n) = 3 2 where the tidal wave number λ = τM n, leading to an approximate solution for the sin M and cos M terms Two important conclusions for tidally synchronized bodies emerge right away: 1) the cos M and sin M terms in libration angle are both negative, and their ratio is defined by the "tidal feedback" coefficient Y ; 2) The classic, triaxiality driven libration term −4Ze sin M is severely suppressed if Y ≫ 1, where the cos M becomes dominant.

SUMMARY AND DISCUSSION
The three dimensionless parameters, X =X∆K, Y , and Z, as well as the orbital eccentricity e, define the spectrum of forced libration of tidally synchronized bodies. Even though our knowledge of the physical parameters listed in Table 1 is imperfect, especially for the exoplanet Kepler 10b, estimating the X, Y, Z parameters for the three objects allows us to outline the range of possible libration spectra. Table 3 gives our best-guess estimates for the Moon, Io, and Kepler 10b.
First, we notice that the "tidal feedback" coefficient Y can indeed be much larger than 1. For the Moon, this coefficient amounts to ∼ 5, which would imply a dominating − cos M term in angle libration per Eq. 30, and a much damped triaxiality driven libration. This, however, is not the case, as both observational data (Table 2) and numerical simulations (Fig. 2) indicate. The reason is explained in Section 4. The dominant term is obviously − sin M, and it is still a triaxiality driven spectrum, because the crucial condition 12 is met. The presence of cosine-terms is only betrayed by the asymmetry of the curve in Fig. 2, which is flatter at the top and sharper at the bottom, indicating that the Moon spends a longer part of the orbital period at a spin rate above the resonance. But for the less viscous Io, condition 19 is satisfied, and the libration curve in angle is almost exactly a cosine (Fig. 3). We conclude that even among the bodies of the Solar system, a large variety of libration spectra should be found, with relatively inviscid, icy satellites such as Io and possibly Titan exhibiting tidally driven, cos M-dominated, libration of much reduced amplitudes due to the tidal feedback.
The case of Kepler 10b stands out because of subtle circumstances. On the one hand, Y is much smaller than unity, so the main term in the spectrum of angle libration is − sin M, just like the Moon. One of the reasons for such small Y is a low effective rigidity Λ2 typical of massive super-Earth exoplanets. This is not a triaxiality driven spectrum, on the other hand, because X ≫ Z, thus the the main action comes from the oscillating component of the tidal torque. Although neither eccentricity nor (B − A)/C are known for this exoplanets, and the numbers given in Table 1 are rather placeholders, the permanent figure of this massive planet is expected to be small, and the eccentricity low too because of the high rate of tidal dissipation . Therefore, the character of libration is likely to be similar to that of the Moon, but for an entirely different reason.
The spectrum of libration of synchronous bodies is important for the tidal energy dissipation (and thus, heating) and the tidal orbital evolution.  surmised that forced libration can boost the energy dissipation in resonance allowing the resonant component of the secular torque to do some friction work. This may be valid for triaxiality driven libration, where the secular tidal torque is in counter-phase with instantaneous velocity. But we have determined that Io's libration must be tidally driven, with the sine-harmonic severely damped, so that the boost to energy dissipation is expected to be negligible.