© ESO, 2014

1. Introduction and context

Tides are one of the key interactions that drive the evolution of planetary systems. Indeed, because of the friction in the host-star and in the planet interiors, a system evolves either to a stable state of minimum energy, where spins are aligned, orbits circularised, and the rotation of each body is synchronised with the orbital motion, or the perturber tends to spiral into the parent body (Hut 1980). Therefore, understanding and modelling the dissipative mechanisms that convert the kinetic energy of tidally excited velocities and displacements into heat is of great importance. These processes, driven by the complex response of a given body (either a star or a planet) to the gravific perturbation by a close companion, strongly depends on its internal structure and its rheology. Indeed, the tidal dissipation in solid (rocky/icy) planetary layers strongly differs from the dissipation in fluid regions in planets and in stars; the one in rocks and ices is often strong with a smooth dependence on the tidal frequency χ, the one in gas and liquids being generally weaker on average and strongly resonant. Therefore, these properties must be taken into account in the study of the dynamical evolution of planetary systems using celestial mechanics.

To reach this objective, the tidal quality factor Q has been introduced in the literature (Goldreich & Soter 1966). Its definition comes from the evaluation of the tidal torque (Kaula 1964) and the analogy with forced damped oscillators: it evaluates the ratio between the maximum energy stored in the tidal distortion during an orbital period and the energy dissipated by the friction. Indeed, a low value of Q corresponds to a strong dissipation, and vice versa. In this framework, Q can be computed from an ab initio resolution of the dissipative dynamical equations for the tidally excited velocities and displacements in fluid and in solid layers of celestial bodies, respectively (e.g. Henning et al. 2009; Efroimsky 2012; Remus et al. 2012b; Zahn 1977; Ogilvie & Lin 2004, 2007; Remus et al. 2012a). This leads to values of Q that varies smoothly as a function of χ in rocks and ices, while numerous and strong resonances are obtained in fluids. However, in celestial mechanics studies, Q is often assumed to be constant or to scale as χ^-1 as a convenient first approach and is evaluated using the scenario for the formation and evolution of planetary systems.

In this work, we show how the dependence of Q on χ impacts this evolution and that it must be taken into account. In Sect. 2, we describe the set-up we studied and the corresponding dynamical equations, which correspond to those adopted by Efroimsky & Lainey (2007), who studied the impact of the rheology of solids on related tidal dissipation and evolution (Sect. 3). In Sect. 4, we study highly resonant tidal dissipation in fluid layers and discuss the strong differences with solids. Finally, we discuss astrophysical consequences for the evolution of planetary systems.

2. Set-up and dynamical equations

2.1. Model

To study the impact of rheology on tidal evolution and of the related variation of Q as a function of χ, we chose to follow Efroimsky & Lainey (2007). We therefore studied a two-body coplanar system with a central extended body A with a mass M_A and a mean radius R_A, and a point-mass tidal perturber B of mass M_B. In a reference-fixed frame ℛ_A:{A,X_A,Y_A,Z_A}, the central body rotates with a spin vector Ω_A and the perturber orbits it. To study a simplified system where we can easily isolate the effect of rheology, this motion was assumed to be circular. Thus, the position of B is directly given by the semi-major axis a, which is the distance that separates B from A in this particular case, and the mean anomaly $\begin{array}{c}{}_{\mathrm{˜}}\\ {\mathit{M}}_{\mathrm{B}}\end{array}\mathrm{=}{\mathit{n}}_{\mathrm{B}}\mathit{t}$, where n_B is the mean motion and t is the time coordinate.

2.2. Dynamical equations

