TIDAL DISSIPATION IN A HOMOGENEOUS SPHERICAL BODY. I. METHODS

The development of the mathematical theory of bodily tides was started by Darwin (1879) who derived a partial sum of the Fourier expansion of the additional potential of a tidally perturbed sphere. A decisive contribution into this theory was offered almost a century later by Kaula (1964), who wrote down a complete series. In a previous paper (Efroimsky & Makarov 2013), we provided a detailed presentation of the Darwin–Kaula expansion, and explained how tidal friction and lagging are built into it. We compared the Darwin–Kaula theory with the one by MacDonald (1964), and demonstrated that the former theory is superior to the latter because it can, in principle, be combined with an arbitrary rheology. Referring the reader to the aforementioned literature for details, we present several central formulae that will be necessary.

An external body of mass M*, located in ${{{\boldsymbol {r}}}}^{\;*} = (r^*,\lambda ^*,\phi ^*)$, generates the following disturbing potential in a point ${{\boldsymbol {{R}}}}= (R,\phi,\lambda)$ on the surface of a sphere of radius R < r*:

$$\begin{eqnarray} W({{\boldsymbol {{R}}}},{{\boldsymbol {r}}}^{ *}) &=& \sum _{{{l}}=2}^{\infty } W_{{l}}({{\boldsymbol {{R}}}}, {{\boldsymbol {r}}}^{ *})\nonumber\\ &=& {-} \frac{G\!M^*}{r^{ \,*}} \sum _{{{l}}=2}^{\infty }\,\left(\frac{R}{r^{ *}}\,\right)^{{{{l}}}}\,P_{{l}}(\cos \gamma) \nonumber\\ &=& {-}\frac{G\!M^*}{r^{\,*}}\sum _{{{l}}=2}^{\infty }\left(\frac{R}{r^{ *}}\right)^{{{{l}}}}\sum _{m=0}^{l} \frac{({l}-m)!}{({l}+m)!}(2-\delta _{0m})P_{lm}\nonumber\\ &&\times\,(\sin \phi)P_{lm}(\sin \phi ^*) \cos m(\lambda -\lambda ^*). \end{eqnarray} \tag{ 1 }$$

Here G denotes Newton's gravity constant, ϕ is the latitude reckoned from the spherical body's equator, λ is the longitude measured from a fixed meridian, and γ is the angular separation between the vectors ${{{\boldsymbol {r}}}}^{\;*}$ and ${{\boldsymbol {{R}}}}$ pointing from the perturbed body's center. The definition and normalization of the Legendre polynomials P_l (cos γ) and the associated Legendre polynomials P_lm (sin ϕ) are given in Appendix A.

While in the above formula the location of the perturber on its trajectory is expressed through the spherical coordinates ${{{\boldsymbol {r}}}}^{\;*} = (r^*,\lambda ^*,\phi ^*)$, a trigonometric transformation (developed by Kaula 1961) enables one to switch to the perturber's orbital elements ${{\boldsymbol {r}}}^{\;*}=(a^*,e^*,{i} ^*,\Omega ^*,\omega ^*,{\cal M}^*)$. Thus, the disturbing potential is expressed as

$$\begin{eqnarray} W({{\boldsymbol {{R}}}},\;{{\boldsymbol {r}}}^{\;*})&=&\sum _{lmpq}\,W_{lmpq}\nonumber\\ &=&{-} \frac{G\!M^*}{a^*}\;\sum _{{l}=2}^{\infty }\left(\frac{R}{a^*}\right)^{{{l}}}\sum _{m=0}^{l}\;\frac{({l} - m)!}{({l} + m)!}(2 -\delta _{0m})P_{{{l}}m}\nonumber\\ &&\times(\sin \phi)\sum _{p=0}^{l}\;F_{{l}mp}({i} ^*)\;\sum _{q= - \infty }^{\infty }\;G_{{l}pq} (e^*)\nonumber\\ &&\times\, \left\lbrace \begin{array}{@{}c@{}} \cos \\ \sin \end{array} \right\rbrace^{{l}-m\, {\rm even}}_{{l}-m\,{\rm odd}} (v_{{l}mpq}^*-m(\lambda +\theta ^*) ), \end{eqnarray} \tag{ 2 }$$

where θ* is the rotation angle of the tidally perturbed body, ³ while F_lmp (i*) and G_lpq (e*) are the inclination functions and the eccentricity polynomials, respectively. The auxiliary linear combinations $v_{lmpq}^*$ are defined by

$$\begin{eqnarray} &&v_{lmpq}^*\;\equiv \;(l-2p)\,\omega ^*+\,(l-2p+q){\cal M}^*+\,m\,\Omega ^*. \end{eqnarray} \tag{ 3 }$$

Conventionally, the letters denoting the elements of the perturber are accompanied with asterisks: $a^*,e^*,{i} ^*,\Omega ^*,\omega ^*,{\cal M}^*$. Following Kaula (1964), the sidereal angle also acquires an asterisk when it appears in a combination $v_{lmpq}^*-m\,\theta ^*$ with the perturber's elements.

The angle θ, however, does not acquire an asterisk when it appears in a linear combination v_lmpq − m θ with the orbital elements of a test body subject to the additional tidal potential of the perturbed body. This strange nomenclature introduced by Kaula (1964)—two different notations for one angle—turns out to be helpful and convenient in the calculation of the back-reaction experienced by the perturber. For comprehensive explanation of this obscure point, see Section 5 in Efroimsky & Makarov (2013).

Over timescales shorter than the apsidal-motion period, the expression in round brackets in the formula (2) can be linearized as

$$\begin{eqnarray} v_{lmpq}^*-m(\lambda +\theta ^*)&=&\omega _{lmpq}\,(t-t_0)\nonumber\\ &&{-} m \lambda + v_{lmpq}^*(t_0) - m \theta ^*(t_0), \end{eqnarray} \tag{ 4 }$$

where the following quantities act as the Fourier tidal modes:

$$\begin{eqnarray} \omega _{lmpq}&\equiv& \dot{v^*_{lmpq}}- m\,\dot{\theta ^*} =\;(l-2p)\;\dot{\omega }^*+\,(l-2p+q)\;{{\bf \dot{\cal {M}}}}^{\,*}\nonumber\\ &&+\,m(\dot{\Omega }^*-\dot{\theta }^*), \end{eqnarray} \tag{ 5 }$$

${{\bf \dot{\cal {M}}}}^{\,*}$ being the perturber's "anomalistic" mean motion (see Section 2.3 below), and t₀ being the time of pericenter passage. (As ever, we set ${\cal {M}}^{\,*}=\,0$ in the pericenter.) The modes ω_lmpq can assume either sign, but the physical forcing frequencies are positive definite:

$$\begin{eqnarray} &&\chi _{lmpq}=|\omega _{lmpq}|. \end{eqnarray} \tag{ 6 }$$

In the preceding subsection, we obeyed the convention by Kaula (1964) and marked with asterisk the orbital elements of the tide-raising body. Kaula introduced this notation because within his model he also considered another exterior body, which was disturbed by the tides generated on the planet by the tide-raising body. This exterior body's elements were denoted by the same letters, but without an asterisk.

When the two outer bodies coincide, the asterisks may be dropped, except on two occasions. The first is writing the masses—while the mass of the planet is denoted with M, the mass of the perturber (the star) will be written as M*. The other occasion requires writing the additional tidal potential of the perturbed body—the potential will have a value $U({{\boldsymbol {r}}},{{\boldsymbol {r}}}^{\,*})$ in a point ${{\boldsymbol {r}}}$, provided the perturber (the star) is located in an exterior point ${{\boldsymbol {r}}}^{\,*}$ (both vectors being planetocentric). The planet's rotation rate θ, as well as the orbital elements of the star as seen from the planet, will hereafter be written without asterisks.

The most important notations employed in this paper are collected in Table 1.

Table 1. Symbol Key

Notation	Description
${{\boldsymbol {r}}}^{\,*}$	The position of the star relative to the center of the planet
r *	The star–planet distance
ϕ *	The declination of the star relative to the equator of the planet
λ *	The right ascension of the star relative to a fixed meridian on the planet
${{\boldsymbol {{R}}}}$	A point on the surface of the planet
${{\boldsymbol {r}}}$	A point outside the planet, located above the surface point ${{\boldsymbol {{R}}}}$
R	The radius of the planet
r	Distance from the center of the planet to an exterior point ${{\boldsymbol {r}}}$
ϕ	The latitude of the point ${{\boldsymbol {{R}}}}$ on the surface of the planet
	(also the declination of an exterior point ${{\boldsymbol {r}}}$ located above the surface point ${{\boldsymbol {{R}}}}$)
λ	The longitude of a point ${{\boldsymbol {{R}}}}$ on the surface of the planet
	(also the right ascension of an exterior point ${{\boldsymbol {r}}}$ located above the surface point ${{\boldsymbol {{R}}}}$)
