Tides between the TRAPPIST-1 planets

TL;DR

This study compares tidal heating from star-planet and planet-planet interactions in the TRAPPIST-1 system, finding stellar tides dominate except under very low viscosities and specific orbital conditions.

Contribution

It provides the first quantitative assessment of the relative importance of planet-planet versus stellar tides in the TRAPPIST-1 system.

Findings

Stellar tides dominate over planet-planet tides in energy dissipation.

Planet-planet tides can contribute up to 20% of total tidal heating for TRAPPIST-1 g.

Orbital eccentricities must be extremely low ($10^{-3}$ to $10^{-4}$) for planet-planet tides to surpass stellar tides.

Abstract

The TRAPPIST-1 system is sufficiently closely packed that tides raised by one planet on another are significant. We investigate whether this source of tidal heating is comparable to eccentricity tides raised by the star. Assuming a homogeneous body with a Maxwell rheology, we find that energy dissipation from stellar tides always dominates over that from planet-planet tides across a range of viscosities. TRAPPIST-1 g may experience the greatest proportion of planet-planet tidal heating, where it can account for between 2% and 20% of the total amount of tidal heating, for high ($10^{21}$ Pa s) and low viscosity ($10^{14}$ Pa s) regimes, respectively. If planet-planet tidal heating is to exceed that from stellar eccentricity tides, orbital eccentricities must be no more than e=$10^{-3}$ to $10^{-4}$ for most of the TRAPPIST-1 planets.

Tables

Table 1. Table 1: Relevant geophysical parameters, orbital parameters, are tidal heating results for the TRAPPIST-1 planets. Eccentricity, e 𝑒 e , semimajor axis, a 𝑎 a , and the masses, m 𝑚 m , are from Grimm et al. ( 2018 ) , while the orbital periods, P 𝑃 P and planetary radii, R 𝑅 R , are from Delrez et al. ( 2018 ) . Uncertainties are shown to 1 ​ σ 1 𝜎 1\sigma , and we include uncertainties for the most unconstrained parameters. We also use a stellar mass of 0.089 ± 0.007 ​ M ☉ plus-or-minus 0.089 0.007 subscript 𝑀 ☉ 0.089\pm 0.007M_{\sun} in this work (Delrez et al., 2018 ) . Tidal heating from planet-planet tides, E ˙ p superscript ˙ 𝐸 𝑝 \dot{E}^{p} (Eq. 22 ), and eccentricity tides, E ˙ e superscript ˙ 𝐸 𝑒 \dot{E}^{e} (Eq. 23 ), are shown for a nominal viscosity of ​ 10 21   Pa   s times E21 times pascal second {10}^{21}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s} . The last two columns are estimates of how small the orbital eccentricity must be if tidal heating due to planet-forcing is to become comparable to eccentricity-forced heating, e 0 subscript 𝑒 0 e_{0} , for both a low ( η l ​ o ​ w < ​ 10 14   Pa   s subscript 𝜂 𝑙 𝑜 𝑤 times E14 times pascal second \eta_{low}<${10}^{14}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s}$ ) and high viscosity ( η h ​ i ​ g ​ h > ​ 10 16   Pa   s subscript 𝜂 ℎ 𝑖 𝑔 ℎ times E16 times pascal second \eta_{high}>${10}^{16}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s}$ ) Maxwell body (Figs. 4 and 5 ).

Planet	$P$ [days]	$e$ [${10}^{-2}$]	$m$ [$M_{\oplus }$]	$R$ [$R_{\oplus }$]	${\dot{E}}^p$ [$W$]	${\dot{E}}^e$ [$W$]	$e_0$ [${10}^{-2}$]
							${\eta }_{low}$	${\eta }_{high}$
b	1.519	$0.622\pm 0.304$	$1.017\begin{array}{c}+0.154\\ -0.143\end{array}$	$1.121\begin{array}{c}+0.031\\ -0.032\end{array}$	$8.1\text{×}{10}^{07}$	$9.7\text{×}{10}^{11}$	0.014	0.006
c	2.435	$0.654\pm 0.188$	$1.156\begin{array}{c}+0.142\\ -0.131\end{array}$	$1.095\begin{array}{c}+0.030\\ -0.031\end{array}$	$5.9\text{×}{10}^{07}$	$1.6\text{×}{10}^{11}$	0.039	0.013
d	4.072	$0.837\pm 0.093$	$0.297\begin{array}{c}+0.039\\ -0.035\end{array}$	$0.784\begin{array}{c}+0.023\\ -0.023\end{array}$	$8.5\text{×}{10}^{05}$	$3.2\text{×}{10}^{09}$	0.041	0.014
e	6.135	$0.510\pm 0.058$	$0.772\begin{array}{c}+0.079\\ -0.075\end{array}$	$0.910\begin{array}{c}+0.026\\ -0.027\end{array}$	$3.3\text{×}{10}^{05}$	$9.4\text{×}{10}^{08}$	0.028	0.009
f	9.261	$1.007\pm 0.068$	$0.934\begin{array}{c}+0.080\\ -0.078\end{array}$	$1.046\begin{array}{c}+0.029\\ -0.030\end{array}$	$8.9\text{×}{10}^{05}$	$1.4\text{×}{10}^{09}$	0.069	0.025
g	12.426	$0.208\pm 0.058$	$1.148\begin{array}{c}+0.098\\ -0.095\end{array}$	$1.148\begin{array}{c}+0.032\\ -0.033\end{array}$	$7.1\text{×}{10}^{05}$	$2.9\text{×}{10}^{07}$	0.096	0.032
h	18.871	$0.567\pm 0.121$	$0.331\begin{array}{c}+0.056\\ -0.049\end{array}$	$0.773\begin{array}{c}+0.026\\ -0.027\end{array}$	$4.1\text{×}{10}^{03}$	$3.6\text{×}{10}^{06}$	0.059	0.019