As recalled in the introduction, the tidal quality factor Q is defined as the ratio between the maximum energy stored in the tidal distortion during an orbital period and the energy dissipated by the friction. It is thus directly related to the rheology of the studied bodies, which generally leads to a dependence of Q as a function of the tidal frequency $\mathit{\chi }\mathrm{=}{\mathit{\chi }}_{\mathrm{2200}}\mathrm{=}\mathrm{2}\hspace{0.17em}\mathrm{|}{\mathit{n}}_{\mathrm{B}}\mathrm{-}{\mathrm{\Omega }}_{\mathrm{A}}\mathrm{|}\mathrm{=}\mathrm{\left|}{\mathit{\omega }}_{\mathrm{2200}}\mathrm{\right|}\mathrm{=}\mathrm{\left|}\mathit{\omega }\mathrm{\right|}$(1)(e.g. Greenberg 2009; Efroimsky 2012), where ω = ω₂₂₀₀ is the principal semi-diurnal Fourier tidal mode that corresponds to the frequency of the perturbation in the frame that rotates with the perturber. This friction induces a geometrical angle δ(χ) between the directions of the tidal bulge and of the line of centres¹$\mathit{\delta }\left(\mathit{\chi }\right)\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}\mathit{\chi }\mathrm{\Delta }\mathit{t}\left(\mathit{\chi }\right)\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}{\sin }^{-1}\mathrm{[}{\mathit{Q}}^{-1}{\left(\mathit{\chi }\right)}^{\mathrm{]}}\mathit{,}$(2)where we have introduced the so-called time lag Δt(χ) (see e.g. Hut 1981).

This lag induces a net torque that modifies the evolution of the spin of body A (Mathis & Le Poncin-Lafitte 2009) $\frac{\mathrm{d}{\mathrm{\Omega }}_{\mathrm{A}}}{\mathrm{d}\mathit{t}}\mathrm{=}\frac{\mathrm{3}}{\mathrm{2}}\frac{{\mathit{k}}_{\mathrm{2}}\mathrm{\left(}\mathit{\chi }\mathrm{\right)}\mathit{G}{\mathit{M}}_{\mathrm{B}}^{\mathrm{2}}{\mathit{R}}_{\mathrm{A}}^{\mathrm{5}}}{{\mathit{I}}_{\mathrm{A}}{\mathit{a}}^{\mathrm{6}}}{\mathit{Q}}^{-1}\left(\mathit{\chi }\right)\mathrm{sgn}\mathrm{(}{\mathit{\omega }}^{\mathrm{)}}\mathit{,}$(3)where I_A is the moment of inertia of body A, k₂(χ) is the Love number, and G is the gravitational constant. The semi-major axis of body B is also modified (Efroimsky & Lainey 2007) $\frac{\mathrm{d}\mathit{a}}{\mathrm{d}\mathit{t}}\mathrm{=}\mathrm{-}\frac{\mathrm{3}{\mathit{k}}_{\mathrm{2}}\mathrm{\left(}\mathit{\chi }\mathrm{\right)}{\mathit{R}}_{\mathrm{A}}^{\mathrm{5}}{\mathit{n}}_{\mathrm{B}}{\mathit{M}}_{\mathrm{B}}}{{\mathit{M}}_{\mathrm{A}}{\mathit{a}}^{\mathrm{4}}}{\mathit{Q}}^{-1}\left(\mathit{\chi }\right)\mathrm{sgn}\mathrm{(}{\mathit{\omega }}^{\mathrm{)}}\mathit{.}$(4)These equations show that Q^-1 has an explicitly linear impact on the evolution of the system; strong variations of Q^-1 thus imply rapid changes for a and Ω_A. Then we took into account the dependence of Q to χ (e.g. Mathis & Le Poncin-Lafitte 2009; Efroimsky & Makarov 2013) to consider the impact of rheological models straightforwardly.

3. Solid tides

To solve our problem, we must close it with the choice of a law giving Q (or more generally, k₂/Q) as a function of χ. It is common to assume Q to be either constant (e.g. MacDonald 1964) or to scale as χ^-1 in the case of a constant tidal time lag (e.g. Hut 1981). However, for solid rocky or icy bodies, Efroimsky & Lainey (2007) suggested to use a power scaling law²: Q = ℰ^αχ^α, which has been experimentaly validated for metals and silicates in the laboratory as well as in seismic and geodetic experiments. Here, the empirical parameters α and ℰ are bound to the rheology; α characterises the frequency dependence and takes values between 0.1 and 0.4; ℰ is an integral relaxation parameter, which has the dimension of time and is determined by the internal-friction mechanism that dominates at the frequency χ. Usually, one or another mechanism or group of mechanisms remains dominant across vast bands of frequencies. Over these bands, ℰ may be regarded as constant or almost constant. This law is particularly interesting since it introduces the frequency dependence with only one more parameter than the constant law and remains close to the realistic physics of solids at the same time.