$W({{\boldsymbol {{R}}}},{{\boldsymbol {r}}}^{\,*})$	Tide-raising potential at a surface point ${{\boldsymbol {{R}}}}$ of the planet
$W_l({{\boldsymbol {{R}}}},{{\boldsymbol {r}}}^{\,*})$	The l-degree part of the tide-raising potential at a surface point ${{\boldsymbol {{R}}}}$ of the planet
$U({{\boldsymbol {r}}},{{\boldsymbol {r}}}^{\,*})$	Additional tidal potential in a point ${{\boldsymbol {r}}}$ outside the planet
$U_l({{\boldsymbol {r}}},{{\boldsymbol {r}}}^{\,*})$	The l-degree part of the additional tidal potential in a point ${{\boldsymbol {r}}}$ outside the planet
V 0	The constant-in-time spherically symmetrical potential of an undeformed planet
V'	The total perturbation of the potential of the planet (V' = W + U)
V	The overall potential of the planet (V = V 0 + V' = V 0 + W + U)
ωlmpq	The Fourier modes of the tide
χlmpq	The physical frequencies of the tidal stresses and strains (χlmpq = |ωlmpq |)
χ	A shorter notation for χlmpq
M	The mass of the planet
M *	The mass of the star
a	The semimajor axis
e	The eccentricity
i	The obliquity (the inclination of the star as seen from the planet)
ω	The argument of the pericenter of the star as seen from the planet
ω	When there is no risk of confusion, ω is also used as a short notation for ωlmpq
Ω	The longitude of the node of the star as seen from the planet
${\cal {M}}$	The mean anomaly of the star as seen from the planet
n	The mean motion
G	Newton's gravitational constant
θ	The rotation angle of the planet
$\dot{\theta \,}$	The spin rate of the planet
kl , hl	The degree-l static Love numbers of the planet
kl (ωlmpq ), hl (ωlmpq )	The degree-l dynamical Love numbers of the planet
l (ωlmpq )	The degree-l tidal phase lag
Δtl (ωlmpq )	The degree-l tidal time lag
Ql (ωlmpq )	The degree-l tidal quality factor defined as Ql (ωlmpq ) ≡ 1/sin |l (ωlmpq )|
Qlmpq	The Q factor, in the notation of Kaula (1964) and Peale & Cassen (1978)
Flmp (i)	The inclination functions
Glpq (e)	The eccentricity polynomials
μ	The mean rigidity of the mantle
J	The mean compliance of the mantle (J = 1/μ)
J	When there is no risk of confusion, J also denotes the Jacobian
ρ	The mean density of the planet
g	The surface gravity on the planet
u	Tidal displacement in the planet
v	The velocity of tidal displacement in the planet
P	The power exerted on the planet by the tidal stresses
〈P〉	The time-averaged power (averaged over one or several cycles of flexure)

Download table as:  ASCIITypeset images: 1 2

At this point, a word of warning is necessary. Deriving Equation (5), we differentiated the expression (3), which gave us the terms with $\dot{\omega }$ and $\dot{\Omega }$. Including these terms in Equation (5) acknowledges the fact that the perturber's trajectory is disturbed, not Keplerian. The disturbance may come solely from tides, as in Kaula (1964, Equation (38)), or from both tides and other sources. One way or another, the perturber's mean anomaly ${\cal {M}}$ is no longer equal to nt, but is now given by $\, {\cal {M}}= {\cal {M}}_0+ \int ^{\,t}_{t_0} n(t) dt$, where $n(t)\equiv \sqrt{G\,(M+\,M^*) a^{-3}(t)\,}.$ Accordingly, the expression for the modes becomes:

$$\begin{eqnarray} \omega _{lmpq}&\equiv& \dot{v}_{lmpq}- m\,\dot{\theta } =\;(l-2p)\;(\dot{\omega }+\,\dot{\cal {M}}_0)\nonumber\\ &&+\,q\,\dot{\cal {M}}_0 +\,(l-2p+q)\;n+\,m\;(\dot{\Omega } - \dot{\theta }). \end{eqnarray} \tag{ 7 }$$

It is, of course, tempting to assume that $\dot{{\cal {M}}_0}\ll n$, thus accepting the approximation

$$\begin{eqnarray} &&\dot{\cal {M}\,}= \dot{{\cal {M}}_0}+ n \approx n, \end{eqnarray} \tag{ 8 }$$

as Kaula (1964) did in his Equations (46) and (47). Within his theory, however, this approximation could not be used. ⁴ This is explained in Appendix B where we consider two examples. One is the case where perturbation of an orbit of a moon is mainly due to the tides the moon creates in the planet. In that situation, $\dot{\omega }$ and $\dot{\cal {M}}_0$ are of the same order but of opposite signs, so they largely compensate one another. This suggests a simultaneous neglect of both rates. The second example is when the dominant perturbation of the orbit comes from the oblateness of the primary. In this case $\dot{\omega }$ and $\dot{\cal {M}}_0$ are of the same order and the same sign—so keeping one of these terms requires keeping the other.

Whether one or both of these rates should be included depends on a particular setting, and each practical case must be examined separately. In general, both rates should be kept.

While keeping $\dot{\omega }$ complicates the formalism, the emergence of $\dot{{\cal {M}}_0}$ complicates the treatment even further. To sidestep this issue, we shall define the mean motion via

$$\begin{eqnarray} &&n\equiv {{\bf \dot{\cal {M}}}}. \end{eqnarray} \tag{ 9 }$$

This, the so-called anomalistic mean motion differs from $\sqrt{G(M+M^*)\,a^{-3}(t)\,}$.

We shall derive the heat-production formulae for two different settings—with a fixed pericenter ω and with ω moving uniformly.

For a static tide, the incremental tidal potential of the perturbed body mimics the perturbation (Equation (2)), except that each term W_l is now equipped with a mitigating multiplier k_l (R/r)^{l + 1}, where k_l is an l-degree Love number. With the star located in ${{\boldsymbol {r}}}^{\,*}$, the additional potential in a point ${{\boldsymbol {r}}}$ will read as

$$\begin{eqnarray} U({{\boldsymbol {r}}},\;{{\boldsymbol {r}}}^{\;*})&=&\sum _{{l}=2}^{\infty } U_{{l}}({{\boldsymbol {r}}}) \nonumber\\ &=& \sum _{{l}=2}^{\infty } k_{l}\;\left(\frac{R}{r}\,\right)^{{l}+1}\;W_{{l}}({{\boldsymbol {{R}}}},\;{{\boldsymbol {r}}}^{\;*}). \end{eqnarray} \tag{ 10 }$$

For time-dependent tides, this expression acquires an extra amendment: the reaction must lag, compared to the action. Naively, this would imply taking each W_l at an earlier instant of time. However, in reality lagging depends on frequency, so each W_l must be first expanded into a Fourier series over tidal modes, whereafter each term of the series should be delayed separately. The magnitude of the tidal reaction is frequency dependent too, so each term of the Fourier series will be multiplied by a dynamical Love number of its own. Symbolically, this may be written in a manner similar to the static expression:

$$\begin{eqnarray} U({{\boldsymbol {r}}},\;{{\boldsymbol {r}}}^{\;*})&=&\sum _{{l}=2}^{\infty } U_{{l}}({{\boldsymbol {r}}},\;{{\boldsymbol {r}}}^{\;*})\nonumber\\ &=& \sum _{l=2}^{\infty } \left(\frac{R}{r}\,\right)^{{l}+1}\;\hat{k}_{l}\;W_{{l}}({{\boldsymbol {{R}}}},\;{{\boldsymbol {r}}}^{\;*}). \end{eqnarray} \tag{ 11 }$$

The hat above $\hat{k}_{l}$ means that this is not a multiplier but a linear operator that mitigates and delays each Fourier mode of W_l differently:

$$\begin{eqnarray} \nonumber U({{\boldsymbol {r}}},\;{{\boldsymbol {r}}}^{\;*})&=& {-} \frac{G\!M^*}{a}\;\sum _{{l}=2}^{\infty }\;\left(\frac{R}{r}\,\right)^{{{l+1}}} \left(\frac{R}{a}\,\right)^{{{l}}}\sum _{m=0}^{l}\;\frac{(l - m)!}{(l + m)!}\nonumber\\ &&\times\, (2 -\delta _{0m})P_{lm}(\sin \phi) \sum _{p=0}^{l}\;F_{lmp}(i) \nonumber \\ &&\times\,\sum _{q= - \infty }^{\infty }\;G_{lpq}(e) k_l(\omega _{lmpq}) \left\lbrace \begin{array}{@{}c@{}}\cos \\ \sin \end{array} \right\rbrace ^{{l} - m\, {\rm even}}_{l - m\, {\rm odd}}\nonumber\\ &&\times\, (v_{lmpq} -m\,(\lambda + \theta) - \epsilon _{l} ), \end{eqnarray} \tag{ 12 }$$

where the Love numbers k_l (ω_lmpq ) and the phase lags _l (ω_lmpq ) are functions of the Fourier modes. The lags emerge as the products

$$\begin{eqnarray} &&\epsilon _l(\omega _{lmpq}) = \omega _{lmpq} \Delta t_l(\omega _{lmpq}), \end{eqnarray} \tag{ 13a }$$

where Δt_l (ω_lmpq ) is the time delay at the mode ω_lmpq . In reality, the time delays are functions not of the Fourier modes (which can assume either sign), but of the actual physical forcing frequencies χ_lmpq = |ω_lmpq |, which are positive definite. Thus it is more accurate to write the delays not as Δt_l (ω_lmpq ) but as Δt_l (χ_lmpq ). Accordingly, the phase lags become

$$\begin{eqnarray} &&\epsilon _l(\omega _{lmpq}) = \omega _{lmpq} \Delta t_l(\chi _{lmpq}). \end{eqnarray} \tag{ 13b }$$

For the reasons of causality, the time delays are positive definite, so the sign of the phase lag always coincides with that of the corresponding Fourier mode. Thus we finally have:

$$\begin{eqnarray} \epsilon _l(\omega _{lmpq}) &=& |\omega _{lmpq}| \,{\rm Sgn}\,(\omega _{lmpq}) \,\Delta t_l(\chi _{lmpq})\nonumber\\ &=& \chi _{lmpq} \,{\rm Sgn}\,(\omega _{lmpq}) \,\Delta t_l(\chi _{lmpq}), \end{eqnarray} \tag{ 13c }$$

where χ_lmpq ≡ |ω_lmpq | are the positive definite forcing frequencies.

The dynamical Love number k_l (ω_lmpq ) and the phase lag _l (ω_lmpq ) are the absolute value and the negative phase of the complex Love number $\bar{k}_l(\omega _{lmpq})$ whose functional dependence upon the Fourier mode is solely determined by l, provided the body is spherical. ⁵

As we saw above, to obtain the decomposition (12) from the Fourier series (2), each lmpq term of the latter had to be endowed with its own mitigating factor k_l = k_l (ω_lmpq ) and phase lag _l = _l (ω_lmpq ). In the past, some authors enquired whether this mitigate-and-lag method is general enough to describe tides. It is, as long as the tides are linear. This is explained in Appendix C.

The expression (12) for the additional tidal potential contains both sines and cosines of the phase lags, and so does the ensuing expression for the surface elevation. However, the resulting expression for the tidal dissipation rate turns out to contain only the combination k_l (ω) sin _l (ω) which is the negative imaginary part of the complex Love number:

$$\begin{eqnarray} &&{k}_{ {l}}(\omega)\;\sin \epsilon _{ {l}}(\omega)=|\bar{k}_{ {l}}(\omega)|\;\sin \epsilon _{ {l}}(\omega)= -{\cal {I}}{{m}}[\bar{k}_{ {l}}(\omega)],\nonumber\\ && {\rm where}\quad \omega =\omega _{lmpq}. \end{eqnarray} \tag{ 14 }$$

This quantity is often denoted as k_l /Q, although it would be more reasonable to employ the notation k_l /Q_l , with the tidal quality factors defined through 1/Q_l ≡ |sin _l |.

A dynamical Love number k_l (ω_lmpq ) is an even function of the tidal mode ω_lmpq , while a phase lag _l (ω_lmpq ) is odd, as can be observed from Equation (13b). Thus the expression for the product k_l  sin _l as a function of the physical frequency χ = χ_lmpq ≡ |ω_lmpq | is:

$$\begin{eqnarray} &&{k}_{ {l}}(\omega)\;\sin \epsilon _{ {l}}(\omega)\;=\;{k}_{ {l}}(\chi)\;\sin \epsilon _{ {l}}(\chi)\;\,{\rm Sgn}\,(\omega), \end{eqnarray} \tag{ 15 }$$

where _l (χ) is non-negative, because the physical frequency is χ.

The frequency dependence of k_l /Q_l = k_l (χ) sin _l (χ) is defined by two major physical circumstances: self-gravitation of the planet and the rheology of its mantle. A rheological law is expressed by a constitutive equation, that is, by an equation interconnecting the strain and the stress. A particular form of this equation is determined by the friction mechanisms present in the considered medium. A realistic rheological law should contain contributions from elasticity, viscosity, and inelastic processes (mainly, dislocation unjamming). Self-gravitation suppresses the tidal bulge. At low frequencies this effectively adds to the mantle's rigidity, whereas at higher frequencies the interplay of rheology and gravity is more complex (Efroimsky 2012b, Figure 2).

The calculation of the frequency dependence k_l (χ) sin _l (χ) for a homogeneous body of a known size, mass, and rheology is presented in detail in Efroimsky (2012a, 2012b). See also the Appendix to Paper II.

While quadrupole (l = 2) terms are sufficient in most problems, exceptions are known. For the orbital evolution of Phobos, the l = 3 and, possibly, even l = 4 terms of the Martian tidal potential may be of relevance (Bills et al. 2005). Studying close binary asteroids, Taylor & Margot (2010) took into account the Love numbers up to l = 6.

The question of how rapidly l > 2 terms fall off with the increase of the degree l is also interesting. Most authors only rely on the geometric factor (R/a)^{2l + 1} to answer this question. As was explained in Efroimsky (2012b), the l-dependence of k_l (ω_lmpq ) sin (ω_lmpq ), too, comes into play and changes the result considerably.

To compare the varying shape of a deformable body against some benchmark configuration, we use ${{\boldsymbol {X}}}$ to denote the initial position occupied by a particle at t = 0. At another time t, the particle finds itself in a new place

$$\begin{eqnarray*} &&{{\boldsymbol {X}}}= {\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), \nonumber \end{eqnarray*}$$

where the function ${\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)$ is a trajectory, that is, a solution to the equation of motion, with the initial condition ${{\boldsymbol {X}}}={{\boldsymbol {X}}}$ set at t = 0.

The current values of all physical and kinematic properties of the medium can be expressed as functions of the instantaneous coordinates ${{\boldsymbol {X}}}$ of a point where these properties are being measured at the present moment t. When referred to the present time and position, such properties are named Eulerian and are equipped with a subscript E; for example: $q_{E}({{\boldsymbol {X}}},t)$. The Eulerian description is fit to answer the question "where," and therefore is convenient in fluid dynamics where the displacement ${{\boldsymbol {X}}}-{{\boldsymbol {X}}}$ can become arbitrarily large and the initial position ${{\boldsymbol {X}}}$ is soon forgotten.

While ${{\boldsymbol {X}}}$ denotes a place in space, the initial condition ${{\boldsymbol {X}}}$ acts as the "number" of a particle presently residing at the place ${{\boldsymbol {X}}}$. Although located at ${{\boldsymbol {X}}}$, the particle originally came from ${{\boldsymbol {X}}}$ and will carry the label ${{\boldsymbol {X}}}$ forever.

Knowing the trajectories of all particles, we can express the properties as functions of the time t and the initial conditions ${{\boldsymbol {X}}}$. To that end, we employ the change of variables x = f  ( X , t). Expressed through the initial conditions, a property q will be termed as Lagrangian and equipped with the subscript L:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}},t) \equiv q_{E}({{\boldsymbol {X}}}, t)\, \end{eqnarray} \tag{ 16a }$$