Equations

$${\Phi }^T=\frac{Gm}{r}{\left(\frac{R}{r}\right)}^2\left(\frac{3{\cos }^2\gamma -1}{2}\right)$$

$$r\approx a{\left(1+e\cos M\right)}^{-1}$$

$$\cos \gamma =\sin \theta \cos (\varphi -{\varphi }_t).$$

$${\Phi }^e(\theta ,\varphi ,t)=\frac{3}{2}\frac{G{m}_{\star }}{a}{\left(\frac{R}{a}\right)}^{2}e\left[2\sin M{\sin }^2\theta \sin (2\varphi )+\cos M\left(3{\sin }^2\theta {\cos }^2\varphi -1\right)\right]$$

$${\vec{r}}_i$$

$${\vec{r}}_j$$

$${\vec{r}}_{ij}={\vec{r}}_j-{\vec{r}}_i=\left[a_j\cos (n_{j}t)-a_i\cos (n_{i}t)\right]\hat{x}+\left[a_j\sin (n_{j}t)-a_i\sin (n_{i}t)\right]\hat{y}.$$

$$r_{ij}={\left[a_i^2+a_j^2-2a_i{a}_j\cos (n_{ij}t)\right]}^{1/2}$$

$$\left(\begin{array}{c}\hat{x}\\ \hat{y}\end{array}\right)=\left(\begin{array}{cc}-\cos f_i& \sin f_i\\ -\sin f_i& -\cos f_i\end{array}\right)\left(\begin{array}{c}{\hat{x}}_r\\ {\hat{y}}_r\end{array}\right).$$

$$\cos \gamma =\sin \theta \cos (\varphi -{\varphi }_{ij})$$

$$\cos {\varphi }_{ij}$$

$$\sin {\varphi }_{ij}$$

$$\cos \gamma =\sin \theta \left(\cos \varphi \cos {\varphi }_{ij}+\sin \varphi \sin {\varphi }_{ij}\right)$$

$${\Phi }_{ij}^p(\theta ,\varphi ,t)$$

$$\frac{\delta {\Phi }_i^p}{\delta {\Phi }_i^e}\sim \left(\frac{m_j}{m_{\star }}\right){\left(\frac{a_i}{r_{ij}}\right)}^3\frac{1}{3e_i}$$

$${\Phi }^p(\theta ,\varphi ,t)=\sum _{l=0}^{\infty }\sum _{m=-l}^l{\Phi }_{lm}^p(t)Y_{lm}(\theta ,\varphi ).$$