To compute the evolution of a and Ω_A with time, we chose to use the numerical code developed by one of us and based on the code ODEX (Hairer et al. 2000). To validate it, we studied the case of the Mars-Phobos system simulated by Efroimsky & Lainey (2007), assuming exactly the same parameters that are summarised in Table 1. The Phobos initial semi-major axis, the initial dissipative time lag, and the constant tidal quality factor are denoted a₀, Δt₀ and Q.

We note that in this case, the spin of A does not change over time compared to a because $\left|\frac{\mathrm{d}{\mathrm{\Omega }}_{\mathrm{A}}}{{\mathrm{\Omega }}_{\mathrm{A}}}\mathrm{·}\frac{\mathit{a}}{\mathit{da}}\right|\mathrm{=}\frac{\mathit{G}{\mathit{M}}_{\mathrm{A}}{\mathit{M}}_{\mathrm{B}}}{\mathrm{2}{\mathrm{\Omega }}_{\mathrm{A}}{\mathit{n}}_{\mathrm{B}}\mathit{a}{\mathit{I}}_{\mathrm{A}}}\mathrm{\approx }{\mathrm{10}}^{-7}\mathrm{\ll }\mathrm{1.}$(5)Our results perfectly reproduce those obtained by Efroimsky & Lainey (2007). Figure 1 shows the evolution of a with time for different values of α = −1 (constant tidal time lag Δt₀; e.g. Singer 1968; Mignard 1979), 0 (constant Q; e.g. Kaula 1964), 0.2, 0.3, and 0.4. These plots highlight the impact of the rheology on the smooth induced evolution of orbital parameters such as the semi-major axis and on the related life-time of the system.

Fig. 1 Temporal evolution of the semi-major axis a for a class of tidal models with k2/Q ∝ χ−α and with different values of the parameter α. The abscissa represents time in years, the vertical axis measures the deviation of the semi-major axis from its initial value a0. Our plots serve to compare realistic rheologies (α = 0.2,0.3,0.4) with less physical ones – those with α = 0 (Kaula 1964) and α = −1 (Singer 1968; Mignard 1979).

4. Tidal dissipation in fluid layers

4.1. Frequency dependence of Q

Revisiting the work by Efroimsky & Lainey (2007) provides a strong basis for exploring the impact of tidal dissipation in fluid bodies. In this section, we therefore chose to study the evolution of a perturber of the mass of Phobos that orbits a hypothetical completely fluid central body with the mass of Mars. Therefore, only the tidal quality factor Q(χ) will change.

In fluid bodies, tidal dissipation is caused by the turbulent viscous friction that acts on the equilibrium tide and on inertial waves, which are driven by the Coriolis acceleration, in convective regions (e.g. Zahn 1977; Ogilvie & Lin 2004, 2007; Remus et al. 2012a), and on thermal and viscous diffusions that act on gravito-inertial waves in stably stratified zones (e.g. Zahn 1977; Ogilvie & Lin 2004, 2007). The excitation of these oscillation eigenmodes by tides then leads to a highly resonant dissipation.