or, in more detail:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}},t) \equiv q_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t). \end{eqnarray} \tag{ 16b }$$

So q_L has the same value as q_E , but has a different functional form, as it is now understood as a function of the initial conditions (the particles' "numbers") ${{\boldsymbol {X}}}$, and not of the present-time coordinates ${{\boldsymbol {X}}}$. Relating the quantities to the initial positions ${{\boldsymbol {X}}}$, the Lagrangian description tells us "which particle", and is thus practical in description of deformable solids.

In anticipation of perturbative treatment, we regard the trajectory x = f ( X , t) as fiducial and equip the appropriate functional dependencies with a superscript 0:

$$\begin{eqnarray} &&q^0_{L}({{\boldsymbol {X}}},t) \equiv q^0_{E}({{\boldsymbol {X}}}, t) \end{eqnarray} \tag{ 17a }$$

which is:

$$\begin{eqnarray} &&q^0_{L}({{\boldsymbol {X}}},t) \equiv q^0_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t). \end{eqnarray} \tag{ 17b }$$

Under disturbance, two changes will take place in a point ${{\boldsymbol {r}}}$ at a time t:

• 1.  
  Properties will now assume different values in this point at this time. So we substitute the unperturbed Eulerian dependencies $q^0_{E}({{\boldsymbol {r}}}, t)$ with

  $$\begin{eqnarray} &&q_{E}({{\boldsymbol {r}}}, t) = q^0_{E}({{\boldsymbol {r}}}, t) + q^{\prime }_{E}({{\boldsymbol {r}}}, t). \end{eqnarray} \tag{ 18 }$$

  This equality, in fact, serves as a definition of the variation $q^{\prime }_{E}({{\boldsymbol {r}}}, t)$: the variation is a change in the functional dependence of a physical property upon the present position ${{\boldsymbol {r}}}$.
• 2.  
  A different particle will now appear in the point ${{\boldsymbol {r}}}$ at the time t. It will not be the same particle as the one expected there at the time t in the absence of perturbation.Accordingly, a particle, which starts in ${{\boldsymbol {X}}}$ at t = 0, will show up at the time t, not in the point ${{\boldsymbol {X}}}={\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)$ but in some other location displaced by u:

  $$\begin{eqnarray} &&{{\boldsymbol {r}}}= {{\boldsymbol {X}}}+ {{\bf u}}= {\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t) + {{\bf u}}({{\boldsymbol {X}}},t). \end{eqnarray} \tag{ 19 }$$

Both of these changes, 1 and 2, will affect the Lagrangian dependencies of the properties upon the initial conditions, so the dependency of each property will acquire a variation $q^{\prime }_{L}({{\boldsymbol {X}}},t)$:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}},t) = q^0_{L}({{\boldsymbol {X}}},t) + q^{\prime }_{L}({{\boldsymbol {X}}},t). \end{eqnarray} \tag{ 20 }$$

In Appendix D, we provide a self-sufficient introduction into the perturbative treatment of a deformable body, both in the Eulerian and Lagrangian languages. There we derive a relation between the perturbations of the Lagrangian and Eulerian quantities:

$$\begin{eqnarray} &&q^{\prime }_{L}({\boldsymbol X},t) = q^{\prime }_{E}({\boldsymbol f}({\boldsymbol X},t), t) + {{\bf u}}({\boldsymbol X},t) \nabla _{x} q_{E} + O({{\bf u}}^2),\quad \end{eqnarray} \tag{ 21 }$$

with the gradient in the second term acting on the unperturbed history: ⁶

$$\begin{eqnarray} &&\nabla _{x} q_{E}\equiv \nabla _{x}q_{E}({\boldsymbol x},t) ,\quad {\rm where}\qquad {\boldsymbol x}={\boldsymbol f}({\boldsymbol X},t). \end{eqnarray} \tag{ 22 }$$

In the formula (21), the first term on the right-hand side, $q^{\prime }_{E}({{\boldsymbol {X}}},t)$, accounts for the change of the final spatial distribution of properties. The other two terms show up because perturbation alters the mapping from ${{\boldsymbol {X}}}$ to the present position.

We need several formulae for density perturbations, which are obtained in Appendix D.

In the Eulerian description:

$$\begin{eqnarray} &&\rho _{E}({{\boldsymbol {r}}},t)= \rho ^{\,0}_{E}({{\boldsymbol {r}}}) + \rho _{E}^{\prime } ({{\boldsymbol {r}}},t), \qquad \quad \end{eqnarray} \tag{ 23 }$$

$$\begin{eqnarray} &&\rho _{E}^{\prime }+\,\nabla _{ {r\,}}\cdot \left(\rho ^{\,0}_{E}\,{{\bf u}}\right) = 0, \qquad \quad \end{eqnarray} \tag{ 24 }$$

Formula (23) renders the interrelation between the functions of the same variable. The unperturbed density $\rho ^{\,0}_{E}$ appears here as a function of the perturbed present positions ${{\boldsymbol {r}}}$, not of the unperturbed reference positions ${{\boldsymbol {X}}}$. This can be traced through the derivation (D21)–(D24). There, the unperturbed density initially shows up as a function of ${{\boldsymbol {X}}}={{\boldsymbol {r}}}-{{\bf u}}$. It then ends up as a function of ${{\boldsymbol {r}}}$, after the Taylor expansion around ${{\boldsymbol {r}}}$ over powers of u is performed.

Accordingly, the symbol ∇_r denotes differentiation with respect to the perturbed position ${{\boldsymbol {r}}}$ upon which $\rho ^{\,0}_{E}$ is set to depend in the above equations. Also remember that in ${{\bf u}}({{\boldsymbol {X}}},t)={{\bf u}}({{\boldsymbol {r}}},t)+\,O({{\bf u}}^2)$ we can neglect O(u²), in the linear approximation. Thus the Lagrangian and Eulerian values of the displacement coincide in the first order. Specifically, in Equation (24), our u can be treated as a function of ${{\boldsymbol {r}}}$. So all entities in that equation are functions of the same variable, the perturbed location.

In the Lagrangian description:

$$\begin{eqnarray} &&\rho _{L}({{\boldsymbol {X}}},t)= \rho ^{\,0}_{E}({{\boldsymbol {X}}}) + \rho _{L}^{\prime }({{\boldsymbol {X}}},t), \quad \quad \end{eqnarray} \tag{ 25 }$$

$$\begin{eqnarray} &&\rho _{L}^{\prime }+\,\rho ^{\,0}_{E}\,\nabla _{ {X\,}}\cdot {{\bf u}} = 0. \qquad \quad \end{eqnarray} \tag{ 26 }$$

Recall that this is an interrelation between functions of the same variable. This time, it is the initial position ${{\boldsymbol {X}}}$. Had we altered the notation from ${{\boldsymbol {X}}}$ to ${{\boldsymbol {r}}}$, nothing would have changed (except that we would write ∇_r instead of ∇_X )—it would still be the same relation between three functions of the same argument.

Relation between the increments $\rho _{L}^{\prime }$ and $\rho _{E}^{\prime }$:

This relation originates from the general formula (21). In our case, the reference trajectory ${{\boldsymbol {X}}}={\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)$ stays identical to the initial position ${{\boldsymbol {X}}}$, so we obtain:

$$\begin{eqnarray} &&\rho _{L}^{\prime }({{\boldsymbol {X}}},t) = \rho _{E}^{\prime }({{\boldsymbol {X}}},t) + {{\bf u}}({{\boldsymbol {X}}},t)\cdot \nabla _{ {x\,}}\rho ^{\,0}_{E}({{\boldsymbol {X}}},t). \end{eqnarray} \tag{ 27 }$$

Once again, we are dealing with a relation between several functions taken all at one and the same point. Here the point is denoted with ${{\boldsymbol {X}}}$. Had we denoted it with ${{\boldsymbol {r}}}$, the only change would be a switch from ∇_x to ∇_r , no matter what meaning we instill into these ${{\boldsymbol {X}}}$ and ${{\boldsymbol {r}}}$.

For an initially homogeneous body, $\nabla \rho ^{\,0}_{E} = 0$, so the forms (24) and (26) of the linearized conservation law coincide and can both be conveniently written as

$$\begin{eqnarray} &&\rho ^{\,0} \nabla _{ {r}}\cdot {{{\bf v}}}+\,\frac{\textstyle \partial \rho }{\textstyle \partial t\,} = 0, \end{eqnarray} \tag{ 28 }$$

where $\rho ^{\,0}\equiv \rho ^{\,0}_{E}$ and the velocity is

$$\begin{eqnarray} &&{{\bf v}}= \frac{\partial {{\bf u}}}{\partial t}. \end{eqnarray} \tag{ 29 }$$

In each point, the density ρ and potential V comprise a mean value and a perturbation:

$$\begin{eqnarray} {\rm density:}\qquad \rho &=&\rho ^{\,0} + \rho ^{\prime }, \end{eqnarray} \tag{ 30a }$$

$$\begin{eqnarray} {\rm potential:}\qquad V&=&V^{\,0}+ V^{\prime }= V^{\,0}+ (W + U), \end{eqnarray} \tag{ 30b }$$

where V⁰ is the constant-in-time spherically symmetrical potential of an undeformed body, whereas V' denotes the potential's perturbation. The perturbation consists of the external tide-raising potential W and the resulting additional potential U of the perturbed body:

$$\begin{eqnarray} &&V^{\prime } = W + U. \end{eqnarray} \tag{ 31 }$$

The potentials and densities will be endowed with a subscript "L" or "E" pointing at the Lagrangian or Eulerian descriptions, accordingly. Owing to the general expression (21), we have:

$$\begin{eqnarray} &&{V_{E}^{\prime }} ({{\boldsymbol {r}}},t) = {V_{L}}^{\prime }({{\boldsymbol {X}}},t) - {{\bf u}}\cdot \nabla _{ {x\,}} V^{\,0}_{E}, \end{eqnarray} \tag{ 32 }$$

the same being valid for ρ, see Equation (27). For unperturbed properties, however, subscripts may be dropped without causing any confusion:

$$\begin{eqnarray} &&V^{\,0} \equiv V^{\,0}_{E},\quad \rho ^{\,0} \equiv \rho ^{\,0}_{E}. \end{eqnarray} \tag{ 33 }$$

In both the perturbed and unperturbed settings, the density and potential are always linked through the Poisson equation:

$$\begin{eqnarray} \nabla _{ {r\,}}^{\,2}\,V_{E}&=&{-} 4\,\pi \,G\,\rho _{E}, \end{eqnarray} \tag{ 34a }$$

$$\begin{eqnarray} \nabla _{ {r\,}}^{\,2}\,V^{\,0}_{E}&=&{-} 4\,\pi \,G\,\rho ^{\,0}_{E}, \end{eqnarray} \tag{ 34b }$$

while the perturbing potential W obeys the Laplace equation outside the perturber:

$$\begin{eqnarray} \nabla _{ {r\,}}^{\,2}\,W_{E}&=&0. \qquad \qquad \end{eqnarray} \tag{ 34c }$$

Subtraction of Equation (34b) from Equation (34a) results in a Poisson equation for the density perturbation:

$$\begin{eqnarray} &&\left. \quad \right. \nabla _{ {r\,}}^{\,2}\,{V_{E}^{\prime} } = - 4 \pi G {\rho _{E}^{\prime }}. \end{eqnarray} \tag{ 35 }$$

The Poisson equation in the Lagrangian description is presented in Appendix D.

The power P exerted on the perturbed body is an integral, over its volume, of the rate of working by tidal forces on displacements. In the Eulerian language, the power reads as

$$\begin{eqnarray} &&P\;=\; \int \,\rho _{E}\;{{{\bf v}}}\,\cdot \,\nabla _{ {r\,}} {V_{E}^{\prime }}\;d^3r, \end{eqnarray} \tag{ 36 }$$

the integration being performed over an instantaneous, deformed volume. Together with

$$\begin{eqnarray} &&\rho _{E} {{{\bf v}}}\cdot \nabla _{ {r\,}} {V_{E}^{\prime }} = \nabla _{ {r\,}}\cdot (\rho _{E}\,{{{\bf v}}} {V_{E}^{\prime }}) - {V_{E}^{\prime }} \nabla _{ {r\,}}\cdot (\rho _{E}\,{{{\bf v}}}), \end{eqnarray} \tag{ 37 }$$

the mass-conservation law

$$\begin{eqnarray} &&\nabla _{ {r\,}}\cdot (\rho _{E}\,{{\bf v}})+\,\frac{\textstyle \partial \rho _{E}}{\textstyle \partial t\,} = 0\, \end{eqnarray} \tag{ 38 }$$

simplifies the expression under the integral to the following form:

$$\begin{eqnarray} &&\rho _{E} {{{\bf v}}}\cdot \nabla _{ {r\,}} {V_{E}^{\prime }} = \nabla _{ {r\,}}\cdot (\rho _{E}\,{{{\bf v}}}\,{V_{E}^{\prime }})+\,{V_{E}}^{\prime } \frac{\,\partial {\rho _{E}}^{\prime }}{\partial t\,} \;\;. \end{eqnarray} \tag{ 39a }$$

Further employment of the Poisson equation in the Eulerian form, (Equation (35)), gives us

$$\begin{eqnarray} \rho _{E} {{{\bf v}}}\cdot \nabla _{ {r\,}} {V_{E}^{\prime }} &=& \nabla _{ {r\,}}\cdot (\rho _{E}\,{{{\bf v}}}\,{V_{E}^{\prime }}) - \frac{1}{4\,\pi \,G} {V_{E}^{\prime }} \frac{\partial }{\partial t}\,\nabla _{ {r\,}}^2 {V_{E}}^{\prime }. \nonumber\\ \end{eqnarray} \tag{ 39b }$$

So the power becomes

$$\begin{eqnarray} P&=& \int \,\nabla _{ {r\,}}\cdot (\rho _{E}\,{{{\bf v}}}\,{V_{E}^{\prime }})\;d^3r + \int \,{V_{E}^{\prime }} \frac{\,\partial {\rho _{E}^{\prime }}}{\partial t\,} d^3r \end{eqnarray} \tag{ 40a }$$

$$\begin{eqnarray} &=& \int _{\Sigma ^{t}}\,\rho _{E} {V_{E}}^{\prime }\;{{{\bf v}}}\cdot d{{\bf S}}^t\nonumber\\ &&- \frac{1}{4\,\pi \,G} \int \,{V_{E}}^{\prime } \frac{\partial }{\partial t\,}\,\nabla _{ {r\,}}^2 {V_{E}}^{\prime } d^3r , \end{eqnarray} \tag{ 40b }$$

where $d{{\bf S}}^t\,\equiv {{\bf \hat{n}}}^{t}\;d\Sigma ^{t}$, with ${{\bf \hat{n}}}^{t}$ and dΣ^t being a unit normal to the deformed surface and an element of area on that surface, both taken at the time t. Correct to the first order in the displacement u, these are related to their unperturbed analogues via

$$\begin{eqnarray} {{\bf \hat{n}}}^{t}&=&(1 - \nabla ^{\Sigma }\otimes {{\bf u}})\,{{\bf \hat{n}}}\qquad {\rm and} \qquad d\Sigma ^{t}=(1 + \nabla ^{\Sigma }\cdot {{\bf u}})\,d\Sigma, \nonumber\\ && \end{eqnarray} \tag{ 41 }$$

where the surface gradient is defined as

$$\begin{eqnarray} &&\nabla ^{\Sigma }\equiv \nabla _{ {x\,}}- {{\bf \hat{n}}} \partial _{{{\bf \hat{n}}}}, \end{eqnarray} \tag{ 42 }$$

so ∇^Σ⊗u is a three-dimensional, second-rank tensor (Dahlen & Tromp 1998). Altogether,

$$\begin{eqnarray} d{{\bf S}}^t &\equiv& {{\bf \hat{n}}}^{t}\;d\Sigma ^{t} = (1+\,\nabla ^{\Sigma }\cdot {{\bf u}})\,{{\bf \hat{n}}} d\Sigma - (\nabla ^\Sigma \otimes {{\bf u}})\,{{\bf \hat{n}}} d\Sigma\nonumber\\ &=& (1+\,\nabla ^{\Sigma }\cdot {{\bf u}})\,d{\bf S} - (\nabla ^\Sigma \otimes {{\bf u}})\,d{\bf S} , \, \end{eqnarray} \tag{ 43a }$$

with $d{\bf S}\equiv {{\bf \hat{n}}}d\Sigma$ pertaining to the unperturbed surface. In a shorter form, the above reads as

$$\begin{eqnarray} &&d{{\bf S}}^t= {\mathbb {J}} d{\bf S}, \end{eqnarray} \tag{ 43b }$$

where the three-dimensional, second-rank tensor

$$\begin{eqnarray} &&{\mathbb {J}} \equiv (1+\,\nabla ^{\Sigma }\cdot {{\bf u}})\,{\mathbb {I}} - \nabla ^\Sigma \otimes {{\bf u}} \end{eqnarray} \tag{ 44 }$$

is, loosely speaking, playing the role of a Jacobian for elements of area. This is fully analogous to the formula

$$\begin{eqnarray} d^3r &=& J\,d^3x = (1+\,\nabla _{ {x}}\cdot {{\bf u}})\,d^3x = [1+\,\nabla _{ {r}}\cdot {{\bf u}}+\,O(u^2)]\,d^3x \nonumber\\ && \end{eqnarray} \tag{ 45 }$$

linking the deformed volume d³ r to the undeformed volume d³ x (see Appendix D.4.1.).

Applied to the density, the general formula (16) renders:

$$\begin{eqnarray} &&\rho _{L}({{\boldsymbol {X}}},t)\equiv \rho _{E}({{\boldsymbol {r}}},t). \end{eqnarray} \tag{ 46 }$$

This, together with the formula (32) for the potential perturbation, enables us to express the power in the Lagrangian description:

$$\begin{eqnarray} &&P\;=\; \int \,\rho _{L}\;{{{\bf v}}}\,\cdot \,\nabla _{ {x\,}}{V_{L}}^{\prime }\;d^3x - \int \,\rho _{L}\;{{{\bf v}}}\,\cdot \,\nabla _{ {x\,}}({{\bf u}}\cdot \nabla _{ {x\,}} V^0)\;d^3x, \nonumber\\ \end{eqnarray} \tag{ 47 }$$

the integral now being taken over the undeformed body. Be mindful that d³ r ∇_r = d³ x ∇_x , so no Jacobian shows up on the right-hand side.

The velocity and displacement being in quadrature, the second term should be dropped after time averaging (denoted with angular brackets):

$$\begin{eqnarray} &&\langle P\rangle \;=\; \int \,\rho _{L}\;{{{\bf v}}}\,\cdot \,\nabla _{ {x\,}}{V_{L}}^{\prime }\;d^3x, \end{eqnarray} \tag{ 48a }$$

For a periodically deformed solid, we set the equilibrium state to play the role of the unperturbed configuration, for which reason ⁷ $\rho _{L}({\boldsymbol X}, t) J = \rho ^{0}_{E}({\boldsymbol x})$. Insertion of this equality into the expression (48a) gives us:

$$\begin{eqnarray} &&\langle P\rangle = \int \rho ^{0}{{{\bf v}}}\,\cdot \,\nabla _{{x}}{V_{L}}^{\prime }\;J^{-1}\,d^3x. \end{eqnarray} \tag{ 48b }$$

The dot-product can be easily rearranged via the formulae analogous to Equations (37)–(39). Due to

$$\begin{eqnarray} &&\rho ^{\,0} {{{\bf v}}}\cdot \nabla _{ {x\,}} {V_{L}}^{\prime } = \nabla _{ {x\,}}\cdot (\rho ^{\,0}\,{{{\bf v}}}\,{V_{L}}^{\prime }) - {V_{L}}^{\prime }\,\nabla _{ {x\,}}\cdot ({{{\bf v}}}\,\rho ^{\,0}) \end{eqnarray} \tag{ 49 }$$

and

$$\begin{eqnarray} &&\rho ^{\,0} \nabla _{ {x\,}}\cdot {{{\bf v}}}+\,\frac{\textstyle \partial \rho }{\textstyle \partial t\,} = 0, \end{eqnarray} \tag{ 50 }$$

the expression under the integral becomes

$$\begin{eqnarray} &&\rho ^{\,0} {{{\bf v}}}\cdot \nabla _{ {x\,}} {V_{L}^{\prime }} = \nabla _{ {x\,}}\cdot (\rho ^{\,0}\,{{{\bf v}}}\,{V_{L}^{\prime }})+\,{V_{L}}^{\prime } \frac{\,\partial \rho _{L}^{\prime }}{\partial t\,} , \end{eqnarray} \tag{ 51 }$$

provided we set ∇_x ρ⁰ = 0, that is, provided we assume that the unperturbed body is homogeneous. ⁸ Then the time-averaged power, for an initially homogeneous body, acquires the form of