$$Y_{lm}(\theta ,\varphi )=\{\begin{array}{ll}{P}_{lm}(\theta )\cos (m\varphi ),& \text{if m≥0}\\ P_{l|m|}(\theta )\sin (|m|\varphi ),& \text{if m<0}\end{array}$$

$${\Phi }_{lm}^p(t)=\frac{1}{2}{a}_{lm0}+\sum _{q=1}^{\infty }\left[a_{lmq}\cos (q{n}_{ij}t)+b_{lmq}\sin (q{n}_{ij}t)\right]$$

$$\delta {\Phi }_{lm}^p=\mathrm{Re}\left\{\sum _{q=0}^{\infty }{k}_l(q{n}_{ij})\left(a_{lmq}-i{b}_{lmq}\right)\left[\cos (q{n}_{ij}t)+i\sin (q{n}_{ij}t)\right]\right\}$$

$$\delta\Phi^{p}_{lm}=\sum_{q=0}^{\infty}\operatorname{Re}(k_{l}(qn_{ij}))\Bigl{[}a_{lmq}\cos(qn_{ij}t)+b_{lmq}\sin(qn_{ij}t)\Bigr{]}+\operatorname{Im}(k_{l}(qn_{ij}))\Bigl{[}b_{lmq}\cos(qn_{ij}t)-a_{lmq}\sin(qn_{ij}t)\Bigr{]}$$

$${\dot{E}}_{ij}^p=-\frac{R_i{n}_{ij}}{8{\pi }^{2}G}\sum _{l,m}(2l+1){\int }_{\Omega }{\int }_0^{T_{ij}}\delta {\Phi }_{lm}^p\left({\partial }_t{\Phi }_{lm}^p\right)dtd\Omega ,$$

$${\dot{E}}_{ij}^p=\frac{R_i{n}_{ij}}{2G}\sum _{l,m}\sum _{q=1}^{\infty }\frac{(l+|m|)!}{(l-|m|)!}\frac{1}{(2-{\delta }_{m0})}q\mathrm{Im}(k_l(q{n}_{ij}))(a_{lmq}^2+b_{lmq}^2)$$

$${\dot{E}}_i^e=-\frac{21}{2}\mathrm{Im}(k_2(n_i))\frac{(n_i{R}_i{)}^2}{G}{e}_i^2$$

$$\mathrm{Im}(k_2(q{n}_{ij}))=-k_f\left[\frac{q\bar{\mu }{n}_{ij}{\tau }_M}{1+{\left(q(1+\bar{\mu })n_{ij}{\tau }_M\right)}^2}\right],$$

$$a_{lmq}$$

$$b_{lmq}$$

$${\int }_{\Omega }{Y}_{lm}{Y}_{l^{\prime }{m}^{\prime }}d\Omega =\frac{4\pi {\delta }_{l{l}^{\prime }}{\delta }_{m{m}^{\prime }}}{2l+1}\frac{(l+|m|)!}{(l-|m|)!}\frac{1}{(2-{\delta }_{m0})}$$

$$\mathrm{Im}(k_2(q{n}_{ij}))=-k_f\left[\frac{{\tau }_{M}q{n}_{ij}\bar{\mu }[1+({\tau }_{M}q{n}_{ij}{)}^{1-\alpha }{\zeta }^{-\alpha }S]}{1+({\tau }_{M}q{n}_{ij}{)}^2(\bar{\mu }+1{)}^2+({\tau }_{M}q{n}_{ij}{)}^{2(1-\alpha )}{\zeta }^{-2\alpha }(\alpha !{)}^2+2({\tau }_{M}q{n}_{ij}{)}^{2-\alpha }{\zeta }^{-\alpha }[S/({\tau }_{M}q{n}_{ij})+(\bar{\mu }+1)C]}\right],$$

$$S$$

$$C$$

Full text

gravitation — planets and satellites: dynamical evolution and stability — planets and satellites: individual (TRAPPIST-1g) — planets and satellites: interiors — planets and satellites: terrestrial planets

††journal: ApJ

\savesymbol

tablenum \restoresymbolSIXtablenum

1 Introduction

Tidal heating occurs due to internal friction as a body deforms in response to a time-varying external gravitational potential. It is known to be a dominant process for a number of solar system bodies, such as the Jovian moon Io (e.g., Peale et al., 1979), the small Saturnian satellite Enceladus (e.g., Squyres et al., 1983; Ross & Schubert, 1989; Roberts & Nimmo, 2008), Triton, a retrograde satellite of Neptune (Nimmo & Spencer, 2015), and Earth. Indeed, tidal heating in the Galilean satellites is likely the reason they remain in a stable orbital configuration today.

We usually consider tides raised by the central body on the orbiting body and vice versa. Orbital eccentricity of the secondary object causes it to pass through a time-varying tidal potential, which induces heating due to periodic deformation. Yet, other secondary objects in a system are also sources of time-varying tidal forces. Such tides are typically negligible because the mass of the central tide raising body is usually far greater than other bodies in the system, and also because the distances between these bodies are vast and the strength of tidal forces decreases with the distance between them cubed. The seven planet extrasolar system, TRAPPIST-1 (Gillon et al., 2016, 2017), is the first system to be discovered where this is not the case. The separation distance at conjunction is small enough that tides raised by neighbouring planets can become significant, and heating must occur as a result. While tidal heating due to orbital eccentricity has been investigated for the TRAPPIST-1 system (e.g., Barr et al., 2018), it has never been addressed in detail for planet-planet tides (Wright, 2018).

There are two main aims of this paper. First, we provide the community with the first theory of planet-planet tides and associated heating (Section 2). Second, we apply our theory to the TRAPPIST-1 planets by assuming homogeneous interior structures and Maxwellian rheology (Section 3) to gain a first order understanding of the effect of planet-planet tides in the system and, in particular, how it compares to tidal heating from orbital eccentricity (Section 4).

2 Tidal theory

The distance to a perturbing object, $r$, and the angle $\gamma$ between a position vector $\vec{P}$ on a planet’s surface and the direction from the planet’s centre to the perturber, are two fundamental properties that determine tidal forces acting on a planet. In this section we derive these for tides raised on one planet by another. First, we revisit these for tides raised by a central object.

The tide raising potential due to any object of mass $m$ at distance $r$ from a planet of radius $R$, limited to spherical harmonic degree-2, is;

[TABLE]

where $\gamma$ is the angle between a point on the planet’s surface and the line of centres between the two bodies (Murray & Dermott, 1999). This expression is sufficiently general that it can be used for tides due to a central object, or a neighbouring planet. The only difference between these two cases is how $r$ and $\gamma$ are derived. In the following we use the superscripts $e$ and $p$ to denote properties of eccentricity- and planet-forced tides, respectively.

2.1 Tides raised by the central object

In the reference frame of a planet orbiting a more massive central body and assuming eccentricity $e\ll 1$, the separation distance between two bodies is

[TABLE]

to first order in eccentricity (Murray & Dermott, 1999). Here, $M=nt$ is the mean anomaly of the planet, $n$ is the mean motion, and $t$ is time. We assume that the planet is synchronously rotating, such that its rotation rate is equal to its mean motion.

It can be shown that $\cos\gamma$ is related to the colatitude $\theta$ and longitude $\phi$ of a perturbed planet with zero obliquity via;

[TABLE]

where $\phi_{t}$ is a libration angle of the perturbing body’s position in longitude that arises from small changes in the planet’s orbital speed due to it’s eccentric orbit. To first order in eccentricity, $\phi_{t}\sim 2e\sin M$ (Murray & Dermott, 1999). Substituting this along with Eqs. 3 and 2 into Equation 1 and ignoring the time-independent terms, we arrive at the dynamic tidal potential due to orbital eccentricity (e.g., Kaula, 1964; Wahr et al., 2006):

[TABLE]

where the superscript $e$ denotes an eccentricity-forced potential. This time-varying tidal potential vanishes when $e=0$. Dynamic tides due to neighbouring planets, which we consider in the next section, are unique in that they exist even for circular orbits.

2.2 Tides raised by neighboring planets

In this section we consider the situation shown in Figure 1, where tides are raised on an inner body, planet $i$, by an outer body, planet $j$. Both planets orbit around the host star with coplanar, circular orbits, and have no obliquity. In the barycentric reference frame of the whole system, the position vectors of planets $i$ and $j$ are;

[TABLE]

where $\hat{x}$ and $\hat{y}$ are the inertial unit vectors in the $x$ and $y$ directions (Figure 1). The vector between the two planets is then the difference in these:

[TABLE]

Taking the magnitude of this expression gives us the planet-planet separation distance as a function of time, a fundamental quantity in determining the magnitude of tidal forces;

[TABLE]

where $n_{ij}=n_{i}-n_{j}$ is the conjunction/closest approach frequency between planets $i$ and $j$. As discussed above, we assume the planet is in synchronous rotation, and define planet $i$’s corotating unit vectors $\hat{x}_{r}$ and $\hat{y}_{r}$ to point towards the star and along the trailing hemisphere, respectively (Figure 1). This corotating reference frame is related to the inertial frame via;

[TABLE]

where $f_{i}$ is the true anomaly of planet $i$. As we assume zero inclination and zero obliquity, $\hat{z}=\hat{z}_{r}$. In a similar fashion to eccentricity-forced tides (Eq. 3), $\cos\gamma$ can be expressed as;

[TABLE]

where now $\phi_{ij}$ is the angle between $\hat{x}_{r}$ and the planet-planet vector $\vec{r}_{ij}$ (7). For eccentricity-forced tides, $\phi_{t}$ remains very small over the entire orbital period, whereas $\phi_{ij}$ varies from [math] to $2\pi$ for planet-forced tides, preventing us from making any small angle approximation. Rewriting $\hat{x}_{r}$ in the inertial reference frame using Eq. 9, then taking its dot product with Eq. 7 and rearranging gives the following relationships for $\phi_{ij}$:

[TABLE]

We can then rewrite $\cos\gamma$ in terms of Equations 11 and 12:

[TABLE]

Substituting this expression into Equation 1 gives the tidal potential on planet $i$ due to an outer planet $j$:

[TABLE]

where the superscript $p$ represents tidal forcing due to another planet. The above expression is sufficiently general that tides raised on the outer planet $j$ by the inner planet $i$ can be computed by swapping the subscripts in Equation 14, in which case $n_{ji}$ becomes negative. Note that this expression contains both the static and time-dependent parts of the planet-forced tidal potential.

The potential arising from body $i$’s tidal deformation, commonly referred to as the response potential, is $\delta\Phi_{i}=k_{2}\Phi_{i}$, where $k_{2}$ is the degree-2 potential tidal Love number (Love, 1911). We can use this expression to get some sense of how planet-planet tidal deformation compares to that from eccentricity tides, as was done by Wright (2018). To order of magnitude, using Equations 4 and 14 evaluated at $\theta=\pi/2$, $\phi=0$ and $M=0$, we find;

[TABLE]

which reaches a maximum at planet $i$’s and $j$’s conjunction, $\text{min}(r_{ij})=|a_{i}-a_{j}|$. Note that in this expression (which includes the factor of $1/3$ missing from Wright (2018)), we have assumed that $k_{2}$ is identical for both planet- and eccentricity-forcing, which is not the case in general. Equation 15 is greatest for TRAPPIST-1 g due to tides from planet f, where the maximum tidal distortion may approach $\sim$90\text{,}\mathrm{\char 37\relax}$$ of that from eccentricity-forcing (Table 1). This is illustrated in Figure 2, which shows the magnitude of each forcing potential on TRAPPIST-1 g due to f as a function of time. For the other planets in the system except TRAPPIST-1 b, $\delta\Phi_{i}^{p}/\delta\Phi_{i}^{e}>0.01$, and is sometimes greater than $0.1$, suggesting that in terms of deformation, planet-planet tides may be an important process in the system. This conclusion does not extend to tidal heating, though, where we also have to consider how the tidal forcing and planetary response changes with time.

In Figure 2 we see that planet-planet and eccentricity tides are temporally very different. The eccentricity-forcing on TRAPPIST-1 g operates at a single (orbital) frequency, while in contrast, tides due to planet f are a complex waveform composed of many frequencies. In order to evaluate tidal heating due to planet-forcing, we must decompose the forcing potential into each frequency component because a planet’s tidal response is inherently frequency-dependent. This is the task of the next section.

2.2.1 General form of the tidal potential

The planet-planet tidal potential in Equation 14 can be expanded and written as an infinite sum of spherical harmonics;

[TABLE]

where $\Phi^{p}_{lm}$ are the degree-$l$ and order-$m$ cosine ($m\geq 0$) and sine ($m<0$) spherical harmonic expansion coefficients of the potential $\Phi^{p}$, and the unnoramalised spherical harmonics are;

[TABLE]

where $P_{lm}(\theta)$ is the unnormalised associated Legendre function (Dahlen & Tromp, 1998). At each degree and order, $\Phi^{p}_{lm}(t)$ can be further expressed as an infinite sum over frequency, $q$;

[TABLE]

where $a_{lmq}$ and $b_{lmq}$ are the Fourier expansion coefficients of $\Phi_{lm}^{p}$ for each frequency $q$, as defined in Appendix A.

For each planet in the TRAPPIST-1 system, we calculate $a_{lmq}$ and $b_{lmq}$ by first evaluating the planet-planet tidal potential (Eq. 14) over the planet’s forcing period and surface. At each point in time, we then find $\Phi_{20}$, $\Phi_{22}$ and $\Phi_{2(-2)}$, giving us a time-series of the spherical harmonic expansion coefficients of the forcing potential (Appendix B). Finally, a Fourier transform is applied to that time series to give us the coefficients of the Fourier series. Of these, the only nonzero Fourier series coefficients are $a_{20q}$, $a_{22q}$, and $b_{2(-2)q}$. We provide these coefficients for each planet in the TRAPPIST-1 system as supplementary material to this manuscript.

As an example, Figure 3 shows the normalised frequency spectrum of these coefficients for the tidal potential on TRAPPIST-1 g due to planet f. What is particularly striking is that neither $a_{22q}$ or $b_{2(-2)q}$ peak at the conjunction frequency, $q=1$, but rather at higher frequencies. Additionally, both $a_{20q}$ and $a_{22q}$ have constant ($q=0$) components. Aside from the constant component, $a_{20q}$ is the only coefficient that peaks at the forcing frequency. Importantly, though, Figure 3 illustrates that planet-planet tides are composed of many frequencies, unlike the low-order eccentricity forcing, and much of the planet-planet forcing operates at frequencies higher than the orbital frequency. Consequently, when investigating the impact of planet-planet tides in the TRAPPIST-1 system, it is essential that we adopt a frequency-dependent approach.

2.2.2 The planetary response to tides

As a planet deforms in response to tides, an additional gravitational potential arises as mass is redistributed around the body (Love, 1911). The potential arising from this deformation is commonly referred to as the response potential, $\delta\Phi$, and the magnitude and phase lag between this and the forcing potential controls the amount of tidal heating. The greater the planetary response, the larger $\delta\Phi$ becomes. For planet-planet tides, the spherical harmonic expansion coefficients of the response potential are;

[TABLE]

where the non-subscript $i$ is the imaginary number and $k_{l}(qn_{ij})$ is the complex degree-$l$ potential Love number (Love, 1911) evaluated at the frequency $qn_{ij}$. The real part of the Love number $\operatorname{Re}(k_{l})$ represents the magnitude of the response potential, while the imaginary part $\operatorname{Im}(k_{l})$ corresponds to the phase lag between the forcing and response potentials. Explicitly writing $k_{l}$ in terms of its real and imaginary components, the response potential can be rewritten as;

[TABLE]

In the following section we use this response potential and the forcing in Equation 14 to derive the tidal heating rate due to planet-planet tides.

2.3 Tidal heating

Tidal heating in a solid body is generated via friction as the body responds through deformation to the external forcing potential. Unless the body is perfectly elastic, the response time is nonzero and lags by some amount behind the forcing potential. It is this lag that results in heating via friction at the microphysical scale. As previously mentioned, the lag is often characterized by the imaginary component of the tidal Love number, $k_{l}$, which depends on the forcing frequency for any anelastic material. In the case of planet-planet tides there are also multiple frequencies in the forcing itself (Figure 3), and the tidal response will be different for each of those frequencies. A frequency-dependent approach is then essential when evaluating planet-planet tidal heating.

The time-averaged rate of energy dissipation over the forcing period $T_{ij}=2\pi/n_{ij}$ is given in Zschau (1978, Eq. 18);

[TABLE]

where $\partial_{t}$ represents a time-derivative, $d\Omega=\sin\theta d\theta d\phi$ is the solid angle, and we have collapsed the summation notation for convenience.

If we substitute Equations 18 and 20 into the above and make use of the orthogonality conditions for the spherical harmonics $Y_{lmp}$ (Appendix B), the time-averaged dissipation rate becomes;

[TABLE]

where $\delta_{lm}$ is the Kronecker delta function. We assume these bodies are incompressible, so tidal heating is unaffected by the static forcing and response potential, $q=0$, and is consequently neglected in the frequency summation. For eccentricity-forcing the equivalent expression is (Segatz et al., 1988);

[TABLE]

which is limited to degree-2 and $k_{2}$ is evaluated at the orbital frequency, $n_{i}$. The Fourier series coefficients in Eq. 22, $a_{lmq}$ and $b_{lmq}$, are computed from the planet-planet tidal potential (Eq. 14) as described in the previous section. The only unknowns in this problem are then the frequency-dependent Love numbers for eccentricity- and planet-forcing, which we describe in the next section.

3 Interior Structures and Rheology

A planet’s internal structure and rheology dictate the response of that planet to tides. There are sufficient uncertainties on the masses and radii of the TRAPPIST-1 planets that there is no great advantage in assuming complex interior structures. We instead choose to model these planets as homogeneous bodies with Maxwell rheologies in order to capture the planet’s frequency-dependent response to tides.

The two properties that control the behaviour of a Maxwell material are its viscosity, $\eta$, and rigidity, $\mu$. The Maxwell time, $\tau_{M}=\eta/\mu$, is a fundamental property of such a material, and gives the transition timescale from elastic to viscous deformation. If a Maxwell material is forced on a timescale less than $\tau_{M}$, the material response is elastic. If the opposite is true, then the response is viscous. For average silicate planetary properties of $\mu=$50\text{,}\mathrm{GPa} and $\eta=${10}^{21}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s} (Henning et al., 2009), $\tau_{M}>$ 600 years. The forcing period from either eccentricity or planet-planet tides in the TRAPPIST-1 system is much shorter than this Maxwell time (Table 1), so we would expect the planet response to be elastic. We note, however, that the forcing period can approach the Maxwell time for low viscosities ($\eta\sim${10}^{14}\text{,}\mathrm{Pa}\text{,}\mathrm{s}$$). We choose a Maxwell rheological model in this work because it is the most simple method to capture the frequency-dependent response of a material.