To illustrate our purpose, we consider from now that the central body is completely convective and rotates rapidly so that 0 ≤ σ ≤ 1, where σ = χ/(2Ω_A), which generates tidally excited inertial waves. There clearly is a strong difference between the tidal quality factor adopted before for solid bodies that scale as a smooth power-law of σ and the one related to inertial waves. Indeed, as demonstrated by Ogilvie & Lin (2004), using a local approach, their viscous dissipation is expressed as a sum of corresponding resonant terms³$\mathit{D}\mathrm{(}{\mathit{\sigma }}^{\mathrm{)}}\mathrm{=}{\mathit{D}}_{\mathrm{0}}\sum _{\left\{\mathit{m,n}\right\}\mathrm{\in }{N}^{\mathrm{\ast }}\mathrm{\times }{N}^{\mathrm{\ast }}}\frac{\mathrm{(}{\mathit{m}}^{\mathrm{2}}\mathrm{+}{\mathit{n}}^{\mathrm{2}}\mathrm{)}\left|\stackrel{}{{{\mathit{\sigma }}^{\mathrm{2}}}_{\mathrm{˜}}}\right|\mathrm{+}{\mathit{n}}^{\mathrm{2}}}{{\left|\mathrm{(}{\mathit{m}}^{\mathrm{2}}\mathrm{+}{\mathit{n}}^{\mathrm{2}}\mathrm{)}\hspace{0.17em}\stackrel{}{{{\mathit{\sigma }}^{\mathrm{2}}}_{\mathrm{˜}}}\mathrm{-}{\mathit{n}}^{\mathrm{2}}\right|}^{\mathrm{2}}}\mathrm{(}{\mathit{m}}^{\mathrm{2}}\mathrm{+}{\mathit{n}}^{\mathrm{2}}\mathrm{)}{\left|\mathit{n}{\mathit{f}}_{\mathit{mn}}\mathrm{-}\mathit{m}{\mathit{h}}_{\mathit{mn}}\right|}^{\mathrm{2}}\mathit{,}$(6)where $\begin{array}{c}{}_{\mathrm{˜}}\\ \mathit{\sigma }\end{array}\mathrm{=}\mathit{\sigma }\mathrm{+}\mathit{iE}{}^{\mathrm{(}}\mathit{m}^{\mathrm{2}}\mathrm{+}{\mathit{n}}^{\mathrm{2}}{}^{\mathrm{)}}$ and E = ν/⁽2Ω_AL²⁾ is the Ekman number of the fluid, ν is the viscosity, and L a characteristic length; m and n are the vertical and horizontal wave-vectors of inertial waves, respectively; finally, f_mn and h_mn are the coefficients of the Fourier series of the excitation. The tidal dissipation is thus a complex set of resonant peaks depending on the viscosity and on the rotation of the fluid. Since $\mathit{Q}{}^{\mathrm{\right(}}\mathit{\sigma }^{\mathrm{\left)}}\mathrm{\propto }{\left[\mathit{D}{}^{\mathrm{\right(}}\mathit{\sigma }^{\mathrm{\left)}}\right]}^{-1}$, we coupled it with the dynamical Eqs. (3), (4) of our model.

Fig. 2 Resonant tidal dissipation spectrum D resulting from inertial modes as a function of the normalised tidal frequency σ = χ/(2ΩA) assuming Nmax = 5.

Fig. 3 Evolution of the semi-major axis a over time with a Q factor proportional to inertial wave dissipation in fluids (green curve), and with a constant Q factor (blue dashed curve). The abscissa represents time in years, the vertical axis measures the evolution of the semi-major axis from the initial value a0.

Fig. 4 Evolution of the tidal dissipation and of the semi-major axis with time for different values of lp and Hp, the width at half-height and the height of the studied single resonant damping peak. The grey dotted line corresponding to α = 0.2 is superposed to the continuous green one except at the position of the peak for Q-1.

4.2. Numerical integration