$$\begin{eqnarray} &&\langle P\rangle \;=\; \int \,\nabla _{ {x\,}}\cdot (\rho ^{\,0}\,{{{\bf v}}}\,{V_{L}}^{\prime })\;d^3x + \int \,{V_{L}}^{\prime } \frac{\,\partial \rho _{L}^{\prime }}{\partial t\,}\;d^3x , \end{eqnarray} \tag{ 52 }$$

where we approximated the Jacobian with unity, thus neglecting higher-order terms.

Similar to Zschau (1978, Equation (2)), we begin with the formula (36) for the power in the Eulerian variables. The next natural step is Equation (40), whereafter integration by parts renders:

$$\begin{eqnarray} P &=& \int \rho _{E}\,{V_{E}}^{\prime }\;{{{\bf v}}}\cdot d{{\bf S}}^t- \frac{1}{4\pi G}\,\int d^3r \nabla _{ {r}}\cdot \left({V_{E}}^{\prime } \,\frac{\partial \nabla _{ {r\,}} {V_{E}}^{\prime }}{\partial t}\,\right) \nonumber\\ &&+ \frac{1}{4\pi G}\,\int d^3r \frac{\partial \nabla _{ {r\,}} {V_{E}}^{\prime }}{\partial t}\,\cdot \,\nabla _{ {r\,}} {V_{E}}^{\prime } \end{eqnarray} \tag{ 53a }$$

$$\begin{eqnarray} &=& \int \rho _{L}\,{V_{E}}^{\prime }\,{{{\bf v}}}\cdot \left({\mathbb {J}} \,d{{\bf S}}\,\right)- \frac{1}{4\pi G}\,\int {V_{E}}^{\prime } \frac{\partial \nabla _{ {r\,}} {V_{E}}^{\prime }}{\partial t}\,\cdot \left({\mathbb {J}}\,d{{\bf S}}\,\right) \nonumber\\ &&+ \frac{1}{8\pi G} \frac{\partial }{\partial t}\int (J\,d^3x) \nabla _{ {r\,}} {V_{E}}^{\prime }\cdot \nabla _{ {r\,}} {V_{E}}^{\prime }. \end{eqnarray} \tag{ 53b }$$

En route from the former expression to the latter, we switch from d S^t and d³ r to ${\mathbb {J}} \,d{{\bf S}}$ and J  d³ x, respectively. Thereby we switch from integration over a deformed body to that over the undeformed one. So ρ_E becomes ρ_L , see Equation (46). A similar switch from V_E ' to V_L ' can be performed using Equation (32), but we prefer to stick to V_E ' for some time, as it will be easier to impose the boundary conditions on the Eulerian potential.

In a leading-order calculation, both the Jacobian and its tensorial analogue may be set unity: ${\mathbb {J}}\approx {\mathbb {I}}$ and J ≈ 1, as evident from the formulae (44) and (45). In the same order, we can substitute ∇_r with ∇_x . In addition, as was explained in Footnote 7, we can substitute ρ_L = ρ⁰/J with ρ⁰, and treat the latter as time-independent. Thus, the time average of the power becomes:

$$\begin{eqnarray} \langle P\rangle &=&\int \, \langle \rho ^{\,0} {V_{E}}^{\prime }\;{{{\bf v}}}\rangle \,\cdot \,d{{\bf S}} - \frac{1}{4\pi G} \int \left\langle {V_{E}}^{\prime } \,\frac{\partial \,\nabla _{ {x}} {V_{E}}^{\prime }}{\partial t}\,\right\rangle \,\cdot \,d{{\bf S}} \nonumber\\ \end{eqnarray} \tag{ 54a }$$

$$\begin{eqnarray} &=&{-} \frac{1}{4\pi G} \int \left\langle {V_{E}}^{\prime } \,\frac{\partial }{\partial t}\,(\nabla _{ {x\,}} {V_{E}}^{\prime } - 4\,\pi \,G\,\rho ^{\,0}\,{{\bf u}})\;\right\rangle \,\cdot \,d{{\bf S}}. \nonumber\\ \end{eqnarray} \tag{ 54b }$$

with the volume integral dropped. ⁹ The potential V_E ' in the above developments was the interior potential, so the above formula should, rigorously speaking, have been written as

$$\begin{eqnarray} \langle P\rangle &=& - \frac{1}{4\pi G} \int \left\langle {{V_{E}}^{\prime }}{}^{{{({\rm interior})}}} \frac{\partial }{\partial t}\,(\nabla _{ {x\,}} {{V_{E}}^{\prime }}{}^{{{({\rm interior})}}} \right. \nonumber\\ &&- \left. 4\,\pi \,G\,\rho ^{\,0}\,{{\bf u}}\,)\;\right\rangle \,\cdot \,d{{\bf S}}. \end{eqnarray} \tag{ 54c }$$

The expression (54c) is somewhat formal. On the one hand, it contains integration over an undeformed surface, an operation appropriate to the Lagrangian description. On the other hand, the quantity under the integral is Eulerian, that is, a function of the perturbed positions. Thus, to employ the expression (54c) in practical calculations, one would first have to express the integrated average product $\langle {{V_{E}}^{\prime }}^{{{\,({\rm interior})}}} ({\textstyle \partial }/{\textstyle \partial t}) \,(\nabla _{ {x\,}} {{V_{E}}^{\prime }}^{{{\,({\rm interior})}}} - 4\,\pi \,G\,\rho ^{\,0}\,{{\bf u}}\,)\;\rangle$ as a function of the unperturbed positions, that is, of the coordinates on the undeformed surface. Simply speaking, one would have to switch from a Eulerian function under the integral to a Lagrangian function, using the formula (32). The reason for our procrastination with this step is the convenience of the Eulerian description for imposing boundary conditions.

Our expression (54c) is equivalent to formula (12) in Zschau (1978). The sole difference is how we justify the substitution of the Lagrangian density ρ_L with the unperturbed ρ⁰. Whereas we approximated the Jacobian with 1 + O(|u|), Zschau (1978, Equation (10)) employed a clever trick that did not rely on the smallness of disturbance. In our notation, the trick looks like this: if in the first term of our expression (54a) we also keep the first-order perturbation $\rho _{L}^{\prime }$ of the density, the time average of the product $\rho _{L}^{\prime }\,{{\bf v}}\,{V_{E}}^{\prime }$  will always be zero, provided all three oscillate at the same frequency. While elegant, Zschau's argument works only for a perturbation at one frequency, not for a spectrum of frequencies.

The treatment by Platzman (1984) contains more inaccuracies. The author's formula (2) looks like our Equation (48b), with the actual density substituted from the beginning by its unperturbed value ρ⁰. Such a start indicates the use of the Lagrangian description. This however comes into contradiction with the way the author writes the conservation law. Platzman's form of that law is equivalent to our Equation (24)—it is written in the Eulerian language. The following Poisson equation is also Eulerian. That the author eventually arrives at the right integral expression (Equation (5) in Platzman (1984)) is more due to luck than to accuracy. In the subsequent derivation, the author's formulae (7) and (10) are incorrect because the fact that the Fourier modes in the Darwin–Kaula theory can be of either sign is neglected. We address this point at the end of Section 5.4.

The Eulerian boundary conditions mimic those from electrostatics (see Appendix E):

$$\begin{eqnarray} &&{{V_{E}}^{\prime }}^{\,\textstyle {^{({\rm interior})}}} = {{V_{E}}^{\prime }}^{\,\textstyle {^{({\rm exterior})}}} \end{eqnarray} \tag{ 55 }$$

and

$$\begin{eqnarray} &&\left[ \frac{\partial }{\partial \hat{{\bf n}}}\, {V_{E}}^{\prime } - 4 \pi G \rho ^{\,0} {\bf u} \right]^{\,\textstyle {^{({\rm exterior})}}}= \left[ \frac{\partial }{\partial \hat{{\bf n}}}\, {V_{E}}^{\prime } - 4 \pi G \rho ^{\,0} {\bf u} \right]^{\,\textstyle {^{({\rm interior})}}}\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!. \end{eqnarray} \tag{ 56 }$$

Insertion thereof into Equation (54c) makes the power look

$$\begin{eqnarray} &&\langle P\rangle = - \frac{1}{4\pi G} \int \left\langle {{V_{E}}^{\prime }}^{{{\,({\rm exterior})}}} \frac{\partial }{\partial t} \nabla _{ {x\,}} {{V_{E}}^{\prime }}^{{{\,({\rm exterior})}}}\;\right\rangle \,\cdot \,d{{\bf S}}. \qquad \end{eqnarray} \tag{ 57 }$$

It is now time to write the expression under the integral (57) as a function of the coordinates on the unperturbed surface, the one over which we integrate. The formula (32) prescribes us to substitute V_E ' with V_L ' − u · ∇_x V₀. As u is zero outside the body, we get ¹⁰

$$\begin{eqnarray} &&\langle P\rangle = - \frac{1}{4\pi G} \int \left\langle {{V_{L}}^{\prime }}^{{{\,({\rm exterior})}}} \frac{\partial }{\partial t} \nabla _{ {x\,}} {{V_{L}}^{\prime }}^{{{\,({\rm exterior})}}}\;\right\rangle \,\cdot \,d{{\bf S}}. \end{eqnarray} \tag{ 58 }$$

To analyze the behavior of V' outside the perturbed body, recall that its two components, U and W, scale differently with the planetocentric radius. As can be seen from Equation (1), the degree-l Legendre component of the perturbing potential changes as ¹¹ W_l  ∝ r^l . According to Equation (11), the degree-l component of the tidal potential obeys U_l  ∝ r^{−(l + 1)}. All in all, the l-degree part of the exterior V' assumes the form of

$$\begin{eqnarray} &&{{V_{L}}^{\prime }}^{\,\textstyle {^{({\rm exterior})}}}= \sum _{l=2}^{\infty }\,\left[\left(\frac{r}{R}\,\right)^{{{l}}}\,W_l(R)+\,\left(\frac{r}{R}\, \right)^{{{-(l+1)}}}\,U_l(R)\,\right], \nonumber\\ \end{eqnarray} \tag{ 59 }$$

while the normal part of its gradient on the free surface is

$$\begin{eqnarray} &&\frac{\partial \,}{\partial r} {{V_{L}}^{\prime }}_{ {l}}^{\,\textstyle {^{({\rm exterior})}}} = R^{-1}\,\sum _{l=2}^{\infty }\,\left[l W_l - (l+1) U_l\,\right]. \end{eqnarray} \tag{ 60 }$$

Plugging it into Equation (58), and benefitting from the orthogonality of surface harmonics, we obtain: ¹²

$$\begin{eqnarray} \langle P\rangle &=&{-} \frac{1}{4\pi G R} \sum _{l=2}^{\infty }\,\int \left\langle l U_l\,\dot{W}_l - (l+1) W_l\,\dot{U}_l \right\rangle dS \qquad \end{eqnarray} \tag{ 61a }$$

$$\begin{eqnarray} &=&\frac{1}{4\pi G R} \sum _{l=2}^{\infty }\,(2\,l+\,1)\,\int \langle W_l \dot{U}_l \rangle dS , \end{eqnarray} \tag{ 61b }$$

which is equivalent to the formulae (18) in Zschau (1978) and (5) in Platzman (1984). This, however, is the last point on which we are still in agreement with our predecessors.

Bringing in the dynamical Love numbers k_l and the phase lags defined in Equation (12), one can express the products $W_l(t)\,\dot{U}_l(t)$ via the spectral components of the disturbance W(t). ¹³

Although the formula