For a homogeneous Maxwell body, the imaginary part of the degree-2 potential Love number at frequency $qn_{ij}$ can be expressed analytically as (e.g., Henning et al., 2009; Renaud & Henning, 2018),

[TABLE]

where $\bar{\mu}=(19/2)\mu/\rho gR$ is the effective rigidity of the homogeneous body with surface gravity $g$, and $k_{f}=3/2$ is the fluid degree-2 Love number (Love, 1911). This expression is sensitive to frequency, and this is illustrated in Appendix C. Using Equation 24, we calculate $\operatorname{Im}(k_{2})$ up to frequency $q=100$ in Eq. 22 to evaluate tidal heating due to planet-planet tides.

4 Results and Discussion

The primary aim of this manuscript is to identify how significant planet-forced tidal heating can be relative to eccentricity forcing. The most significant unknown in this problem is the bulk viscosity of the material, $\eta$, which we vary over several orders of magnitude. Figure 4 shows the tidal dissipated power for each TRAPPIST-1 planet as a function of their bulk viscosities. The last panel shows the ratio of dissipated power from planet-planet tides to eccentricity tides for each planet.

Most notably, Figure 4 shows that tidal heating from planet forcing is always less than eccentricity forcing over the explored viscosity parameter space. This fact does not change over a range of reasonable rigidities, or with different rheologies (Appendix D). We can then conclude, using the eccentricities in Table 1 from Grimm et al. (2018), that planet-planet tidal heating is likely secondary to stellar eccentricity tides. TRAPPIST-1 g is the planet that appears to experience the most significant amount of planet-planet tidal heating, relative to eccentricity forcing.