To evaluate the effects of these resonances on dynamics, we computed the evolution of the semi-major axis of the orbit with the same parameters as for solid tides, but giving as input a synthetic ${\mathit{Q}}^{-1}{}^{\mathrm{\right(}}\mathit{\sigma }^{\mathrm{\left)}}$ factor written as $\mathit{D}{}^{\mathrm{\right(}}\mathit{\sigma }^{\mathrm{\left)}}$ given in Eq. (6). Our fluid is characterised by its Ekman number, E = 10^-5, which is a value often adopted in the literature for planetary convective layers, and which allows one to derive a peaked dissipation (see Fig. 2)⁴. The maximal rank of the sum (N_max) was chosen to be relatively low, with N_max = 5, to increase the computation speed. Following Ogilvie & Lin (2004), we describe the excitation with the coefficients $\begin{array}{c}\\ {\mathit{f}}_{\mathit{mn}}\mathrm{=}\frac{\mathrm{1}}{\mathit{m}{\mathit{n}}^{\mathrm{2}}}\mathit{,}& {\mathit{g}}_{\mathit{mn}}\mathrm{=}\mathrm{0}\mathit{,}& \mathrm{and}& {\mathit{h}}_{\mathit{mn}}\mathrm{=}\mathrm{0.}\end{array}$(7)The simulation clearly shows that, in contrast to solid tides, where the power scaling law implies a smooth evolution of a, a contrasted Q factor, that strongly depends on the tidal frequency, enables abrupt changes of a (see Fig. 3). As the perturber comes nearer to the central body, its mean motion increases. That is why dissipation varies strongly during the evolution of the system and, at each time it meets a resonance, there is a jump of a, which is bounded to the properties of the peak: the higher and wider the peak, the higher the amplitude of the jump. This is the resonance-locking identified by Witte & Savonije (1999) in the stellar case. However, we note that when σ > 1, at the end of the simulation, the evolution of a becomes smooth again. The reason for this behaviour is that we are beyond the range of frequencies where inertial waves are excited. Then, the tidal dissipation is that of the equilibrium tide that corresponds to the non-resonant background of D observed in Fig. 2. Finally, as demonstrated in Fig. 3, the evolution of a system where the tidal dissipation is caused by resonant eigenmodes (here the inertial waves) cannot be described properly using models where Q or Δt are assumed to be constant (α = 0 and α = −1, respectively).

4.3. Scaling law

The resonant properties of the tidal dissipation in fluids (see Fig. 2) generate jumps of the value of the semi-major axis a during the evolution of the system (see Fig. 3). In this framework, the link between the orbital dynamics and the rheology of the fluid is the shape of resonances of the viscous dissipation D (see Eq. (6)). A resonance occurs when a term of the sum in Eq. (6) exceeds all the others. In this section, our goal is therefore to obtain a scaling law that relates a jump of a to the height, H_p, and to the width at half-height, l_p, of the corresponding single resonant damping peak (see Fig. 4) defined as ${\mathit{Q}}_{\mathrm{p}}^{-1}\mathrm{(}{\mathit{\sigma }}^{\mathrm{)}}\mathrm{=}\frac{{\mathit{H}}_{\mathrm{p}}}{{\left[\mathrm{4}\mathrm{(}\sqrt{\mathrm{2}}\mathrm{-}{\mathrm{1}}^{\mathrm{)}}{\left(\frac{\mathit{\sigma }\mathrm{-}{\mathit{\sigma }}_{\mathrm{p}}}{{\mathit{l}}_{\mathrm{p}}}\right)}^{\mathrm{2}}\mathrm{+}\mathrm{1}\right]}^{\mathrm{2}}}\mathit{,}$(8)where σ_p is the resonant frequency. Then, the dissipation Q^-1 was chosen to be the sum of a smooth background denoted ${\mathit{Q}}_{\mathrm{0}}^{-1}$ that corresponds to the one studied in Sect. 3, and of a resonant one ${\mathit{Q}}_{\mathrm{p}}^{-1}$ (Eq. (8)), which leads to the following equation for a using Eq. (4): $\frac{\mathrm{d}\mathit{a}}{\mathrm{d}\mathit{t}}\mathrm{=}\mathrm{-}\frac{\mathrm{3}{\mathit{k}}_{\mathrm{2}}{\mathit{R}}_{\mathrm{A}}^{\mathrm{5}}{\mathit{n}}_{\mathrm{B}}{\mathit{M}}_{\mathrm{B}}}{{\mathit{M}}_{\mathrm{A}}{\mathit{a}}^{\mathrm{4}}}\mathrm{[}{\mathit{Q}}_{\mathrm{0}}^{-1}\mathrm{(}{\mathit{\sigma }}^{\mathrm{)}}\mathrm{+}{\mathit{Q}}_{\mathrm{p}}^{-1}{\mathrm{(}{\mathit{\sigma }}^{\mathrm{)}}}^{\mathrm{]}}\mathrm{sgn}\mathrm{\left(}\mathit{\omega }\mathrm{\right)}\mathit{.}$(9)Assuming that the peak influences the system when the condition ${\mathit{Q}}_{\mathrm{p}}^{-1}\mathrm{\ge }{\mathit{Q}}_{\mathrm{0}}^{-1}$ is fullfilled, and that the resulting variation is rapid compared with the mean evolution, we can derive the amplitude of the jump, $\frac{\mathrm{\Delta }\mathit{a}}{\mathit{a}}\mathrm{\approx }\frac{\mathrm{2}{\mathit{l}}_{\mathrm{p}}}{\mathrm{3}\sqrt{\sqrt{\mathrm{2}}\mathrm{-}\mathrm{1}}\mathrm{(}\mathrm{1}\mathrm{+}{\mathit{\sigma }}_{\mathrm{p}}\mathrm{)}}{\left[\sqrt{\frac{{\mathit{H}}_{\mathrm{p}}}{{\mathit{Q}}_{\mathrm{0}}^{-1}\mathrm{(}{\mathit{\sigma }}_{\mathrm{p}}\mathrm{)}}}\mathrm{-}\mathrm{1}\right]}^{\frac{\mathrm{1}}{\mathrm{2}}}\mathrm{·}$(10)In Fig. 4, we plot the evolution of the semi-major axis for different values of H_p and l_p and the corresponding dissipation. These graphs illustrate the scaling law (Eq. (10)) by showing that the width of a peak has a greater impact on a than its height. Moreover, the values of Δa/a obtained using direct numerical simulations perfectly match those predicted by Eq. (10). Finally, as σ_p, H_p, l_p are directly related to the value of the Eckman number E = ν/⁽2Ω_AL²⁾ (cf. Ogilvie & Lin 2004), we see how