$$\begin{eqnarray} && \sum _{l=2}^{\infty }(2l+1)\langle W_l(t)\,\dot{U}_l(t)\rangle \nonumber\\ &&\quad = \sum _{\omega }\,(2l+1)\,\frac{\omega }{2}\,W^{\,2}_l(\omega)\,k_l(\omega) \sin \epsilon _l(\omega) \end{eqnarray} \tag{ 62 }$$

is often used in the literature (Zschau 1978; Platzman 1984; Segatz et al. 1988), ¹⁴ accurate examination demonstrates that it is incorrect. To appreciate this, one simply has to insert the expansions (2) and (12) into the formula (61b) and see what happens.

That the answer differs from Equation (62) was noticed by Peale & Cassen (1978). However, their development also needs correction. Later, we examine this matter in great detail and provide a full inventory of the terms emerging in the spectral expansion for damping rate. At this point, we only mention the two key circumstances:

• 1.  
  The conventional expression (62) ignores the degeneracy of modes, that is, a situation where several modes ω_lmpq with different sets lmpq take the same numerical value ω. As will be demonstrated in Section 6, the sum over modes ω in Equation (62) should be substituted with a sum over distinct values of the modes:

  $$\begin{eqnarray*} &&{\rm instead\; of} \sum _\omega W_l^2(\omega)\,\, \hbox{in Equation (62) use this} : \nonumber\\ &&\sum _\omega \left[\sum _{{\omega _{ {{lmpq}}}= \omega }} W_{l}(\omega _{lmpq})\,\right]^2 \end{eqnarray*}$$

  $$\begin{eqnarray*} && {\rm where} \textstyle\sum _{{\omega _{ {{lmpq}}}= \omega }} W_{l}(\omega _{lmpq})\,\, {\rm denotes\; a\; sum\; of\; all\; terms\; for\; which } \nonumber\\ && \omega _{lmpq}\, {\rm takes\; a\; value }\,\omega. \end{eqnarray*}$$

  In short: first sum all the terms corresponding to one value of ω, then square the sum, and only then sum over all the values of ω.
• 2.  
  Much less intuitive is the fact that the spectral sum will contain extra terms that are missing completely in the expression (62). As we shall see in Appendix F, these terms look (up to some caveat) as W_l (ω) W_l (− ω). They show up because two modes of opposite values, ω and −ω, correspond to the same physical frequency |ω|.

For the time being, we use the notation $\sum ^{\,\textstyle \sharp }$:

$$\begin{eqnarray} && \sum _{l=2}^{\infty }(2l+1)\,\langle W_l(t)\,\dot{U}_l(t) \rangle \nonumber\\ &&\quad = {\sum _{\omega }}^{\textstyle {\sharp } }(2l+1)\,\frac{\omega }{2}\,W^{\,2}_l(\omega)\,k_l(\omega) \sin \epsilon _l(\omega), \end{eqnarray} \tag{ 63 }$$

where the superscript $^{\textstyle {\sharp } }$ reminds the reader that the spectral sum needs to be amended down the road.

Insertion of Equation (63) into Equation (61b) results in: ¹⁵

$$\begin{eqnarray} &&\langle P \rangle = \frac{1}{8\pi G R} {\sum _{\omega }}^{\textstyle {\sharp } }\,(2\,l+\,1) \omega \,k_l(\omega) \sin \epsilon _l(\omega)\int W^{{\,2}}_l(\omega) dS . \nonumber\\ \end{eqnarray} \tag{ 64 }$$

If not for the superscript $^{\textstyle {\sharp } }$, this expression would coincide with the results by Zschau (1978) and Platzman (1984). ¹⁶ The superscript reminds us of the important caveat in the evaluation of the sum: the factors $W^{\,2}_l(\omega)$ should be substituted with more complicated expressions, whereas the sum should be carried not over all modes ω = ω_lmpq , but over all distinct values of ω, see Section 6.

Consider a low-inclined perturber. From the tides it creates, the perturber gets predominantly transversal orbital disturbance. We need two planetary equations in the Gauss form (Brouwer & Clemence 1961, p. 301, Equation (33)): ²²

$$\begin{eqnarray*} \dot{\omega } &=& \frac{\sqrt{1 - e^2}}{n a e} \left[- R \cos f + \left(1+\,\frac{r}{p}\,\right) T \sin f\,\right] \nonumber\\ &&- \frac{\,\sin (\omega +\,f) \,\cot i\,}{n a \sqrt{1 - e^2}} \,\frac{\,r\,}{a} W, \end{eqnarray*}$$

$$\begin{eqnarray*} &&\dot{{\cal {M}}_0} = \frac{1 - e^2}{n a e} \left[\left(\cos f - 2 \,\frac{\,r\,}{p}\, e\,\right)\,R - \left(1+\,\frac{r}{p}\,\right) T \sin f\,\right], \nonumber \end{eqnarray*}$$

where f is the true anomaly, p ≡ a (1 − e²) is the semilatus rectum, while R, T, and W are the radial, transversal, and normal-to-orbit forces, respectively. In a situation where the perturbation is predominantly transversal and the terms with R and W may be neglected, we obtain:

$$\begin{eqnarray} &&\dot{\omega }+ \dot{{\cal {M}}_0} \approx \frac{e}{1+\,\sqrt{1 - e^2}} \frac{\sqrt{1 - e^2}}{n a} \left(1+\,\frac{r}{p}\,\right) T \sin f. \qquad \end{eqnarray} \tag{ B4 }$$

A low-inclined moon gets predominantly transversal orbital disturbance from the tides it creates in the planet. Inserting the latter expression in the formula (B2) for the Fourier modes, we see that for the modes with a zero q (like the semidiurnal tide parameterized with lmpq = 2200) the input from the pericenter rate may be omitted if the eccentricity is not too large. Indeed, for q = 0, the term $q\dot{\cal {M}}_0$ vanishes, while the term $(l-2p)(\dot{\omega }+\dot{\cal {M}}_0)$ is now approximated with (l − 2p) multiplied by the expression (B4). Although $\dot{\omega }$ and $\dot{\cal {M}}_0$ can, separately, be substantial, their sum (B4) is smaller by the order of e. Being (in this particular case) of the same order but of opposite sign, $\dot{\omega }$ and $\dot{\cal {M}}_0$ largely compensate one another. Therefore, if we choose to drop $\dot{\cal {M}}_0$, we should also drop $\dot{\omega }$. In this special case, dropping both will be legitimate.

As a useful aside, we would remind that the mean longitude is defined through $L\equiv {\cal {M}}+\,\omega +\,\Omega$, its rate being $\dot{L}\equiv \sqrt{G\,(M+\,M^*) a^{-3}(t)\,}+\,\dot{\cal {M}}_0+\,\dot{\omega }+\,\dot{\Omega }$. As we have just seen, the rates $\dot{\cal {M}}_0$ and $\dot{\omega }$ largely compensate one another and may both be neglected in the considered case. If, above that, the rate of the node happens to be negligible, then the mean motion from the Kepler law will be close to the mean longitude rate.

The situation is different where the principal perturbation is due to the oblateness of the tidally perturbed primary. The mean rates (Vallado 2007, pp. 647–648)

$$\begin{eqnarray*} \dot{{\cal {M}}_0} &=& \frac{\,3\,}{\,4\,} n J_2 \frac{R^2}{a^2} \frac{\,2 - 3 \sin ^2i\,}{(1 - e^2)^{3/2}} \nonumber\\ &&{\rm and}\nonumber\\ \dot{\omega } &=& \frac{\,3\,}{\,4\,} n J_2 \frac{R^2}{a^2} \frac{\,4 - 5 \sin ^2i\,}{(1 - e^2)^{2}}. \end{eqnarray*}$$

are of the same order and sign. Therefore, when keeping $\dot{\omega }$, we must also include $\dot{{\cal {M}}_0}$.

The easiest way to get rid of $\dot{{\cal {M}}_0}$ is to define the mean motion as $n\equiv {{\bf \dot{\cal {M}}}}$. This, the so-called anomalistic mean motion will, however, differ from $\,\sqrt{G\,(M+\,M^*) a^{-3}(t)\,}$.

Under perturbation, two changes will happen in a point ${{\boldsymbol {r}}}$ at a time t:

• 1.  
  Physical fields will now acquire different values in this point at this time.For example, a moon fixed in a certain position relative to the planetary surface will render some distribution of its potential over the volume of the host planet. The same moon fixed in a different position will generate a different distribution of its potential. This, in its turn, will yield a different deformation of the planet and therefore a different spatial distribution of its tidal-response potential and of all other quantities.Thus, instead of the unperturbed Eulerian dependencies $q^0_{E}({{\boldsymbol {r}}}, t)$, we now have

  $$\begin{eqnarray} &&q_{E}({{\boldsymbol {r}}}, t) = q^0_{E}({{\boldsymbol {r}}}, t) + q^{\prime }_{E}({{\boldsymbol {r}}}, t). \end{eqnarray} \tag{ D1 }$$
• 2.  
  A different particle will now arrive in the point ${{\boldsymbol {r}}}$ at the time t. It will not be the same particle as the one expected there at the time t in the absence of perturbation.On the other hand, a particle that starts in ${{\boldsymbol {X}}}$ at t = 0, will appear, at the time t, not in the point ${{\boldsymbol {X}}}={\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)$ but in some other place

  $$\begin{eqnarray} &&{{\boldsymbol {r}}}= {{\boldsymbol {X}}}+ {{\bf u}}= {\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t) + {{\bf u}}({{\boldsymbol {X}}},t). \end{eqnarray} \tag{ D2 }$$

These two changes will influence the Lagrangian dependencies on the initial conditions. The dependency of each field will obtain a variation $q^{\prime }_{L}({{\boldsymbol {X}}},t)$:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}},t) = q^0_{L}({{\boldsymbol {X}}},t) + q^{\prime }_{L}({{\boldsymbol {X}}},t). \end{eqnarray} \tag{ D3 }$$

In the absence of perturbation, the particle ${{\boldsymbol {X}}}$ was destined to arrive in ${{\boldsymbol {X}}}$, wherefore $q^0_{L}({{\boldsymbol {X}}},t)$ was defined through Equation (17). Under perturbation, the same particle ${{\boldsymbol {X}}}$ is expected to end up in ${{\boldsymbol {r}}}$, so the Lagrangian dependency becomes

$$\begin{eqnarray} q_{L}({{\boldsymbol {X}}},t)&\equiv &q_{E}({{\boldsymbol {r}}}, t) \end{eqnarray} \tag{ D4a }$$

$$\begin{eqnarray} &&\nonumber \\ &=&q_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)+\,{{\bf u}}({{\boldsymbol {X}}},t), t) \end{eqnarray} \tag{ D4b }$$

$$\begin{eqnarray} &&\nonumber \\ &=&q_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t) + {{\bf u}}({{\boldsymbol {X}}},t) \nabla _{ {x\,}} q_{E} + O({{\bf u}}^2) \end{eqnarray} \tag{ D4c }$$

$$\begin{eqnarray} &&\nonumber \\ &=&q^0_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t) + q^{\prime }_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t) + {{\bf u}}({{\boldsymbol {X}}},t) \nabla _{ {x\,}} q_{E} + O({{\bf u}}^2). \nonumber\\ \end{eqnarray} \tag{ D4d }$$

Insertion of Equation (D3) into the left-hand side of Equation (D4) will give us:

$$\begin{eqnarray*} q^0_{L}({{\boldsymbol {X}}},t) + q^{\prime }_{L}({{\boldsymbol {X}}},t) &=& q^0_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t) + q^{\prime }_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t) \nonumber\\ && + {{\bf u}}({{\boldsymbol {X}}},t) \nabla _{ {x\,}} q_{E} + O({{\bf u}}^2). \end{eqnarray*}$$