Interestingly, there are two main solutions to the tidal heating ratio, one in a high viscosity regime ($\eta\gtrsim$ ${10}^{16}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s}$) and the other in a much lower viscosity regime ($\eta\lesssim${10}^{14}\text{,}\mathrm{Pa}\text{,}\mathrm{s}), as highlighted in the bottom right panel of Figure [4](#S4.F4). For all planets, transition from the high to low viscosity regime increases the ratio of planet to eccentricity tidal heating by around one order of magnitude. For TRAPPIST-1 g, this corresponds to planet-planet tides accounting for as low as $2\text{\,}\mathrm{\char 37\relax}$ and up to $20\text{\,}\mathrm{\char 37\relax}$ of the total amount of tidal heating from the high to low viscosity regime, respectively. TRAPPIST-1 h has the second most significant amount of planet-planet tidal heating, which accounts for up to $1\text{\,}\mathrm{\char 37\relax}$ of the total amount of tidal heating, but only for the low viscosity regime. Planet-planet tides are negligible for all other planets, where they account for $<$1\text{\,}\mathrm{\char 37\relax} of the total amount of tidal heating.

Based on our comparison of the maximum tidal deformation in Section 2 (Eq. 15), it is surprising that tidal heating from planet-planet tides is in general so much less than that from stellar eccentricity tides. This is especially true on TRAPPIST-1 g, where the maximum tidal deformation from planet f is nearly equal to that from eccentricity tides. The reason for this difference is because the planet responds far more elastically to planet-forcing than it does to eccentricity-forcing. Much of the power in the frequency spectrum for planet-forcing is at frequencies higher than the orbital frequency (Fig. 3). At these high frequencies, the imaginary part of $k_{2}$ becomes very small, meaning that the planet responds more elastically to the forcing (Appendix C). The more elastic the response is, the less internal friction there is and consequently tidal dissipation drops.

The peaks in tidal heating shown in Figure 4 all occur because, for those viscosities, the Maxwell time becomes comparable to the forcing period. Eccentricity and planet-planet tides have different forcing periods, meaning they each have different viscosities where the Maxwell time approaches the period of forcing. This is why for each planet, the maximum tidal heating occurs at different viscosities for eccentricity and planet-forcing. As a consequence, there is a transition in the tidal heating ratio as the viscosity changes. As discussed above, a planet responds more elastically to planet-planet tides than it does to eccentricity tides in general, which is why the peak in tidal heating occurs at smaller viscosities for planet-forcing than it does for eccentricity-forcing. We observe a similar tidal heating ratio transition for other rheologies (Appendix D).

4.1 Importance of orbital eccentricity

The orbital eccentricities of the TRAPPIST-1 planets are a fundamental property in this problem (Eq. 23). While the most current eccentricities are reasonably well constrained (Grimm et al., 2018), we caution that the eccentricity and mass of a planet can be correlated (e.g., Lithwick et al., 2012). Additionally, the eccentricities of the planets are likely to change significantly over short timescales (e.g., Luger et al., 2017) due to both tidal damping and perturbations from the other planets. Due to this, we also explore how changes in eccentricity affect the planet to eccentricity-forced tidal heating ratio, shown in Figure 5. The dashed lines in Figure 5 represent possible solutions for the planet to eccentricity-forced tidal heating ratio. We explore a low (left) and high (right) viscosity solution.

None of the TRAPPIST-1 planets lie inside the gray region, which is where planet-planet heating is comparable to or greater than eccentricity heating. For a low viscosity of $\eta=${10}^{13}\text{,}\mathrm{Pa}\text{,}\mathrm{s}, TRAPPIST-1 g need only halve its eccentricity to create equal parts planet and eccentricity forced tidal heating. The eccentricities of the other planets must be reduced from the Grimm et al. ([2018](#bib.bib8)) values (Table [1](#S3.T1)) by an order of magnitude or more for this to happen. Planet-planet heating is negligible relative to eccentricity heating for a viscosity of $\eta=${10}^{21}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s}. We conclude that for the majority of the TRAPPIST-1 planets, their eccentricities need to be $12$ orders of magnitude lower than those predicted in Grimm et al. (2018) for planet-planet tidal heating to become at least as significant as eccentricity-driven stellar tidal heating. The exact eccentricities that correspond to a tidal heating ratio of unity, $e_{0}$, are given in Table 1. For most planets this corresponds to $e_{0}=$ ${10}^{-3}{10}^{-4}$. TRAPPIST-1 g is the only planet where planet-forced tidal heating accounts for $>$1\text{,}\mathrm{\char 37\relax} of the total amount of tidal heating throughout the explored viscosity parameter space. Up to $20\text{\,}\mathrm{\char 37\relax}$ of all tidal heating on TRAPPIST-1 g can come from planet-forced tides in the low viscosity regime (assuming $e=$0.0051 (Grimm et al., 2018)).

4.2 Orbital evolution

Tidal heating due to orbital eccentricity reduces the eccentricity of a planet as well as slightly shrinking its semimajor axis (Murray & Dermott, 1999). Planet-planet tides, which exist even in circular orbits, can only shrink the semimajor axis of the tidally heated body. Depending on the relative speed of inward migration between two planets tidally heating each other, an unstable “pile-up” may occur. While possible, we find that the characteristic semimajor axis decay timescale for all TRAPPIST-1 planets, $\tau_{i}=a_{i}/\dot{a}_{i}=Gm_{\star}m_{i}/2a_{i}\dot{E}_{ij}^{P}$ (Murray & Dermott, 1999), is $>$100\text{,}\mathrm{G}\mathrm{y}\mathrm{r}$$ because these planets are so close to their host star. We therefore conclude that planet-planet tides likely have a negligible effect on the orbital evolution of the TRAPPIST-1 planets. We also note that torques from neighboring planets may impact the rotation rate of the tidally distorted planet, but that is beyond the scope of this work.

4.3 Caveats

The absolute amount of tidal dissipation due to planet-forcing is very much of interest, but we caution that this strongly depends on the assumed rheological model and interior structure of the body. Maxwellian rheologies can greatly underestimate the amount of possible tidal dissipation, as has been demonstrated for the Jovian satellite Io (Bierson & Nimmo, 2016; Renaud & Henning, 2018). For that reason we do not focus on the absolute amount of tidal heating, but rather the heating ratio between these two modes of tidal forcing. In Appendix D we present an equivalent version of Figure 4 that uses an Andrade rheology (Andrade, 1910; Jackson & Faul, 2010), which produces much greater amounts of tidal heating from both eccentricity- and planet-forcing in the high viscosity regime. The ratio of planet to eccentricity forced tidal heating in the high viscosity regime also increases slightly for all planets, but this is rather parameter dependent. Furthermore, most rheological models and their associated parameters are derived from laboratory experiments under constant or simple periodic forcing conditions (e.g., Jackson et al., 2004; Jackson & Faul, 2010). While certainly periodic, planet-forcing contains multiple frequency components (Figure 3). Ideally, a rheological model applied to planet-planet tides should therefore be based on laboratory work under the same forcing conditions, which to our knowledge has not yet been performed.

An implicit assumption in our derived tidal heating rate (Eq. 22) is that the system responds linearly to the forcing. In other words, the total response of a planet to tides can be given by the summation of the response at each individual frequency. For the inner two TRAPPIST-1 planets, where planet-planet tides are the most extreme, this assumption should be taken with caution.

Finally, we ignore fluid tides in this work for simplicity. It has been suggested that some of the TRAPPIST-1 planets may have low density envelopes comprising up to $5\text{\,}\mathrm{\char 37\relax}$ their mass (Grimm et al., 2018). If these envelopes are primarily liquid, they may respond more strongly to planet-forced tides than the solid body because dynamical ocean tides have higher natural frequencies in their response (e.g., Kamata et al., 2015; Hay & Matsuyama, 2019).

5 Conclusions

The importance of tidal heating due to tides raised by one planet on another in the TRAPPIST-1 system is investigated in this manuscript. To do this, we derive the tidal potential on a planet due to the gravitational attraction of its neighbours. The potential, which contains may high frequency components, is decomposed first into spherical harmonics and then into Fourier coefficients. We calculate tidal heating on each planet due to tides from every other planet in the TRAPPIST-1 system using these Fourier coefficients and assuming homogeneous Maxwell material interiors. The amount of tidal heating due to neighboring planets is compared to eccentricity-forced tidal heating from the star.

Planet-planet tidal heating is found to always be less than that due to eccentricity tides (assuming orbital eccentricities from Grimm et al. (2018)), because planets respond more elastically to planet-forcing than they do to eccentricity-forcing due to high frequency components in the planet-planet tidal potential. For low viscosities ($<${10}^{14}\text{,}\mathrm{Pa}\text{,}\mathrm{s}), planet-planet tidal heating is around an order of magnitude more significant than if the viscosity were high ($>${10}^{16}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s}), although eccentricity-forced tides still dominate. This transition in tidal heating ratio from a high to low viscosity is because planet-forced tides operate at the planetary conjunction frequency (and higher), while eccentricity-forcing only occurs at the orbital frequency.

Viscosities lower than ${10}^{14}\text{\,}\mathrm{Pa}\text{\,}\mathrm{s}$ seem unlikely unless the body is largely icy or has an extremely high degree of partial melt (e.g., Bierson & Nimmo, 2016; Barr et al., 2018; Renaud & Henning, 2018). We therefore prefer a high viscosity scenario for these planets, which lowers the effectiveness of planet-forced tidal heating relative to eccentricity-forcing. Adopting different rheological models, such as the Andrade rheology, may improve this outlook (Appendix D).

Changes in orbital eccentricity also have a strong effect on the relative importance of planet-forced tidal heating. For TRAPPIST-1 g, its eccentricity from Grimm et al. (2018) must halve for tidal heating from eccentricity- and planet-forcing to become equal. For the other planets, an eccentricity drop of between $12$ orders of magnitude is needed for planet-forced heating to become equal to eccentricity-forced heating.

Overall, planet-planet tidal heating is not found to be significant in the TRAPPIST-1 system when compared to tidal heating due to orbital eccentricity, except for TRAPPIST-1 g. This conclusion is based on our assumption of a homogeneous solid Maxwell body. We neglect any dynamical fluid tides in our model, which may be a source of additional tidal heating and an avenue for future research.

6 Acknowledgments

We thank Renu Malhotra for insightful and encouraging discussions, and Jessie Brown for detailed comments on the manuscript. Some of the ideas presented in this manuscript were formulated through discussion with the participants of the 2018 Keck Institute for Space Studies workshop on tidal heating. This work was supported by the National Aeronautical and Space Agency (NASA) through the Habitable Worlds program (NNH15ZDA001N-HW).

Appendix A Fourier Expansion

We decompose each spherical harmonic coefficient of the tidal potential, $\Phi^{p}_{lm}$ (Eq. 16), into Fourier series coefficients (Section 2). We define the degree-$l$ and order-$m$ Fourier series coefficients at frequency $q$ as;

[TABLE]

where $T_{ij}=2\pi/n_{ij}$ is the forcing/conjunction period. These integrals are evaluated using a Discrete Fast Fourier Transform from the Python library numpy.

Appendix B Spherical Harmonics

B.1 Spherical harmonic coefficients

The tidal potential due to planet forcing (Eq. 14) is evaluated over the reference radius of the deformed body and is then decomposed into spherical harmonics at each point in time. To find the spherical harmonic expansion coefficients at each degree-$l$ and order-$m$, we use the FORTRAN-95 library shtools (Wieczorek & Meschede, 2018).

B.2 Orthogonality

Spherical harmonics are orthogonal over all degrees $l$ and orders $m$ with the normalisation (e.g., Wieczorek, 2015);

[TABLE]

where $d\Omega=\sin\theta d\theta d\phi$ is the solid angle and $\delta_{lm}$ is the Kronecker delta function. This condition is used to derive the energy dissipation expression in Eq. 22.

Appendix C Love numbers of a Maxwell material

The imaginary component of the degree-2 tidal Love number, $\operatorname{Im}(k_{2})$, varies strongly with frequency. We show the $\operatorname{Im}(k_{2})$ of a Maxwell material (Eq. 24) for each TRAPPIST-1 planet as a function of frequency in Figure 6, where the frequency shown is relative to the conjunction frequency with each planet’s nearest neighbour.

The imaginary part of $k_{2}$ decreases rapidly with increasing frequency. For planet-planet tides, which contain high frequency components, this means that the $\operatorname{Im}(k_{2})$ that corresponds to the average frequency in the forcing is much smaller than it is for the conjunction frequency. For this reason, planet-planet tidal heating suffers because the dominant frequency in the forcing corresponds to a smaller $\operatorname{Im}(k_{2})$.

Appendix D Andrade Rheology

Most anelastic materials do not behave in a Maxwellian fashion. A more advanced rheological model that is better suited for tides is the Andrade rheology (e.g., Andrade, 1910; Jackson & Faul, 2010). Recently, Renaud & Henning (2018) showed that tidal dissipation could be much stronger in high viscosity materials when using an Andrade rheology instead of a Maxwell approach. Unfortunately, modeling Andrade materials uses several more free parameters than a Maxwell material does, and some of these are difficult properties to measure in laboratory experiments. The imaginary part of the degree-2 Love number for a homogeneous Andrade body is (e.g., Renaud & Henning, 2018, Table 3);

[TABLE]

where $\alpha$ and $\zeta$ are the Andrade empirical exponent and timescale, respectively, and the Andrade constants are;

[TABLE]

Using nominal values for the empirical constants from Renaud & Henning (2018, Table 1), $\alpha=0.8$ and $\zeta=1$, we calculate $\operatorname{Im}(k_{2})$ using D1 for planet-planet and eccentricity tides to estimate the resulting tidal dissipation (Eqs. 22, 23), shown in Figure 7.

We note two significant differences between the Maxwell (Fig. 4) and Andrade models. Most significantly, the amount of tidal dissipation has increased by several orders of magnitude for both planet- and eccentricity-forcing. Secondly, the different shape of the dissipation curves results in slightly higher ratios between planet-planet and eccentricity tidal heating when the viscosity is $>${10}^{17}\text{,}\mathrm{Pa}\text{,}\mathrm{s}$$. These results are sensitive to both $\alpha$ and $\zeta$, and this could be explored in future work.