the orbital dynamics is directly impacted by the fluid rheology and resonances.

5. Conclusions

We examined the impact of the frequency dependence of tidal dissipation in solids and fluids on the orbital evolution of a coplanar two-body system. We showed the strongly different evolutions induced by tides in rocks and by tides exerted on fluid layers where eigenmodes are resonantly excited. A smooth dependence of the tidal dissipation on the tidal frequency drives a smooth orbital evolution, while a peaked dissipation induces an erratic one. In each case, we pointed out the direct impact of the rheology properties on the dynamics of the system. Finally, we showed that it is important to take the dependence of the tidal dissipation on the tidal frequency into account and the important consequences it may have for the evolution of star-planet(s) and planet-moon(s) systems in the solar and exoplanetary systems. In this context, the impact of the frequency dependence of the tidal torque on resonances will be examined in a forthcoming work.

Acknowledgments

The authors are grateful to the referee, M. Efroimsky, for his detailed review, which has allowed us to improve the paper. This work was supported by the Programme National de Planétologie (CNRS/INSU), the GRAM specific action (CNRS/INSU-INP, CNES), the Paris Observatory, the Campus Spatial de l’Université Paris Diderot, and the Emergence-UPMC grant (contract number: EME0911). C.L.P.L. and S.M. dedicate this article to Dr. M. Le Poncin.

References

All Tables

All Figures

Fig. 1 Temporal evolution of the semi-major axis a for a class of tidal models with k2/Q ∝ χ−α and with different values of the parameter α. The abscissa represents time in years, the vertical axis measures the deviation of the semi-major axis from its initial value a0. Our plots serve to compare realistic rheologies (α = 0.2,0.3,0.4) with less physical ones – those with α = 0 (Kaula 1964) and α = −1 (Singer 1968; Mignard 1979).

In the text

	Fig. 2 Resonant tidal dissipation spectrum D resulting from inertial modes as a function of the normalised tidal frequency σ = χ/(2ΩA) assuming Nmax = 5.
In the text

	Fig. 3 Evolution of the semi-major axis a over time with a Q factor proportional to inertial wave dissipation in fluids (green curve), and with a constant Q factor (blue dashed curve). The abscissa represents time in years, the vertical axis measures the evolution of the semi-major axis from the initial value a0.
In the text

	Fig. 4 Evolution of the tidal dissipation and of the semi-major axis with time for different values of lp and Hp, the width at half-height and the height of the studied single resonant damping peak. The grey dotted line corresponding to α = 0.2 is superposed to the continuous green one except at the position of the peak for Q-1.
In the text