Subtracting Equation (17b) from this formula, we arrive at a relation between the perturbations of the Lagrangian and Eulerian quantities:

$$\begin{eqnarray} &&q^{\prime }_{L}({{\boldsymbol {X}}},t) = q^{\prime }_{E}({\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t), t) + {{\bf u}}({{\boldsymbol {X}}},t) \nabla _{ {x\,}} q_{E} + O({{\bf u}}^2), \end{eqnarray} \tag{ D5 }$$

where the first term on the right-hand side, $q^{\prime }_{E}({{\boldsymbol {X}}},t)$, expresses the change in the final spatial distribution of the field q. The other two terms show up because perturbation changes the mapping from ${{\boldsymbol {X}}}$ to the current location.

A slightly different, although equally valid viewpoint is possible. In a reference setting at time t, an observer located in ${{\boldsymbol {X}}}$ will see the arrival of a particle that started from ${{\boldsymbol {X}}}$:

$$\begin{eqnarray} &&q^0_{L}({{\boldsymbol {X}}},t) \equiv q^0_{E}({{\boldsymbol {X}}}, t). \end{eqnarray} \tag{ D6 }$$

In a perturbed situation, the same observer in ${{\boldsymbol {X}}}$ will register, at the time t, the arrival of a different particle, one that started from ${{\boldsymbol {X}}}-{{\bf U}}$:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}}-{{\bf U}},t) \equiv q_{E}({{\boldsymbol {X}}}, t), \end{eqnarray} \tag{ D7a }$$

which is:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}},t) - {{\bf U}}\nabla _{X}q_{L} + O(U^2) = q_{E}({{\boldsymbol {X}}}, t). \end{eqnarray} \tag{ D7b }$$

Subtraction of Equation (D6) from Equation (D7b) gives us the variations:

$$\begin{eqnarray} &&q_{L}({{\boldsymbol {X}}},t) - q^0_{L}({{\boldsymbol {X}}},t) - {{\bf U}}\nabla _{X}q_{L} + O(U^2) = q_{E}({{\boldsymbol {X}}}, t) - q^0_{E}({{\boldsymbol {X}}}, t) \nonumber\\ \end{eqnarray} \tag{ D8a }$$

or, simply,

$$\begin{eqnarray} &&q^{\prime }_{L}({{\boldsymbol {X}}},t) = q^{\prime }_{E}({{\boldsymbol {X}}}, t) + {{\bf U}}\nabla _{X}q_{L} + O(U^2). \end{eqnarray} \tag{ D8b }$$

Introducing the Jacobian $\;J\equiv \,{\textstyle dV^t}/{\textstyle dV^0} = \;{\rm det}\,({\textstyle \partial x_i}/{\textstyle \partial X_j})$, we write:

$$\begin{eqnarray} &&{{\bf U}}\nabla _{X}q_{L} = {{\bf U}}J \nabla _{ {x\,}} q_{E} = {{\bf u}}\nabla _{ {x\,}} q_{E}, \end{eqnarray} \tag{ D9 }$$

with u ≡ U J. Thus Equation (D8b) and Equation (D5) are equivalent insofar as O(U²) = O(u²).

While the language of Equation (D5) is more conventional than that of Equation (D8), the latter description is easier for physical interpretation. Suppose we are observing the gradual cooling of a flow. In an unperturbed setting, a particle that started in ${{\boldsymbol {X}}}$ at the time t = 0, will show up in ${{\boldsymbol {X}}}$ at the time t. Accordingly, a measurement of the temperature in ${{\boldsymbol {X}}}$ at the time t will render, in the absence of perturbation, a value to which the particle ${{\boldsymbol {X}}}$ has cooled down by this time—see the equality (Equation (D6)).

Under perturbation, the rate of cooling of each particle will change. In addition, owing to the change of trajectories, a different particle will show up in ${{\boldsymbol {X}}}$ at the time t. Now this will be a particle that started its movement at t = 0 from some point ${{\boldsymbol {X}}}-{{\bf U}}$. So a measurement of the temperature in the point ${{\boldsymbol {X}}}$ at the time t will now give us a temperature value to which the particle ${{\boldsymbol {X}}}-{{\bf U}}$ has cooled down—see Equation (D7a).

The difference between the measurements performed in ${{\boldsymbol {X}}}$ in the perturbed and unperturbed cases will, according to Equation (D8b), read as $q^{\prime }_{E}({{\boldsymbol {X}}}, t) = q^{\prime }_{L}({{\boldsymbol {X}}},t) - {{\bf U}}\nabla _{X}q_{L} + O(U^2)$. The first term on the right renders the cooling down of the particle arriving in ${{\boldsymbol {X}}}$ at the time t, while the second and third terms reflect the fact that, under disturbance, we register a particle arriving from a point displaced by U, compared to the particle that would be brought to ${{\boldsymbol {X}}}$ by an unperturbed flow.

With aid of Equation (D9), the expression (D8b) can be equivalently rewritten as

$$\begin{eqnarray} &&D_t = \partial _t + {{\bf v}}\,\nabla _{ {x}}, \end{eqnarray} \tag{ D10 }$$

where v ≡ ∂u/∂t, while D ≡ d/dt is the comoving derivative. The physical interpretation of Equation (D10) is obvious: the rate of cooling of a moving particle, dq/dt, can be measured by a quiescent observer. The observer, however, must amend his result, ∂q/∂t, with a correction taking into account the fact that, being quiescent, he is measuring the difference between the temperature of different particles passing by, not of the same particle.

Hereafter, we shall restrict our consideration to the case of a periodically deformed solid. It is natural to associate the reference trajectory x = f  ( X , t) with the equilibrium configuration. In this configuration, the particles stay idle, so ${{\boldsymbol {X}}}$ coincides with the initial value ${{\boldsymbol {X}}}$:

$$\begin{eqnarray} &&{{\boldsymbol {X}}}= {\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t),\quad {\rm where}\; {\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t) = {{\boldsymbol {X}}}\quad {\rm for\; all }\; t, \quad \end{eqnarray} \tag{ D11 }$$

while all properties keep in time their fiducial values:

$$\begin{eqnarray} &&q^0_{L}({{\boldsymbol {X}}},t) = q^0_{E}({{\boldsymbol {X}}},t) = q^0. \end{eqnarray} \tag{ D12 }$$

Be mindful that the superscript 0 did not originally mark a value fixed in time, but a trajectory chosen to be reference. It is only now that the role of a reference configuration is played by an equilibrium body, that the superscript 0 begins to denote an unchanging value.

For a particle originally located in ${{\boldsymbol {X}}}$, its perturbed trajectory ${{\boldsymbol {r}}}$ differs from its reference trajectory ${{\boldsymbol {X}}}$ by some u:

$$\begin{eqnarray} &&{{\boldsymbol {r}}}= {{\boldsymbol {X}}}+ {{\bf u}}= {\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t) + {{\bf u}}({{\boldsymbol {X}}},t). \end{eqnarray} \tag{ D13 }$$

When the reference trajectory is the equilibrium, insertion of Equation (D11) into Equation (D13) results in

$$\begin{eqnarray} &&{{\boldsymbol {r}}}= {{\boldsymbol {X}}}+ {{\bf u}}({{\boldsymbol {X}}},t), \end{eqnarray} \tag{ D14a }$$

which can also be written as

$$\begin{eqnarray} &&{{\boldsymbol {r}}}= {{\boldsymbol {X}}}+ {{\bf u}}({{\boldsymbol {X}}},t), \end{eqnarray} \tag{ D14b }$$

because, in this case, the unperturbed trajectory ${{\boldsymbol {X}}}$ always coincides with the initial value ${{\boldsymbol {X}}}$.

We work in a linearized approximation, neglecting the term O(u²) in Equation (D5) and writing all expansions up to terms linear in the displacement u or velocity v ≡ d u/dt.

For a short and simple explanation of the linearized Lagrangian and Eulerian descriptions of tides, see Wang (1997). A more comprehensive treatment is offered in the book by Dahlen & Tromp (1998, Section 3.1.1).

Denote the Eulerian value of the mass density with $\rho _{E}({{\boldsymbol {r}}},t)$. As mass cannot be destroyed or created, so its amount in a comoving volume V^t of a flow stays constant:

$$\begin{eqnarray} &&\frac{d}{dt}\,\int _{V^t}\rho _{E}\,dV^t = 0. \end{eqnarray} \tag{ D15 }$$

For the reference history ${{\boldsymbol {X}}}={\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)$, this would imply:

$$\begin{eqnarray} &&\int _{V^{t}_{{\rm ref}}}\rho ^{\,0}_{E}({{\boldsymbol {X}}}, t)\,d^4{{\boldsymbol {X}}}= \int _{V^0}\rho ^{\,0}_{E}({{\boldsymbol {X}}}, 0)\,d^4{{\boldsymbol {X}}}. \end{eqnarray} \tag{ D16a }$$

For a perturbed history ${{\boldsymbol {r}}}={\boldsymbol {{\boldsymbol {f}}}}({{\boldsymbol {X}}},t)+{{\bf u}}({{\boldsymbol {X}}},t)$, we have:

$$\begin{eqnarray} &&\int _{V^{t}_{{\rm pert}}}\rho _{E}({{\boldsymbol {r}}}, t)\,d^4{{\boldsymbol {r}}}= \int _{V^0}\rho _{E}({{\boldsymbol {X}}}, 0)\,d^4{{\boldsymbol {X}}}. \end{eqnarray} \tag{ D16b }$$

For each individual particle, its perturbed trajectory ${{\boldsymbol {r}}}$ stems from the same initial position ${{\boldsymbol {X}}}$ as the appropriate reference trajectory ${{\boldsymbol {X}}}$, so the initial densities are the same:

$$\begin{eqnarray} &&\rho _{E}({{\boldsymbol {X}}}, 0)\;=\;\rho ^{\,0}_{E}({{\boldsymbol {X}}}, 0). \end{eqnarray} \tag{ D17 }$$

At later times, however, $\rho _{E}({{\boldsymbol {r}}},t)$ and $\rho ^{\,0}_{E}({{\boldsymbol {X}}},t)$ have different functional forms.

D.4.1. The Continuity Law in the Eulerian Description

The right-hand sides of the formulae (D16a) and (D16b) coincide, as they render the mass of the same initial distribution $\rho _{E}({{\boldsymbol {X}}},0)=\rho ^{\,0}_{E}({{\boldsymbol {X}}},0)$. Thus the left-hand sides of the two formulae also coincide:

$$\begin{eqnarray} &&0 = \int _{V^{t}_{{\rm pert}}}\rho _{E}({{\boldsymbol {r}}}, t)\,d^4{{\boldsymbol {r}}}- \int _{V^t_{{\rm ref}}}\rho ^{\,0}_{E}({{\boldsymbol {X}}}, t)\,d^4{{\boldsymbol {X}}}, \end{eqnarray} \tag{ D18 }$$

as the mass stays unchanged, no matter whether the system follows the reference history or a perturbed one. Now switch from the perturbed coordinates, ${{\boldsymbol {r}}}$, to the reference ones, ${{\boldsymbol {X}}}$:

$$\begin{eqnarray} 0&=&\int \,d^4{{\boldsymbol {X}}} [J \rho _{E}({{\boldsymbol {r}}}, t) - \,\rho ^{\,0}_{E}({{\boldsymbol {X}}}, t)] \end{eqnarray} \tag{ D19a }$$

$$\begin{eqnarray} &&\nonumber \\ &=&\int \,d^4{{\boldsymbol {X}}} [\left(1+\,\nabla _{ {{{x}}}}\cdot {{\bf u}}\,\right) \rho _{E}({{\boldsymbol {r}}}, t) - \,\rho ^{\,0}_{E}({{\boldsymbol {X}}}, t)], \qquad \end{eqnarray} \tag{ D19b }$$

where the Jacobian is:

$$\begin{eqnarray} J \equiv \frac{dV^{t}_{{\rm pert}}}{dV^{t}_{{\rm ref}}} &=& {\rm det}\; \frac{\partial r_i}{\partial x_j} = 1 + \nabla _{ {{{x}}}}\cdot {{\bf u}}+ O({{\bf u}}^2) \nonumber\\ &=& 1 + J \nabla _{ {r}}\cdot {{\bf u}}+ O({{\bf u}}^2). \end{eqnarray} \tag{ D20a }$$

From this, we see that the Jacobian can also be written as ²³

$$\begin{eqnarray} &&J = \frac{\textstyle 1 + O({{\bf u}}^2)}{\textstyle 1 - \nabla _{ {r}}\cdot {{\bf u}}} = 1 + \nabla _{ {r}}\cdot {{\bf u}}+ O({{\bf u}}^2). \end{eqnarray} \tag{ D20b }$$

From Equation (D19a), we obtain the exact equality

$$\begin{eqnarray} &&\rho _{E}({{\boldsymbol {r}}}, t) J = \rho ^{\,0}_{E}({{\boldsymbol {X}}}, t), \end{eqnarray} \tag{ D21 }$$

a linearized version thereof being

$$\begin{eqnarray} 0&=&\rho _{E}({{\boldsymbol {r}}}, t) (1+\,\nabla _{ {r}}\cdot {{\bf u}}) - \rho ^{\,0}_{E}({{\boldsymbol {r}}}-{{\bf u}}, t) + O({{\bf u}}^2) \qquad \end{eqnarray} \tag{ D22a }$$

$$\begin{eqnarray} &=& \rho _{E}({{\boldsymbol {r}}}, t) - \rho ^{\,0}_{E}({{\boldsymbol {r}}}, t) + \rho _{E}({{\boldsymbol {r}}}, t) \nabla _{ {r}}\cdot {{\bf u}}+ {{\bf u}}\nabla _{ {r\,}}\rho ^{\,0}_{E}({{\boldsymbol {r}}}, t) + O({{\bf u}}^2). \nonumber\\ \end{eqnarray} \tag{ D22b }$$

For a small t, the deviation u between the two trajectories is linear in time, and so is the difference between the perturbed and reference density functions. Thus, we may change $\rho _{E}({{\boldsymbol {r}}},t)\,\nabla _{ {r}}\cdot {{\bf u}}$ to $\rho ^{\,0}_{E}({{\boldsymbol {r}}},t)\,\nabla _{ {r}}\cdot {{\bf u}}$, to obtain an expression correct to first order in u:

$$\begin{eqnarray} &&\rho _{E}^{\prime } + \nabla _{ {r}}\cdot \left(\rho ^{\,0}_{E}\,{{\bf u}}\right) = 0, \end{eqnarray} \tag{ D23 }$$

where the finite variation is

$$\begin{eqnarray} &&{\rho _{E}}^{\prime }({{\boldsymbol {r}}},t) \equiv \rho _{E}({{\boldsymbol {r}}},t) - \rho ^{\,0}_{E}({{\boldsymbol {r}}},t). \end{eqnarray} \tag{ D24 }$$

We would reiterate that the perturbative approach to Eulerian quantities implies a comparison between their present spatial distributions. So the two histories are compared in the same point ${{\boldsymbol {r}}}$ and at the same time t.

In Equation (D22b), we could also have changed ${{\bf u}}\nabla _{ {r\,}}\rho ^{\,0}_{E}({{\boldsymbol {r}}},t)$ to ${{\bf u}}\nabla _{ {r\,}}\rho _{E}({{\boldsymbol {r}}},t)$. Then, instead of Equation (D23), we would have obtained ρ_E ' + ∇_r · (ρ_E   u) = 0, without the superscript 0 in the second term. In the linear approximation, however, this would be no better than Equation (D23). Traditionally, the form (D23) is preferred in the literature.

However, when switching to a differential form of the conservation law, we no longer need to keep the superscript 0, because the difference between $\rho _{E}({{\boldsymbol {r}}},t)$ and $\rho ^{\,0}_{E}({{\boldsymbol {r}}},t)$ becomes infinitesimally small. So the differential law reads as:

$$\begin{eqnarray} &&\frac{\partial \rho _{E}}{\partial t} + \nabla _{ {r}}\cdot \left(\rho _{E}\,{{\bf v}}\right) = 0, \end{eqnarray} \tag{ D25 }$$

where v ≡ ∂u/∂t, the partial derivative giving the rate of change with coordinates fixed.

Employment of the perturbative formula (D23) near a deformable free boundary requires some care. On the one hand, the reference density $\rho ^{\,0}_{E}({{\boldsymbol {r}}})$ makes an abrupt step there. On the other hand, due to deformation of the boundary, we may get a finite present density in a point where the reference density used to be zero, and vice versa.

D.4.2. The Continuity Law in the Lagrangian Description

The Lagrangian density is introduced in the standard way (Equation (D4a)):

$$\begin{eqnarray} &&\rho _{L}({{\boldsymbol {X}}},t)\equiv \rho _{E}({{\boldsymbol {r}}},t), \end{eqnarray} \tag{ D26 }$$

so the formula (D21) becomes:

$$\begin{eqnarray} &&\rho _{L}({{\boldsymbol {X}}}, t) J = \rho ^{\,0}_{E}({{\boldsymbol {X}}}, t). \end{eqnarray} \tag{ D27a }$$

For a periodically deformed solid, the reference density $\rho ^{\,0}_{E}({{\boldsymbol {X}}},t)$ is the density of the undeformed, stable configuration. So $\rho ^{\,0}_{E}({{\boldsymbol {X}}},t) = \rho ^{\,0}_{E}({{\boldsymbol {X}}},0) = \rho ^{\,0}({{\boldsymbol {X}}})$ is time-independent, and the equality (Equation (D27a)) becomes simply

$$\begin{eqnarray} &&\rho _{L}({{\boldsymbol {X}}}, t) J = \rho ^{\,0}_{E}({{\boldsymbol {X}}}). \end{eqnarray} \tag{ D27b }$$

In accordance with the general formula (D5), we interrelate the density variations as

$$\begin{eqnarray} &&\rho _{L}^{\prime }({{\boldsymbol {X}}}, t) = \rho _{E}^{\prime }({{\boldsymbol {X}}}, t) + {{\bf u}}\nabla _{ {x}} \rho ^{\,0}_{E}({{\boldsymbol {X}}}, t). \end{eqnarray} \tag{ D28a }$$

In the considered case of small periodic variations, the reference trajectory is simply ${{\boldsymbol {X}}}={{\boldsymbol {X}}}$ at all times; so on the right-hand side of the above formula we have a gradient of a constant-in-time stationary distribution: $\,\nabla _{ {x\,}} \rho ^{\,0}({{\boldsymbol {X}}},t) = \nabla _{ {X\,}} \rho ^{\,0}({{\boldsymbol {X}}})$. Thence we obtain:

$$\begin{eqnarray} &&\rho _{L}^{\prime }({{\boldsymbol {X}}}, t) = {\rho _{E}}^{\prime }({{\boldsymbol {X}}}, t) + {{\bf u}}\nabla _{ {X\,}} \rho ^{\,0}_{E}({{\boldsymbol {X}}}). \end{eqnarray} \tag{ D28b }$$

Combining this formula with Equation (D23), we arrive at ²⁴

$$\begin{eqnarray} &&\rho _{L}^{\prime } + \rho ^{\,0}_{E} \nabla _{ {X}}\cdot {{\bf u}}= 0, \end{eqnarray} \tag{ D29 }$$

where $\rho _{L}^{\prime }=\rho _{L}^{\prime }({{\boldsymbol {X}}},t)$, while $\rho ^{\,0}_{E}=\rho ^{\,0}_{E}({{\boldsymbol {X}}})$.

D.5.1. In the Eulerian Description

Perturbed or not, the density always obeys the Poisson equation, while the perturbing potential W obeys the Laplace equation outside the perturber:

$$\begin{eqnarray} \nabla _{ {r\,}}^{\,2}\,V_{E}&=&{-} 4\,\pi \,G\,\rho _{E}, \end{eqnarray} \tag{ D30a }$$

$$\begin{eqnarray} &&\nonumber \\ \nabla _{ {r\,}}^{\,2}\,V^{\,0}_{E}&=&{-} 4\,\pi \,G\,\rho ^{\,0}_{E}, \end{eqnarray} \tag{ D30b }$$

$$\begin{eqnarray} &&\nonumber \\ \nabla _{ {r\,}}^{\,2}\,W_{E}&=&0. \end{eqnarray} \tag{ D30c }$$

Subtraction of Equation (D30b) from Equation (D30a) results in a Poisson equation for the density perturbation:

$$\begin{eqnarray} &&\left. \quad \right. \nabla _{ {r\,}}^{\,2}\,{V_{E}}^{\prime } = - 4 \pi G {\rho _{E}}^{\prime } \end{eqnarray} \tag{ D31 }$$

or, equivalently:

$$\begin{eqnarray} &&\left. \quad \right. \nabla _{ {r\,}}^{\,2}\,{U_{E}} = - 4 \pi G {\rho _{E}}^{\prime }, \end{eqnarray} \tag{ D32 }$$

where we took into account the relations (31) and (D30c).

D.5.2. In the Lagrangian Description

Insertion of the formulae (27) and (32) into the Eulerian version of the Poisson equation, Equation (D31), results in the Lagrangian version of this equation:

$$\begin{eqnarray} &&\nabla _{ {r\,}}^{\,2}\, ({V_{L}}^{\prime }- {{\bf u}}\cdot \nabla _{ {x\,}} V^{\,0}) = - 4 \pi G (\rho _{L}^{\prime }- {{\bf u}}\cdot \nabla _{ {x\,}} \rho ^{\,0}). \qquad \end{eqnarray} \tag{ D33a }$$

A switch to differentiation over the initial position, ∇_x , would entail corrections of the order of O(u²). In neglect of those, the equation may be written as

$$\begin{eqnarray} &&\nabla _{ {{x\,}}}^{\,2} ({V_{L}}^{\prime }- {{\bf u}}\cdot \nabla _{ {{x\,}}} V^{\,0}) = - 4 \pi G (\rho _{L}^{\prime }- {{\bf u}}\cdot \nabla _{ {{x\,}}} \rho ^{\,0}). \qquad \end{eqnarray} \tag{ D33b }$$

For an initially homogeneous body, the above formulae simplify to:

$$\begin{eqnarray} &&\nabla _{ {r\,}}^{\,2} ({V_{L}}^{\prime }- {{\bf u}}\cdot \nabla _{ {x\,}} V^{\,0}) = - 4 \pi G \rho _{L}^{\prime } \end{eqnarray} \tag{ D34a }$$

and

$$\begin{eqnarray} &&\nabla _{ {x\,}}^{\,2} ({V_{L}}^{\prime }- {{\bf u}}\cdot \nabla _{ {{x}\,}} V^{\,0}) = - 4 \pi G \rho _{L}^{\prime }. \end{eqnarray} \tag{ D34b }$$