# Tidal dissipation in stably stratified and semi-convective regions of
rotating giant planets: incorporating Coriolis forces

Christina M. Pontin  Department of Applied Mathematics, School of Mathematics,
University of Leeds, Leeds, LS2 9JT, UK  Adrian J. Barker Department of
Applied Mathematics, School of Mathematics, University of Leeds, Leeds, LS2
9JT, UK Adrian Barker<EMAIL_ADDRESS>Rainer Hollerbach Department of
Applied Mathematics, School of Mathematics, University of Leeds, Leeds, LS2
9JT, UK<EMAIL_ADDRESS>
(Received Aug 15, 2023; Revised Oct 25, 2023; Accepted Nov 6, 2023)

###### Abstract

We study how stably stratified or semi-convective layers alter tidal
dissipation rates associated with the generation of inertial, gravito-
inertial, interfacial and surface gravity waves in rotating giant planets. We
explore scenarios in which stable (non-convective) layers contribute to the
high rates of tidal dissipation observed for Jupiter and Saturn in our solar
system. Our model is an idealised spherical Boussinesq system incorporating
Coriolis forces to study effects of stable stratification and semi-convective
layers on tidal dissipation. Our detailed numerical calculations consider
realistic tidal forcing and compute the resulting viscous and thermal
dissipation rates. The presence of an extended stably stratified fluid core
significantly enhances tidal wave excitation of both inertial waves (due to
rotation) in the convective envelope and gravito-inertial waves in the dilute
core. We show that a sufficiently strongly stratified fluid core enhances
inertial wave dissipation in a convective envelope much like a solid core
does. We demonstrate that efficient tidal dissipation rates (and associated
tidal quality factors $Q^{\prime}$) – sufficient to explain the observed
migration rates of Saturn’s moons – are predicted at the frequencies of the
orbiting moons due to the excitation of inertial or gravito-inertial waves in
our models with stable layers (without requiring resonance-locking). Stable
layers could also be important for tidal evolution of hot and warm Jupiters,
and hot Neptunes, providing efficient tidal circularisation rates. Future work
should study more sophisticated planetary models that also account for
magnetism and differential rotation, as well as the interaction of inertial
waves with turbulent convection.

Tides (1702) — Astrophysical fluid dynamics (101) — Solar system gas giant
planets (1191) — Extrasolar gaseous giant planets (2172) — Internal waves
(819) — Exoplanet tides (497)

††journal: ApJ

##  1 Introduction

Jupiter and Saturn in our solar system are the best studied giant planets in
the Universe. They have been explored by recent and ongoing space missions
such as Cassini and Juno, thereby allowing us to probe their gravity and
magnetic fields, and to observe fascinating features in their ring systems
(e.g. Bolton et al., 2017; Wahl et al., 2017; Guillot et al., 2018; Durante et
al., 2020; Ingersoll, 2020; Connerney et al., 2022). This wealth of data has
constrained models of the interior structures of these bodies, and we have
found that classical models consisting of a solid rocky core with two
distinct, chemically homogeneous layers above, are inconsistent with
observations. For Jupiter, gravity field measurements from Juno have inferred
a dilute fluid core (which may or may not be stably stratified/non-convective)
containing heavy elements extending out to 40-50% of Jupiter’s radius (Wahl et
al., 2017). Saturn has had an extended stably-stratified layer inferred in its
interior, probably produced by compositional gradients, and possibly extending
out to 60% of the planetary radius, based on features observed in its rings
thought to be produced by resonances with global oscillation modes within the
planet (Marley & Porco, 1993; Hedman & Nicholson, 2013; Fuller, 2014;
Mankovich & Fuller, 2021). These recent findings strongly motivate new studies
of giant planets containing extended stably-stratified fluid layers in their
interiors (e.g. Vazan et al., 2018). Such studies could also shed light on the
structures and properties of extrasolar planets including hot and warm
Jupiters, as well as hot Neptunes (e.g. Guillot et al., 2022).

The dissipative tidal responses of Jupiter and Saturn have been probed using
astrometric observations over the past century or so that measure the orbital
migration of the Galilean and Saturnian moons (Lainey et al., 2009, 2012,
2017, 2020). To explain the observed rates of outward migration of these
moons, we require much more efficient tidal dissipation in these planets than
previously estimated (Goldreich & Soter, 1966). Such efficient rates of tidal
dissipation are not currently understood theoretically, but have motivated an
increasing number of works to explore tides in giant planets. Possible
dissipative mechanisms that have been proposed include: inertial waves
(restored by Coriolis forces) in convective regions (Ogilvie & Lin, 2004; Wu,
2005; Goodman & Lackner, 2009), gravity or gravito-inertial waves in stable
layers (restored by buoyancy forces, and also by rotation) – which might be
locked in resonance (Fuller et al., 2016) – and interactions of equilibrium
tides with turbulent convection, though the latter mechanism is not widely
believed to be important and is more uncertain (Goldreich & Nicholson, 1977;
Duguid et al., 2020; Terquem, 2021; Barker & Astoul, 2021; de Vries et al.,
2023). Finally, dissipation in the visco-elastic rocky/icy cores of giant
planets have also been proposed to be important (e.g. Remus et al., 2012;
Storch & Lai, 2014; Remus et al., 2015; Storch & Lai, 2015; Lainey et al.,
2017), and are worthy of future study, even if fluid mechanisms are typically
favoured currently due to the instability and mixing of such cores in the high
temperature and pressure environments near the centres of giant planets (e.g.
Mazevet et al., 2015).

The efficiency of tidal dissipation is often quantified by the modified tidal
quality factor

$Q^{\prime}=\frac{3}{2k_{2}}\frac{2\pi E_{0}}{\int D\mathrm{d}t},$ (1)

where $k_{2}$ is the quadrupolar Love number (essentially a measure of the
degree of central concentration of mass), $E_{0}$ is the stored tidal energy,
$D$ is the energy dissipation rate, and the integral represents the energy
dissipated over one tidal period. Observations require
$Q^{\prime}\approx(0.94\pm 0.44)\times 10^{4}$ for Saturn and
$Q^{\prime}\approx(1.59\pm 0.25)\times 10^{5}$ for Jupiter, for most tidal
frequencies probed (Lainey et al., 2009, 2012, 2017).

Stabilising compositional gradients can lead to the inhibition of ordinary
convection if they can compete with the destabilising thermal gradients. Since
planetary interiors have more efficient thermal diffusion of heat than viscous
diffusion of momentum (i.e. they have thermal Prandtl numbers smaller than
one, Guillot et al., 2004), convection instead takes the form of double-
diffusive convection (an oscillatory linear instability or “overstability”
involving excitation of internal gravity waves, e.g. Garaud, 2018), whose
nonlinear evolution can produce layering or density staircases, involving
well-mixed steps separated by sharp, diffusive interfaces (Wood et al., 2013).
Whether or not this outcome is expected for conditions relevant in giant
planets is an open question (Fuentes et al., 2022). Nevertheless, the possible
presence and effects of stably-stratified layers in giant planets, which may
or may not involve density staircases, are important to study and are likely
to strongly influence the planet’s tidal response. Indeed, the influence of
interior stably-stratified layers on planetary tidal flows is currently poorly
understood, and has motivated our prior work on this problem (André et al.,
2017, 2019; Pontin et al., 2020, 2023).

We present new theoretical models of giant planets containing stable layers
similar to those constrained observationally for Saturn, and hypothesised for
Jupiter, to explore dissipation of tidal flows inside these giant planets. Our
idealised (Boussinesq) model of a rotating and tidally-forced planet analyses
the dissipative fluid response in a spherical shell using both linear
theoretical analysis and numerical calculations obtained with high-resolution
spectral methods. In this paper, we build upon Pontin et al. (2023) (hereafter
paper 1) to study rotating and stratified planets, by fully incorporating
Coriolis forces to study the effects of rotation on the tidal response. Most
notably, the incorporation of Coriolis forces allows for inertial (restored by
rotation) and gravito-inertial (restored by both buoyancy and rotation) waves
to be excited by tidal forcing and subsequently dissipated. Since Coriolis
forces are important for planetary applications, including the tidal responses
of Jupiter and Saturn at the forcing frequencies of their moons, our model is
more realistic than the one considered in paper 1. However, we continue to
neglect centrifugal deformations and to consider a spherical body for
simplicity. Incorporating centrifugal deformations considerably complicates
the analysis (e.g. Braviner & Ogilvie, 2014, 2015; Barker et al., 2016;
Dewberry & Lai, 2022) but is unlikely to change our results substantially, and
we leave its study to future work.

Our work is complementary to the recent studies of Lin (2023) and Dewberry
(2023) who study more sophisticated compressible and self-gravitating
planetary models numerically, even incorporating centrifugal deformations in
the case of the latter, which are important advances. We have chosen to focus
on a simpler Boussinesq model here to allow a much wider exploration of
parameter space, and also to permit a deeper understanding of the effects of
stable layers, and semi-convective ones, on the tidal response. Our Boussinesq
(incompressible) model incorporates buoyancy forces, and hence allows us to
study gravito-inertial and inertial wave excitation and dissipation, but it
neglects the density variation expected in giant planets. The overall features
of inertial wave excitation in spherical shells are however broadly similar in
incompressible and compressible polytropic models (e.g. Ogilvie, 2013), so
while the particular linear predictions for resonant frequencies are likely to
depend strongly on the particular model, much of the behaviour we observe here
is likely to be similar in more complex models. Here we study the tidal
response of semi-convective layers for the first time in global models with
rotation, and perform a detailed exploration of parameter space for interior
stable layers.

The structure of this paper is as follows. § 2 presents our model, including
the governing equations, energetics, planetary density/entropy profiles and
numerical methods adopted (see paper 1 for further details). § 3 presents an
overview of our numerical exploration of parameter space, analysing the
frequency-dependent tidal response and its dissipative properties. In § 4 we
vary the parameters (rotation rate, core sizes, smooth vs layered density
profiles, diffusivities) and also analyse an integrated measure of the
response, given by the frequency-averaged dissipation rate, to determine how
this correspondingly varies. In § 5, we compare various models of dilute
stably-stratified fluid cores with a rigid core model. Finally, we apply our
results to a Saturn-like model in § 6 and conclude in § 7.

##  2 Model

###  2.1 Governing equations

We briefly recap the model adopted in paper 1 and highlight our extensions to
incorporate rotation. We consider the linearised momentum equation in the
Boussinesq approximation, in the frame rotating at the rate
$\boldsymbol{\Omega}=\Omega\boldsymbol{e}_{z}$,

$\frac{\partial\boldsymbol{u}}{\partial
t}+2\boldsymbol{\Omega}\times\boldsymbol{u}=-\frac{1}{\rho_{0}}\boldsymbol{\nabla}p-b\boldsymbol{g}-\boldsymbol{\nabla}\psi+\nu\boldsymbol{\nabla}^{2}\boldsymbol{u},$
(2)

where $\boldsymbol{u}$, $p$ are velocity and Eulerian pressure perturbations,
$\rho_{0}$ is the reference density and $\psi$ is the tidal potential. We use
spherical polar coordinates $(r,\theta,\phi)$ centred on the planet, and the
rotation axis (along $z$) corresponds to $\theta=0$. We define a (non-
standard) dimensionless buoyancy variable $b=-\rho/\rho_{0}$, where $\rho$ is
the Eulerian density perturbation, and adopt a gravitational acceleration
$\boldsymbol{g}=-g(r)\boldsymbol{\hat{r}}$. The quantity $b$ is proportional
to the entropy perturbation in our Boussinesq model. We focus predominantly on
a homogeneous body with constant density $\rho_{0}$, for which $g=g_{0}r$,
where $g_{0}$ is the surface gravity with $r$ measured in units of the
planetary radius (see paper 1 for exploration of different profiles). Self-
gravity is neglected (i.e. we adopt the Cowling approximation) to allow us to
develop more detailed understanding and because it is only likely to lead to a
moderate linear effect on any quantitative results. The perturbation of the
gravitational potential due to surface deformations of our homogeneous model
would also lead to unrealistically large effects compared with more realistic
centrally-condensed models, further motivating our neglect of its effects
here. The flow is assumed to be incompressible, so

$\boldsymbol{\nabla}\cdot\boldsymbol{u}=0.$ (3)

The heat equation written in terms of $b$ is

$\frac{\partial b}{\partial
t}+\frac{u_{r}}{g}N^{2}=\kappa\boldsymbol{\nabla}^{2}b,$ (4)

where $N^{2}$ is the (radially varying) squared Brunt-Väisälä or buoyancy
frequency, defined to be,

$N^{2}=g\left(\frac{1}{\Gamma_{1}}\frac{\mathrm{d}\ln
p_{0}}{\mathrm{d}r}-\frac{\mathrm{d}\ln\rho_{0}}{\mathrm{d}r}\right)\propto\frac{\mathrm{d}s}{\mathrm{d}r},$
(5)

where $\Gamma_{1}=\left(\frac{\partial\ln
p_{0}}{\partial\ln\rho_{0}}\right)_{\mathrm{ad}}$ is the first adiabatic
exponent and $s(r)$ is the specific entropy.

We have adopted a single scalar field for the density ($b$) with a single
“thermal diffusivity” ($\kappa$) to model our density/entropy stratification.
If compositional gradients are important then the density should strictly
consist of two distinct components (thermal and compositional) with two
associated (unequal) diffusivities. This is necessary for double-diffusive
convection to operate in giant planet interiors and potentially produce
layered density structures (e.g. Garaud, 2018), but we leave exploration of
double-diffusive effects here to future work.

We introduce tidal forcing by considering the dominant component of the tidal
potential

$\psi=\psi_{0}~{}r^{2}Y^{2}_{2}(\theta,\phi)e^{-i\omega t},$ (6)

where
$\psi_{0}=\sqrt{\frac{6\pi}{5}}\frac{M_{2}}{M}\omega_{d}^{2}\big{(}\frac{r_{0}}{a}\big{)}^{3}$.
We take the $l=m=2$ dimensional tidal amplitude for a circular orbit, as
defined in Ogilvie (2014), where $M$ is the mass of the planet, $M_{2}$ is the
companion mass, $a$ is the orbital semi-major axis, and the dynamical
frequency $\omega_{d}=\sqrt{GM/r_{0}^{3}}$, where $G$ is the gravitational
constant. The forcing frequency is $\omega=2(\Omega_{o}-\Omega)$, where
$\Omega_{o}$ is the orbital frequency of the satellite and $\Omega$ is the
spin frequency of the planet – the most relevant one for a circular and
aligned orbit of a non-synchronously orbiting moon. Note that we solve for the
entire linear tidal response to $\psi$ directly; we do not split up the tide
into an equilibrium and a dynamical tide, though only the dynamical/wavelike
response is typically important for dissipation in our models.

We expand perturbations using spherical harmonics with a harmonic time-
dependence such that,

$u_{r}(r,\theta,\phi,t)=\sum_{l=m}^{\infty}\tilde{u}_{r}^{l}(r)Y^{m}_{l}(\theta,\phi)e^{-i\omega
t},$ (7)
$u_{\theta}(r,\theta,\phi,t)=r\sum_{l=m}^{\infty}\bigg{[}\tilde{u}_{b}^{l}(r)\frac{\partial}{\partial\theta}+\frac{\tilde{u}_{c}^{l}(r)}{\sin\theta}\frac{\partial}{\partial\phi}\bigg{]}Y^{m}_{l}(\theta,\phi)e^{-i\omega
t},$ (8)
$u_{\phi}(r,\theta,\phi,t)=r\sum_{l=m}^{\infty}\bigg{[}\frac{\tilde{u}_{b}^{l}(r)}{\sin\theta}\frac{\partial}{\partial\phi}-\tilde{u}_{c}^{l}(r)\frac{\partial}{\partial\theta}\bigg{]}Y^{m}_{l}(\theta,\phi)e^{-i\omega
t},$ (9)
$p(r,\theta,\phi,t)=\sum_{l=m}^{\infty}\tilde{p}^{l}(r)Y^{m}_{l}(\theta,\phi)e^{-i\omega
t},$ (10)
$b(r,\theta,\phi,t)=\sum_{l=m}^{\infty}\tilde{b}^{l}(r)Y^{m}_{l}(\theta,\phi)e^{-i\omega
t},$ (11)

where we adopt the standard normalisation

$\int_{0}^{2\pi}\int_{0}^{\pi}[Y_{l^{\prime}}^{m^{\prime}}(\theta,\phi)]^{*}Y_{l}^{m}(\theta,\phi)\sin^{2}\theta\,\mathrm{d}\theta\,\mathrm{d}\phi=\delta_{l,l^{\prime}}\,\delta_{m,m^{\prime}}.$

The resulting differential equations in radius for each $l$ and $m$ are
(dropping tildes)

$\displaystyle(-i\omega)u_{r}^{l}+2\Omega
r\Big{(}-imu_{b}^{l}+(l-1)q_{l}u_{c}^{l-1}-(l+2)q_{l+1}u_{c}^{l+1}\Big{)}$
$\displaystyle=-\frac{1}{\rho_{0}}\frac{\mathrm{d}p^{l}}{\mathrm{d}r}+gb^{l}-\frac{\mathrm{d}\psi^{l}}{\mathrm{d}r}\delta_{l,2}-\nu\frac{l(l+1)}{r^{2}}\Big{[}u_{r}^{l}-\frac{\mathrm{d}}{\mathrm{d}r}(r^{2}u_{b}^{l})\Big{]},$
(12) $\displaystyle(-i\omega)r^{2}u_{b}^{l}+2\Omega
r^{2}\left(\frac{-im}{l(l+1)}\left(\frac{u_{r}^{l}}{r}+u_{b}^{l}\right)\right.$
$\displaystyle\left.\hskip
71.13188pt+\frac{l-1}{l}q_{l}u_{c}^{l-1}+\frac{l+2}{l+1}q_{l+1}u_{c}^{l+1}\right)$
$\displaystyle\hskip
56.9055pt=-\frac{p^{l}}{\rho_{0}}-\psi^{l}\delta_{l,2}+\nu\,\left[\frac{2u_{r}^{l}}{r}+\frac{1}{r^{2}}\frac{\mathrm{d}}{\mathrm{d}r}\bigg{(}r^{4}\frac{\mathrm{d}u_{b}^{l}}{\mathrm{d}r}\bigg{)}\right.$
$\displaystyle\left.\hskip 85.35826pt-(l-1)(l+2)u_{b}^{l}\right],$ (13)
$\displaystyle(-i\omega)r^{2}u_{c}^{l}+2\Omega
r^{2}\bigg{(}\frac{-im}{l(l+1)}u_{c}+\frac{q_{l}}{l}\frac{u_{r}^{l-1}}{r}-\frac{q_{l+1}}{l+1}\frac{u_{r}^{l+1}}{r}$
$\displaystyle\hskip
85.35826pt+\frac{l-1}{l}q_{l}u_{b}^{l-1}-\frac{l+2}{l+1}q_{l+1}u_{b}^{l+1}\bigg{)}$
$\displaystyle\hskip
34.14322pt=\nu\Bigg{[}\frac{1}{r^{2}}\frac{\mathrm{d}}{\mathrm{d}r}\bigg{(}r^{4}\frac{\mathrm{d}u_{c}^{l}}{\mathrm{d}r}\bigg{)}-(l-1)(l+2)u_{c}^{l}\Bigg{]},$
(14)
$\frac{1}{r^{2}}\frac{\mathrm{d}}{\mathrm{d}r}(r^{2}u_{r}^{l})-l(l+1)u_{b}^{l}=0,$
(15)
$(-i\omega)b^{l}+N^{2}\frac{u_{r}^{l}}{g}=\kappa\,\Bigg{[}\frac{1}{r^{2}}\frac{\mathrm{d}}{\mathrm{d}r}\bigg{(}r^{2}\frac{\mathrm{d}b^{l}}{\mathrm{d}r}\bigg{)}-\frac{l(l+1)}{r^{2}}b^{l}\Bigg{]},$
(16)

where $q_{l}=\big{(}\frac{l^{2}-m^{2}}{4l^{2}-1}\big{)}^{\frac{1}{2}}$. These
are consistent with the equations in Ogilvie (2009) when $b=0$. The equations
for each $m$ are uncoupled due to the axisymmetric basic state, but the
Coriolis term couples components with different angular wavenumbers $l$.
Hence, we solve our system for $m=2$ and $l\in[2,l_{max}]$, where $l_{max}$ is
taken to be sufficiently large such that the solution converges. Typically,
$l_{max}=100$, though larger values are required for the smallest
diffusivities explored.

We non-dimensionalise our system with the planetary radius $r_{0}$ as the unit
of length, the reference density $\rho_{0}$ for density, and $\omega_{d}^{-1}$
as the time unit. Therefore, we introduce the following dimensionless
variables: $r=r_{0}\hat{r}$, $u_{r}=r_{0}\omega_{d}\hat{u}_{r}$,
$u_{b}=\omega_{d}\hat{u}_{b}$, $u_{c}=\omega_{d}\hat{u}_{c}$,
$\psi=\psi_{0}\hat{\psi}$, $b=\hat{b}$, $g=g_{0}\hat{g}$,
$p=\rho_{0}r_{0}g\hat{p}$, $\omega=\sqrt{\frac{g_{0}}{r_{0}}}\hat{\omega}$,
where $g_{0}=\omega_{d}^{2}r_{0}$ is the surface gravity. We henceforth drop
hats on all variables. This leads to four dimensionless parameters for each
frequency $\omega$ (in units of $\omega_{d}$),

$\small\frac{\Omega}{\omega_{d}},~{}~{}~{}\frac{N^{2}}{\omega_{d}^{2}},~{}~{}~{}\frac{\nu}{r_{0}^{2}\omega_{d}}~{}~{}~{}\textrm{and}~{}~{}~{}\frac{\kappa}{r_{0}^{2}\omega_{d}}~{}~{}~{}\bigg{(}\textrm{alternatively
}\mathrm{Pr}=\frac{\nu}{\kappa}\bigg{)}.$

These are varied in later analysis, and are referred to simply as $\Omega$,
$N^{2}$ (which can depend on $r$ with mean value $\bar{N}^{2}$), $\nu$ and Pr,
since we use units such that $r_{0}=\omega_{d}=1$.

The molecular values of $\nu$ and $\kappa$ expected in giant planets, as well
as the values of Pr, are somewhat uncertain and are likely to depend on
radius. It is believed that Pr ranges from 0.01 to 0.1 in deep giant planet
interiors (e.g. Guillot et al., 2004; French et al., 2012, albeit for
Jupiter), with smaller values in the thin outermost atmosphere. Values of the
Ekman number $\nu/(r_{0}^{2}\Omega)$ are likely to be very small, on the order
of $10^{-15}$ or smaller using molecular values. It is not clear whether
molecular values are appropriate though or whether tidal waves should be
damped by turbulent diffusivities instead, which may be expected to be much
larger (though still small, as discussed in e.g. de Vries et al., 2023). We
assume $\nu$ and $\kappa$, and hence Pr, are constant in radius, and explore
the widest range of these parameters that is computationally accessible,
though we are unable to study values as small as the microscopic ones in
planetary interiors – common with many other problems in astrophysics, thereby
requiring us to understand how $\nu$ and $\kappa$ affect our results before we
can potentially extrapolate to real planets. Another uncertainty is if there
are semi-convective layers, turbulent motions within convective layers may
lead to $\mathrm{Pr}\sim 1$, whereas the diffusive interfaces between them may
have molecular values of $\mathrm{Pr}\sim 10^{-2}$. This matter is very
speculative however, so we prefer to study models with constant diffusivities
in this work.

###  2.2 Boundary Conditions

We assume an idealised, perfectly rigid, solid core that is impermeable at
$r=\alpha r_{0}$, thus

$u_{r}^{l}=0\quad\mathrm{at}\quad r=\alpha r_{0}.$ (17)

At the planetary radius $r=r_{0}$ there is a perturbed free surface with
vanishing normal stress, requiring (paper 1)

$W^{l}-\frac{g_{0}}{(-i\omega)}u_{r}^{l}-2\nu\frac{\mathrm{d}u_{r}^{l}}{\mathrm{d}r}=\psi_{0}\delta_{l2}\quad\mathrm{at}\quad
r=r_{0},$ (18)

where $W^{l}=\frac{p^{l}}{\rho_{0}}+\psi$. At both boundaries, we apply
stress-free conditions (no tangential stress). This is an approximation at
$r=\alpha r_{0}$, but is less computationally costly than applying a no-slip
condition, and it is unlikely to produce significant differences in our
results. Hence, we require

$\frac{\mathrm{d}u_{b}^{l}}{\mathrm{d}r}+\frac{u_{r}^{l}}{r^{2}}=\frac{\mathrm{d}u_{c}^{l}}{\mathrm{d}r}=0\quad\mathrm{at}\quad
r=\alpha r_{0}\quad\mathrm{and}\quad r=r_{0}.$ (19)

The inner core has fixed entropy, such that

$b^{l}=0\quad\mathrm{at}\quad r=\alpha r_{0},$ (20)

and we adopt vanishing perturbations to the buoyancy flux through the surface,
i.e.

$\frac{\partial b^{l}}{\partial r}=0\quad\mathrm{at}\quad r=r_{0}.$ (21)

Our results are not particularly sensitive to choices of thermal boundary
conditions though.

###  2.3 Energetics

We define the mean rate of energy injection by the tidal forcing

$I=\int_{V}\rho_{0}\boldsymbol{u}\cdot(-\nabla\psi)\mathrm{d}V=-\oint_{S}\rho_{0}\psi\boldsymbol{u}\cdot\mathrm{d}\boldsymbol{S}.$
(22)

By taking the scalar product of equation 2 with $\rho_{0}\boldsymbol{u}$,
using equation 4 and integrating over the volume, we find

$I=\frac{\mathrm{d}E_{K}}{\mathrm{d}t}+\frac{\mathrm{d}E_{PE}}{\mathrm{d}t}+\frac{1}{V}\oint_{S}p\boldsymbol{u}\cdot\mathrm{d}\boldsymbol{S}+D_{ther}+D_{visc},$
(23)

where

$E_{K}=\int_{V}\frac{\rho_{0}}{2}|\boldsymbol{u}|^{2}\mathrm{d}V,$ (24)
$E_{PE}=\int_{V}\frac{g^{2}}{2N^{2}}b^{2}\mathrm{d}V,$ (25)

are the kinetic and potential energies of the system, respectively, except
when $N=0$, in which case $E_{PE}=0$ and
$\frac{\mathrm{d}E_{PE}}{\mathrm{d}t}=0$. The volume integrated viscous
dissipation rate is given by

$D_{visc}=-\int_{V}\rho_{0}\nu\boldsymbol{u}\cdot\nabla^{2}\boldsymbol{u}\,\mathrm{d}V,$
(26)

and volume integrated thermal dissipation rate is written

$D_{ther}=-\int_{V}\rho_{0}\kappa\frac{g^{2}}{N^{2}}b\nabla^{2}b\,\mathrm{d}V,$
(27)

and when $N=0$, $D_{ther}=0$.

$N$ $s$$r$$0$$\alpha r_{0}$$\beta r_{0}$$r_{0}$Inner CoreOuter CoreConvective
RegionPlanet Radius$N$ $s$$r$$0$$\alpha r_{0}$$\beta r_{0}$$r_{0}$

Figure 1: Left panels: illustrative entropy profile and Brunt-Väisälä
frequency ($N^{2}$) profiles for a continuously stratified layer and for semi-
convective layers. A solid inner core extends to $\alpha r_{0}$ and the stable
layer (“dilute core”) extends to $\beta r_{0}$, outside of which is a
convective envelope to planetary radius $r_{0}$. Right panels: examples of the
numerical profiles adopted.

In paper 1, we showed how the viscous dissipation can be separated into bulk
and boundary contributions. The former can be written (e.g. Ogilvie, 2009),

$D_{interior}=\frac{\rho_{0}\nu}{2}\sum_{l}\Bigg{(}l(l+1)\Bigg{(}\bigg{|}\frac{u_{r}^{l}}{r}+r\frac{\mathrm{d}u_{b}^{l}}{\mathrm{d}r}\bigg{|}^{2}+\bigg{|}r\frac{\mathrm{d}u_{c}^{l}}{\mathrm{d}r}\bigg{|}^{2}\Bigg{)}\\\
\hskip
14.22636pt+3\bigg{|}\frac{\mathrm{d}u_{r}^{l}}{\mathrm{d}r}\bigg{|}^{2}+(l-1)l(l+1)(l+2)(|u_{b}^{l}|^{2}+|u_{c}^{l}|^{2})\Bigg{)},$
(28)

and is typically by far the dominant contribution for the parameters we study.
We subsequently refer to the viscous dissipation using equation 26, but the
differences are negligible for our purposes.

In a steady state, averaged over the tidal period $2\pi/\omega$, the injection
rate $I$ balances the total viscous and thermal dissipation
$D_{total}=D_{visc}+D_{ther}$ (ignoring a small boundary term discussed in
paper 1). In this study we will analyse how the dissipation rates $D_{visc}$
and $D_{ther}$ depend on the properties of the system.

###  2.4 Modelling stratified layers

We consider stably-stratified regions in a giant planet, or ones with semi-
convective layers (“density staircases”), as might be produced by double-
diffusive convection. We model stable layers both near a small solid core to
model a dilute stably-stratified fluid core, and ones further out to represent
a layer near a potential H/He molecular to metallic transition where helium
rain may occur. Firstly, we model a continuous stably stratified region, with
a constant $N$ which sits above the solid core extending a defined distance
into the planetary envelope. This is shown on the left-hand side of Figure 1
and is described by a step function in $N^{2}(r)$ that is non-zero from
$\alpha r_{0}$ to $\beta r_{0}$. When $\beta\neq 1$, we consider a smoothly
varying profile,

$N^{2}(r)=\frac{\bar{N}^{2}}{2}\Big{(}\tanh\big{(}\Delta(\beta
r_{0}-r)\big{)}+1\Big{)},$ (29)

with $\Delta=100$ unless otherwise specified, and set $N^{2}=0$ when
$N^{2}(r)<\frac{\bar{N}^{2}}{10^{7}}$. For cases where $\beta=1$ we simply set
$N^{2}(r)=\bar{N}^{2}$.

To model a semi-convective “dilute fluid core” with $n_{max}$ steps within the
layer we consider a series of $\delta$-like-functions for $N^{2}$ to give a
staircase-like density profile, as shown on the right-hand side of Figure 1.
We use

$N^{2}(r)=\begin{cases}\frac{N_{0}^{2}}{2}\bigg{(}1+\cos\Big{(}2\pi\frac{r-r_{n}}{\delta
r}\Big{)}\bigg{)}\quad&|r-r_{n}|<\frac{\delta r}{2},\\\
0\quad&\text{otherwise},\end{cases}$ (30)

for $1<n<n_{max}$, where $r_{n}=\alpha r_{0}+nd$,
$d=\frac{(\beta-\alpha)r_{0}}{n_{max}}$, and $\epsilon=\delta r/d$.
$N_{0}^{2}$ is set to a value which gives a mean stratification equivalent to
a constant stratification $\bar{N}^{2}$,
$N_{0}^{2}=\bar{N}^{2}\frac{(\beta-\alpha)r_{0}}{n_{max}~{}\delta r}.$ This
allows for comparisons between stratified layers and semi-convective regions.
If $\beta=1$, we omit the final step at the planetary radius (i.e.
$0<n<n_{max}-1$), and $N_{0}$ is altered accordingly to maintain the mean
stratification.

To model a stratified layer at the metallic/molecular transition zone, we
consider one wide “interface”:

$N^{2}(r)=\begin{cases}\frac{N_{0}^{2}}{2}\bigg{(}1+\cos\Big{(}2\pi\frac{r-\beta}{\delta
r}\Big{)}\bigg{)}\quad&|r-\beta|<\frac{\delta r}{2},\\\
0\quad&\text{otherwise},\end{cases}$ (31)

and $N_{0}^{2}=\bar{N}^{2}\frac{(r_{0}-\alpha)r_{0}}{\delta r}$. This models
an isolated stable layer embedded within a convective medium.

###  2.5 Frequency averaged dissipation

The tidal response is expected to be strongly frequency-dependent when
wavelike tides are excited, but we have found it helpful to define a single
quantitative measure of the tidal dissipation that can be compared as we vary
our parameters. To do so, we define a frequency-averaged dissipation measure
for a given dissipation rate $D$

$\bar{D}=\int^{\omega_{max}}_{\omega_{min}}\frac{D(\omega)}{\omega}~{}\mathrm{d}\omega,$
(32)

which gives more emphasis to lower frequencies. Numerically, we adopt a small
non-zero lower bound $\omega_{min}$, and unless otherwise stated,
$\omega_{max}=\bar{N}$. This allows analysis of the low frequency regime while
slightly reducing the contribution of the surface gravity modes (f-modes) that
are not our primary focus. We also use the different weighting

$\bar{D}_{2}=\int^{\omega_{max}}_{\omega_{min}}\frac{D(\omega)}{\omega^{2}}~{}\mathrm{d}\omega.$
(33)

This quantity is directly relevant for comparison with Ogilvie (2013), where
in the low frequency limit, for the unstratified case with a solid core, they
used impulsive forcing to calculate the associated energy transfer into
inertial modes $\hat{E}$ analytically. This quantity is related to tidal
dissipation rates $D$ by

$\hat{E}=\frac{1}{2\pi}\int^{\infty}_{-\infty}\frac{D}{\omega^{2}}~{}\mathrm{d}\omega,$
(34)

which provides motivation for using an $\omega^{-2}$ weighting. Other
weighting factors could be used instead, and it is not at all clear what is
the most useful one for astrophysical purposes, so we will explore various
possibilities.

###  2.6 Numerical method for the forced problem

We solve the system of ordinary differential equations in radius, equations 12
to 16 for each $l$, using a Chebyshev collocation method, where the ordinary
differential equations in $r$ are converted into a linear system of equations
on a Chebyshev extrema grid. We consider points in radius as a set of
$(n_{cheb}+1)$ Gauss-Lobatto-Chebyshev points, which are defined as,

$x_{j}=\cos\bigg{(}\frac{j\pi}{N_{c}}\bigg{)},\quad j=0,\cdots,N_{c},$ (35)

where $N_{c}$ is taken as the (smallest) appropriate value for which numerical
convergence is found. This value varies depending on the parameters, but we
typically take $N_{c}=100$ to $N_{c}=400$. This method is well-suited for many
non-periodic problems as it is a spectral method free of the Runge phenomenon,
which converges exponentially fast with resolution $N_{c}$ for smooth
solutions (Boyd, 2001). We then have a linear algebra problem which we solve
using the inbuilt MATLAB routine “mldivide”, where matrices are stored in
sparse form to reduce memory requirements. The solutions for $\boldsymbol{u}$,
$p$ and $b$ can then be used in equations 26 and 27.

###  2.7 Eigenvalue problem

In paper 1 (see also Pontin et al., 2020) we studied analytically the non-
diffusive free gravity modes in our model without rotation, and we also
explored numerically the eigenvalues of the dissipative system. In the
presence of rotation, analytical progress is severely hampered by the
$l$-coupling caused by Coriolis forces, so in almost all cases we must resort
to numerical calculations to obtain eigenfrequencies. This is particularly the
case when inertial modes are excited, i.e. for $|\omega|<2\Omega$ in
convective regions (or for gravito-inertial modes in stable layers), where the
inviscid (non-diffusive) problem becomes mathematically ill-posed. To compute
the eigenvalues (free modes) of the dissipative problem numerically, we set
$\psi_{0}=0$, and manipulate our system of equations into a linear generalised
eigenvalue problem with eigenvalue $(-i\omega)$, i.e.,

$(-i\omega)\mathrm{A}\begin{bmatrix}\boldsymbol{u}\\\ p\\\
b\end{bmatrix}=\mathrm{B}\begin{bmatrix}\boldsymbol{u}\\\ p\\\
b\end{bmatrix},$ (36)

where $\mathrm{A}$ and $\mathrm{B}$ are matrices that describe equations 12 to
16 with the same boundary conditions as above. This is solved numerically to
obtain eigenvalues ($-i\omega$), and corresponding eigenvectors for
$\boldsymbol{u}$, $p$, $b$. We use the iterative eigs solver in MATLAB to scan
the relevant frequency range, as a non-iterative method would be prohibitive
in its memory requirements for even modest spatial resolutions, since this is
a two-dimensional eigenvalue problem in $r$ and $\theta$. The frequency of a
mode is $\operatorname{Re}[\omega]$ and its corresponding damping rate is
$\operatorname{Im}[\omega]$ (there are no unstable modes).

##  3 Overview of the tidal response of rotating planets

In paper 1 we studied in detail the response of non-rotating planets
containing stable layers (or semi-convective layers), which consists of
internal, surface and interfacial gravity modes. When incorporating rotation,
the additional excitation of low frequency inertial waves in convective
(neutrally stratified) regions, and the modification of internal gravity modes
to become gravito-inertial modes, means we see new resonances at low
frequencies that align with these modes compared with calculations without
rotation.

We illustrate the dissipative response of a rotating planet as a function of
frequency in Figure 2(a). These plots show the viscous ($D_{visc}$), thermal
($D_{therm}$), and total ($D_{total}$) dissipation rates. In this example, we
show results for a rotating body with a large solid core of size $\alpha=0.5$
surrounded by a well-mixed convective fluid envelope with $N=0$, thereby
allowing pure inertial waves to be excited (hence $D_{therm}=0$). We adopt
$\Omega=0.4$ here to approximately model Saturn’s rotation rate, and take
$\nu=\kappa=10^{-6}$ for this illustration (for molecular diffusion we expect
much smaller values in reality but this value is readily accessible
computationally). For this example we only plot positive frequencies, but note
that with rotation the response at a frequency $\omega$ will in general differ
from that at $-\omega$.

Pure inertial waves are excited at low frequencies for which
$|\omega|<2\Omega=0.8$, and surface gravity modes associated with the free
surface are excited at higher frequencies $|\omega|\gtrsim 1$. We observe
peaks of enhanced dissipation close to resonances with inertial modes for low
frequencies, which produce the irregular frequency dependence characteristic
of these modes. Here we use a large core size to enhance the appearance of
inertial waves in the tidal response, to more clearly illustrate their
properties (Ogilvie, 2009; Rieutord, 2009; Rieutord & Valdettaro, 2010),
though this can be thought to model strongly stably stratified dilute cores,
as we will illustrate later. The numerically-obtained eigenvalues (purple
squares) agree well with the locations and heights of the peaks, where the
least damped modes typically correspond to the tallest resonances, with the
strongest dissipation. We have plotted only the least damped modes, with a
cut-off value of $\mathrm{Im}\left[\omega\right]$ chosen for aesthetic
reasons.

(a) $\alpha=0.5$, $\bar{N}^{2}=0$

(b) $\alpha=0.5$, $\bar{N}^{2}=0$

(c) $\alpha=0.1$, $\beta=0.5$, $\bar{N}^{2}=1$

Figure 2: Illustrative examples of tidal dissipation as a function of
frequency in two rotating models, restricted to positive frequencies in panels
(a) and (c), with numerically-obtained eigenvalues over-plotted for both
models. The $y$-axis for the eigenvalues indicates their decay rates
(imaginary parts). One example of an entirely convective fluid envelope
outside a solid inner core (panels 2(a) and 2(b), where the latter shows
negative frequencies) and one with a stably stratified layer extending to half
the planetary radius and a tiny solid core (panel 2(c)). In both cases we set
$\Omega=0.4$ and $\nu=\kappa=10^{-6}$.

As in the non-rotating cases in paper 1, we observe a strong surface gravity
mode resonance around $\omega\sim 1$, but here it is shifted to a
significantly lower frequency. We show in Appendix A (see also Lebovitz, 1961;
Braviner & Ogilvie, 2014; Barker et al., 2016), that rotation causes the
splitting of this mode for a given $l$ for different $m$ values. In the limit
of slow rotation, the $l=2$ surface gravity mode is split into five modes,
with frequencies $\omega=\sqrt{2}\omega_{d}-\frac{m}{2}\Omega$, where
$m=-2,-1,0,1,2$ instead of the non-rotating frequency $\sqrt{2}\omega_{d}$ (in
agreement with Lebovitz, 1961, for a Maclaurin spheroid in this limit). As we
consider $m=2$ forcing, we only observe the $m=2$ surface gravity mode, with a
frequency $\omega\approx 1.07$ that is nicely predicted by (A10) (the small
$\Omega$ expression predicts $1.01$). There is an additional resonance close
to $\omega=1.8$, caused by the coupling of different harmonic degrees
appearing close to a resonance with the $l=4$, $m=2$, surface gravity mode
(again shifted due to splitting). Around $\omega=1.5$ there is an eigenvalue
solution that does not align with any peak, corresponding to the $l=3$ surface
gravity mode. However, this mode is not actually excited, even though it is a
solution to the unforced eigenvalue problem. The reason for this is that while
all $l$’s are coupled by the Coriolis force, the equatorial symmetry imposed
by the $l=2$ tidal forcing means that only modes with this symmetry
(equatorially symmetric for $u_{r}$, therefore only with even $l$’s for $m=2$)
are excited.

Figure 2(c) shows a similar but potentially more realistic case where instead
of a large solid core there is a stably stratified layer extending to half the
planetary radius, also with $\Omega=0.4$ and $\nu=\kappa=10^{-6}$. We consider
a buoyancy profile as in equation (29) with $\alpha=0.1$, $\beta=0.5$, and
$\bar{N}=1$. This case is the same as Fig. 2(b) in paper 1 but with the
addition of rotation. In this case, as well as exciting inertial waves in the
outer convective envelope, as observed in the solid core case, we also excite
gravito-inertial waves in the extended stably stratified layer over a larger
range of frequencies111The range of frequencies for gravito-inertial wave
propagation is explained in e.g. Rieutord (2009). Alternatively, please see
André et al. (2017) equations 2.14 (or 2.33) for this upper bound in the case
of plane waves when the rotation axis is perpendicular to the local gravity
vector. $|\omega|<\sqrt{\bar{N}^{2}+4\Omega^{2}}=1.28$.

(a) $\omega=0.181$

(b) $\omega=0.561$

Figure 3: Illustrative examples of $|\boldsymbol{u}|$ for two different
forcing frequencies of the forced response in a convective envelope with
$\alpha=0.5$, $\bar{N}^{2}=0$, $\Omega=0.4$, and $\nu=\kappa=10^{-6}$. Both
panels show examples of tidally-forced inertial waves.

In Figures 3 and 4 we show the spatial structure of the forced response, which
displays features similar to those reported in previous studies for the free
and forced modes in rotating and stratified planets (e.g. Rieutord &
Valdettaro, 1997; Dintrans et al., 1999; Ogilvie & Lin, 2004; Rieutord, 2009).
In Figure 3, we show solutions at two different forcing frequencies for the
case shown in Figure 2(a). We observe inertial wave beams in the convective
envelope that propagate with the angle $\theta_{i}$ between the rotation axis
along $z$ and the wavevector $\boldsymbol{\hat{k}}$ as predicted by
$\cos{\theta_{i}}=\omega/(2\Omega)$, thus satisfying the dispersion relation
$\omega^{2}=4\Omega^{2}\cos\theta_{i}$, with the group velocity (along which
energy travels) propagating perpendicular to this direction and lying along
the visible wave beams.

(a) $\omega=0.514$

(b) $\omega=0.95$

Figure 4: Illustrative examples of $|\boldsymbol{u}|$ for the forced response
within a planet with a stably stratified core and convective envelope, with
$\alpha=0.1$, $\beta=0.5$, $\bar{N}^{2}=1$, $\Omega=0.4$,
$\nu=\kappa=10^{-6}$. Panel 4(a) shows a gravito-inertial wave response in the
stratified layer and inertial wave response in the convective envelope. Panel
4(b) shows gravito-inertial wave in the stratified layer only.

In Figure 4(a) we consider the spatial structure of the forced solution at a
frequency of $\omega=0.514$ for the case in Figure 2(c), a frequency which is
within the range of both gravito-inertial waves in the stratified region and
inertial waves in the convective region. We observe both these waves to be
excited, the gravito-inertial waves within the stratified core, and the
inertial waves within the convective envelope. At this point we note the
similarity between the spatial structure in the outer envelope with the
example shown in Figure 3(b), where the forcing frequency is similar. When we
consider a higher frequency $\omega=0.95$, outside the inertial wave range but
within the gravito-inertial wave range, shown in Figure 4(b), we can see that
only gravito-inertial waves within the stably stratified layer are excited and
the response is evanescent above. This demonstrates how in a region with
stable stratification, instead of the straight lines characteristic of an
inertial wave we observe the curved lines characteristic of gravito-inertial
waves (e.g. Dintrans et al., 1999, though their profile has $N^{2}\propto r$).

One significant difference between the tidal response in rotating and non-
rotating models is that rotation breaks the symmetry between positive and
negative forcing frequencies. Figure 2(b) shows the dissipation rate for
negative forcing frequencies in a convective envelope with a solid core for
the same model as Figure 2(a), demonstrating that $D(\omega)\neq D(-\omega)$
in general. At small negative frequencies we have additional Rossby mode
resonances (also known as planetary waves). Rossby modes are a subset of
inertial modes obtained by considering conservation of vorticity, and in the
absence of a background flow and stratification, have a dispersion relation
(Papaloizou & Pringle, 1978; Zaqarashvili et al., 2021),

$\omega=-\frac{2m\Omega}{l(l+1)},$ (37)

which has strictly the opposite sign to rotation, indicating that they
propagate in the retrograde direction. These are more evident in thin shells
(large $\alpha$) and tend to excite modes with $l=3$ and $m=2$, for which
$\omega=-1/3$. At higher frequency magnitudes, the direction of the shift in
frequency of the surface gravity mode discussed above depends on the sign of
$\omega$, as predicted in Appendix A, and which we confirm numerically here.

##  4 Variation of the parameters

We now begin to explore how the dissipation depends on the parameters of our
model. To do so, we follow paper 1 and employ a frequency-averaged measure of
the dissipation, but now using the two different weightings introduced in §
2.5. In many of our parameter ranges, the highest frequency gravito-inertial
and surface gravity modes have comparable frequencies, therefore it is not
always possible to separate their behaviour using these integrated measures.
Without rotation, we can always take an upper limit on our integrals to be
$\omega_{max}=\bar{N}$ to include internal gravity modes and exclude surface
gravity modes, but we are unable to do the same here because gravito-inertial
waves propagate up to $\omega_{max}=\sqrt{\bar{N}^{2}+4\Omega^{2}}$, which
overlaps with surface gravity modes (which are themselves shifted to lower
frequencies for positive $\omega$). The frequency-averaged measures are
helpful to explore how the dissipative properties depend on the various
parameters, but we emphasise that results involving these quantities by
themselves should be interpreted with caution. This is partly because they use
a single number to represent an entire complicated spectrum, and partly
because this measure is less robust to the frequency interval considered than
in the non-rotating cases explored in paper 1.

###  4.1 Variation of rotation rate $\Omega$

First, we vary $\Omega$ to determine how the planetary rotation rate affects
its dissipative properties. We focus on two models: one with an entirely
stably stratified planet $\beta=1.0$ and one with an extended stably
stratified “dilute core” extending to half the planetary radius $\beta=0.5$,
with a well-mixed convective region above.

(a) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$

(b) $\bar{D}=\int D/\omega^{2}~{}\mathrm{d}\omega$

(c) $\bar{D}=(1/\Delta\omega)\int D~{}\mathrm{d}\omega$

Figure 5: Comparison of different frequency weightings for the frequency-
averaged dissipation as a function of rotation rate. Other parameters are
fixed at $\alpha=0.1$, $\beta=0.5$, $\bar{N}=1$, $\nu=\kappa=10^{-6}$. In all
cases the integration limits have been taken to be $\omega_{min}=0$ and
$\omega_{max}=\sqrt{\bar{N}^{2}+4\Omega^{2}}$.

(a) $\beta=0.5$

(b) $\beta=1$

(c) $\beta=1$

Figure 6: Dissipation rate as a function of $\omega$ for varying rotation
rates $\Omega$, with $\alpha=0.1$, $\bar{N}=1$, $\nu=\kappa=10^{-6}$ in all
cases. Panel 6(a) shows a stratified layer extending to half the planetary
radius, panels 6(b) and 6(c) show a uniform layer extending to the planetary
radius, with 6(c) showing just the low frequency regime (travelling wave).

In Figure 5 we consider $\beta=0.5$, and first compare the outcomes of three
different definitions of the frequency-averaged dissipation. We consider
weightings with $\frac{1}{\omega}$, $\frac{1}{\omega^{2}}$ and a “mean
average”, which we define to be
$\bar{D}=\frac{1}{\omega_{max}-\omega_{min}}\int D~{}\mathrm{d}\omega$. In all
three cases we have kept the integration limits the same as the range for
gravito-inertial waves, i.e. $\omega_{min}=0$ and
$\omega_{max}=\sqrt{\bar{N}^{2}+4\Omega^{2}}$; the same quantity for the
corresponding non-rotating case has been plotted for reference in each case.
There is a slight quantitative difference between these three definitions but
no qualitative differences. The similarity in all three measures implies that
we can focus on only one of these for an overview of the dissipative response.

In all three cases the total dissipation is larger than obtained without
rotation and there is a large jump between $\Omega=0.2$ and $\Omega=0.4$.
Since $\omega_{max}=\sqrt{N^{2}+4\Omega^{2}}$ is not a robust limit to
separate the surface and internal modes, it is also instructive to consider
the frequency-dependent solution, which we plot in Figure 6 for all of these
cases. Considering the vertical lines in Figure 6, which mark the gravito-
inertial wave limit, we can see that for the rotation rates
$\Omega=0,0.1,0.2$, the surface gravity mode is well above $\omega_{max}$,
whereas for $\Omega=0.4,0.6$, it is below this limit. In fact, we see that in
the case of $\Omega=0.6$, the surface gravity mode no longer appears as an
isolated peak. It is therefore likely that the surface gravity mode is to a
large part dictating the trends observed in Figure 5.

In Figure 7 we consider different integration limits, but show results only
for the $1/\omega^{2}$ frequency weighting. In the top panel, we consider
$\omega_{max}=1$, which removes the surface gravity mode from most cases (for
$\Omega=0.6$ it is still included). The frequency-averaged quantities increase
away from the non-rotating baseline as rotation rate increases and the
inertial wave response is enhanced. In the bottom panel, we consider a higher
integration limit $\omega_{max}=1.8$, which incorporates the surface gravity
mode behaviour in all cases. We observe that more rapid rotation leads to more
efficient dissipation according to this measure. Using the higher limit, we
find that there are some differences between the different weightings (not
shown) and those that give the most emphasis to inertial waves exhibit an
increase as rotation rate is increased.

(a) $\bar{D}=\int D/\omega^{2}~{}\mathrm{d}\omega$, $\omega_{min}=0$,
$\omega_{max}=1$

(b) $\bar{D}=\int D/\omega^{2}~{}\mathrm{d}\omega$, $\omega_{min}=0$,
$\omega_{max}=1.8$

Figure 7: Comparison of different integration limits for frequency-averaged
dissipation with $1/\omega^{2}$ weighting as a function of rotation rate.
Other parameters are fixed at $\alpha=0.1$, $\beta=0.5$, $\bar{N}=1$,
$\nu=\kappa=10^{-6}$.

We now turn to explore a case with a fully stratified planet with $\beta=1$ in
Figure 8, with all other parameters fixed at $\alpha=0.1$, $\bar{N}=1$,
$\nu=\kappa=10^{-6}$. Due to the similarities between the weightings, we show
here just the frequency-averaged measure used in the non-rotating case in
paper 1 ($1/\omega$), but we have included all three limits of integration
discussed so far in this section. We see that trends observed for $\beta=0.5$
also hold here. The only key difference is that, for data points in which the
surface gravity mode is not included in the integration, the overall
dissipation is lower than for $\beta=0.5$, due to the absence of a convective
envelope permitting inertial mode excitation.

(a) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$, $\omega_{min}=0$,
$\omega_{max}=1$

(b) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$, $\omega_{min}=0$,
$\omega_{max}=\sqrt{N^{2}+4\Omega^{2}}$

(c) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$, $\omega_{min}=0$,
$\omega_{max}=1.8$

Figure 8: Comparison of integration limits for the frequency-averaged
dissipation as a function of rotation rate for a fully stratified planet.
Other parameters at fixed at $\alpha=0.1$, $\beta=1$, $\bar{N}=1$,
$\nu=\kappa=10^{-6}$.

(a) $\Omega=0.2$, $\omega=0.25$

(b) $\Omega=0.2$, $\omega=0.95$

(c) $\Omega=0.2$, $\omega=1.1$

(d) $\Omega=0.4$, $\omega=0.25$

(e) $\Omega=0.4$, $\omega=0.95$

(f) $\Omega=0.4$, $\omega=1.1$

(g) $\Omega=0.6$, $\omega=0.25$

(h) $\Omega=0.6$, $\omega=0.95$

(i) $\Omega=0.6$, $\omega=1.1$

Figure 9: Examples of the spatial structure for different rotation rates and
forcing frequencies; in all cases we consider a stably stratified layer
extending to half the planetary radius $\alpha=0.1$, $\beta=0.5$, $\bar{N}=1$,
$\nu=\kappa=10^{-6}$.

Returning to the frequency-dependent dissipation in Figures 6(a) and 6(b),
which show the total dissipation ($D_{total}$) for four different rotation
rates, as well as the corresponding non-rotating cases from paper 1 (black-
dashed line); the cases shown are for $\beta=1$ and $\beta=0.5$, respectively.
We can see that there is a non-trivial balance between the roles of the
buoyancy and Coriolis forces. In all rotating cases we see gravito-inertial
waves in the expected range $\omega<\sqrt{\bar{N}^{2}+4\Omega^{2}}$. As well
as observing the increasing range of modes we also see that the irregular
pattern characteristic of inertial modes for $\beta=0.5$ is more pronounced as
rotation rate increases. Indeed, at a low rotation rate, $\Omega=0.1$, the
profile is similar to that of $\Omega=0$, as buoyancy forces appear to
dominate. When comparing the two figures, we note that for a partially
stratified planet, varying rotation has a more pronounced effect. This is
because the stratified layer acts similarly to a large solid core in enhancing
the excitation of inertial waves, which are not excited in a homogeneous full
sphere (e.g. Ogilvie, 2009).

Figure 6(c) shows the low frequency range for $\beta=1$. Based on our analysis
in paper 1 we interpret this regime as being the travelling wave regime, for
which tidally-excited inwardly-propagating gravito-inertial waves are excited
in the stably stratified core by perturbations (by non-wavelike tides plus
inertial waves) at $r\sim\beta r_{0}$, which are then subsequently fully
damped by viscous and thermal diffusion before they reflect from the solid
core and return to their launching sites. We can see that at low frequencies
there is a clear $O(1)$ dependence on $\Omega$, though it does not typically
lead to order of magnitude variations. It is possible that a similar analysis
to that carried out in the travelling wave regime in paper 1 (and considered
by Papaloizou & Savonije, 1997; Ogilvie & Lin, 2004; Ivanov et al., 2013, in
different models) could explain this dependence (e.g. adopting the
“traditional approximation”) but we leave this to future work as it only leads
to $O(1)$ differences in the dissipation.

Turning our attention to the spatial structure of the velocity magnitude in
the $\beta=0.5$ case in Figure 9, which shows the forced solution for four
values of $\Omega$ increasing in each row, with increasing forcing frequency
in the different columns. We see clear agreement in properties with the
predicted ranges for inertial and gravito-inertial waves in the convective
envelope and outer core, respectively, as we move between the different
regimes. Additionally we observe larger magnitude inertial waves as the
rotation rate $\Omega$ increases, corresponding to the larger dissipation
rates. Note the changing colour bars between the different cases. Solutions in
the first column are visually similar to the non-rotating solutions in
Fig.3(a) of paper 1, except that inertial waves are excited in the convective
envelope. In the middle column, showing a frequency outside the inertial range
in the envelope in all panels except for $\Omega=0.6$, the solution varies
substantially for the different rotation rates, from g-mode like in the top
panel, to being strongly modified by rotation in the core in the bottom two
panels. The rightmost column shows a higher frequency for which the free
surface perturbations become more important, for which $\omega=1.1$ is outside
the inertial range except for $\Omega=0.6$. Panel (l) shows a mode with a
complicated structure that exhibits properties of each of surface gravity,
inertial and gravito-inertial modes near the surface, in the envelope and
core, respectively.

Varying the rotation rate can therefore modify the modes and tidal dissipation
rates substantially, depending on the relevant tidal frequency, and whether it
lies relative to the free inertial, gravito-inertial and surface gravity
modes. Overall, the dissipation is typically enhanced for faster rotation
rates.

###  4.2 Variation of core sizes $\alpha$ and $\beta$ for uniform $N$

A large uncertainty is how deep any stable layers extend throughout planetary
interiors, with current estimations and observational constraints varying
significantly. We know from previous studies that the size of a solid core can
significantly enhance the excitation of inertial waves in incompressible
models (Ogilvie, 2009; Goodman & Lackner, 2009; Rieutord & Valdettaro, 2010).
Here we compare two scenarios, cases with a solid core with a convective
envelope above ($\bar{N}=0$, varying $\alpha$) and cases where a small solid
core is surrounded by a stably stratified (“dilute core”) layer again beneath
a convective envelope ($\bar{N}=1$, $\alpha=0.1$, varying $\beta$).

(a) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$, $\omega_{min}=0$,
$\omega_{max}=1$ Figure 10: Frequency-averaged dissipation as a function of
core size for both cases where a convective layer sits above a solid core
($\bar{N}=0$ and varying $\alpha$) and those where a stably stratified layer
extends to the same radius ($\bar{N}=1$, $\alpha=0.1$, varying $\beta$). Other
parameters kept constant at $\Omega=0.4$ and $\nu=\kappa=10^{-6}$.

(a) $\bar{N}=0$ and varying $\alpha$

(b) $\bar{N}=1$, $\alpha=0.1$, varying $\beta$

(c) Core extending to $0.25r_{0}$

(d) Core extending to $0.55r_{0}$

(e) Core extending to $0.85r_{0}$

Figure 11: Frequency dependence of the dissipation rate for different core
sizes, comparing both a solid core and stably stratified core, with
$\Omega=0.4$, $\nu=\kappa=10^{-6}$ in all cases. Panels 11(a) and 11(b) show a
solid core and stratified layer, respectively, for different core sizes.
Panels 11(c) to 11(e) compare a solid core to a stratified layer for three
different radii.

Figure 10 shows the frequency-averaged profiles for these two types of model.
The darker lines show neutrally stratified cases and the dependence on the
radius of the solid core $\alpha$, whilst the lighter coloured lines show the
cases varying the size of the stably stratified core $\beta$. The black line
shows the equivalent trend found in the non-rotating cases, where we consider
a small core with a stably stratified region above. We can see that, in both
cases, the frequency-averaged dissipation shows a strong dependence on core
size, but this trend varies between the two core types. For core sizes less
than $0.5r_{0}$ there are larger dissipation rates where there is a stably
stratified core; for larger core sizes (larger than would be expected in the
case of Jupiter and Saturn), we see the opposite trend. This suggests that for
the smaller core sizes, the stratified layer acts as a solid core for the
excitation of inertial waves, with additional contribution to the overall
dissipation rate arising from the excitation of gravito-inertial waves in the
stratified layer. We consider the trend found for the largest core sizes with
caution however, as when considering the frequency-dependent dissipation in
Figure 11, we see that the surface gravity mode has shifted to frequencies
less than the integration limit $\omega=1$ when considering a solid core,
which also contributes to the dissipation rate. We again see the enhancement
of the dissipation in solid core cases with rotation compared with non-
rotating results.

Figure 11 shows the frequency-dependent dissipation rates, which allow us to
understand more about the contribution of a stratified core. In Figures 11(a)
and 11(b) we compare cases with differently-sized cores for a solid core and a
stratified core, respectively. In both cases, we see that the amplitudes of
the resonances increase as the core sizes increase, but we note some
differences. For cases with a stratified core, we observe the increased
frequency range of resonances matching gravito-inertial waves rather than
inertial waves. Additionally, when a stably stratified layer extends to
$0.85r_{0}$, there are regular discrete peaks that are characteristic of
internal gravity waves that dominate around a forcing frequency of $0.2$. In
Figures 11(c) to 11(e), we compare the total dissipation rates for cases with
a solid core and stably stratified core for three different radii, $0.25$,
$0.55$ and $0.85$. We see clearly that for the smallest core size, the
inertial wave response is very similar in both models, departing from each
other at frequencies around $0.6-0.9$ only. However, as we increase the core
size we increase the contribution of the gravito-inertial waves so this
agreement between a solid core model and a dilute core is worse. Indeed, when
the core extends as far as $0.85r_{0}$, we see at frequencies less than $0.6$
that the dissipation in a stably stratified core resembles the equivalent non-
rotating system more than the equivalent solid core case. We note also that
the shift in the frequency of the surface gravity mode varies between these
cases as the stratified layer also affects the frequency of this mode (as
predicted analytically without rotation in paper 1). Finally, we note that
when considering these thin shell convective regions (large cores) we expect
negative frequency Rossby modes to be more significant than the positive
frequency modes we show here, which require further investigation to explore.

###  4.3 Variation of step number for a density staircase

(a) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$, $\omega_{min}=0.1$,
$\omega_{max}=1$

(b) $\bar{D}=\int D/\omega~{}\mathrm{d}\omega$, $\omega_{min}=0.1$,
$\omega_{max}=1$

Figure 12: Frequency-averaged dissipation as a function of the number of
steps. Other parameters kept constant at $\alpha=0.1$, $\beta=1.0$,
$\bar{N}=1$, $\Omega=0.4$, $\nu=\kappa=10^{-6}$. Panel 12(a) compares total,
viscous and thermal dissipation with the total dissipation of the uniformly
stratified case and the non-rotating case. Panel 12(b) compares two different
viscosities/diffusivities for the same case.

(a) $\nu=\kappa=10^{-7}$

(b) $\nu=\kappa=10^{-6}$

(c) $\nu=\kappa=10^{-4}$

Figure 13: Frequency dependence of dissipation for different numbers of steps
in the staircase density profile. In all cases $\alpha=0.1$, $\beta=1$,
$\bar{N}=1$. In panel 13(b), we fix $\nu=\kappa=10^{-6}$ and in panel 13(c) we
fix $\nu=\kappa=10^{-4}$.

(a) steps $=0$

(b) steps $=1$

(c) steps $=5$

(d) steps $=9$

(e) steps $=13$

Figure 14: Examples of the spatial structure of the forced response for a
density staircase with 0 (i.e. unstratified),1,5,9 and 13 steps, in all cases
with $\alpha=0.1$, $\beta=1.0$, $\bar{N}=1$, $\nu=\kappa=10^{-6}$, and a
forcing frequency $\omega=0.25$.

We now turn to explore the consequences of semi-convective layers. As in non-
rotating cases in paper 1, we consider a staircase-like density profile as
defined by equation 30, and we consider cases in which $\alpha=0.1$,
$\beta=1.0$, $\bar{N}=1$, $\delta r=0.03$, $\Omega=0.4$ and
$\nu=\kappa=10^{-6}$, with varying step numbers. We consider $\beta=1.0$ here
to allow us to explore the largest number of steps possible, allowing us to
resolve modes with variations on the scale of the steps while minimising
computational costs. In Figure 12(a) we show the frequency-averaged
dissipation as the step number is increased; we have in this case shown just
the $\frac{1}{\omega}$ weighting and used an upper integration limit of
$\omega_{max}=1$ to exclude the surface gravity waves. Note, the lower limit
of $\omega_{min}=0.1$ has been chosen for purely numerical reasons; this is
because the lowest frequency waves require the highest resolution and
contribute little to the overall dissipation using this measure. For
comparison, we also show the equivalent non-rotating case as a function of
step number (black-dashed line), as well as the uniformly stratified layer
with equivalent mean stratification (black solid line).

We see that the dissipation for all step numbers considered is higher than
both a uniformly stratified rotating medium and cases without rotation
according to this measure. Although there is a trend towards the uniformly
stratified case as we increase the number of steps, our results do not
converge to this result as quickly as we had found in the non-rotating cases
in Figure 11 of paper 1. The dissipation is initially nearly 8 times larger
with one step as the case with a constant $N$. For large step numbers though,
there is a trend towards the uniformly stratified case, implying that if the
planetary interior contains a large number of steps it will behave on average
similar to an equivalent continuously stratified medium. Here the staircase
dissipation remains slightly larger than the uniform case even for 13 steps,
though it is very similar.

In Figure 12(b) we explore the effects of varying viscosity and thermal
diffusivity according to the same frequency-averaged measure, showing total
dissipation for $\nu=\kappa=10^{-6}$ compared with $\nu=\kappa=10^{-4}$.
Although it is hard to draw robust conclusions from two values of the
diffusivities, these initial results suggest that our findings may be robust
to varying the viscosity and thermal diffusivity to approach planetary values.
We also expect that varying the Prandtl number may also be unimportant for the
total dissipation, and that, as found in paper 1, it will primarily alter the
balance between viscous and thermal dissipation rates rather than the total
dissipation.

Figure 13 shows the corresponding frequency-dependent dissipation rates for
various cases with different step numbers for $\nu=\kappa=10^{-6}$ in panel
(a) and $\nu=\kappa=10^{-4}$ in (b). At low frequencies the behaviour varies
significantly for different numbers of steps. Note that the frequency limit
for propagation of gravito-inertial waves is
$\omega<\sqrt{N_{max}^{2}+4\Omega^{2}}$ and when considering a staircase-like
structure $N_{max}\geq\bar{N}$, thereby increasing the range slightly (but
with no significant implications for the most astrophysically relevant
parameter values, so we do not focus on this aspect).

We first compare the cases with zero steps, one step and five steps, noting
that zero steps consists of a small solid core with a convective envelope to
the outer edge. At low frequencies, although inertial waves are excited in the
convective envelopes in all three cases, these are barely visible for the case
of zero steps, where there is only the small solid core to launch inertial
waves from. As we increase the number of steps to one and then five, the
excitation of inertial waves occurs from boundaries at additional and
increasingly larger radii, leading to enhanced dissipation. Following this, we
now compare the cases with five, nine and 13 steps and find there are
significantly smaller differences between these cases. At this point the outer
radius does not vary significantly as the number of steps is increased. The
individual resonant peaks do shift as the number of steps is increased, but
overall the dissipation spectrum looks similar when there are “enough steps”
(here this means more than five).

Figure 14 shows the forced solutions at $\omega=0.25$ for all five step
numbers explored for $\nu=\kappa=10^{-6}$. We see inertial wave beams excited
in the convective layers, which have a very small amplitude response dominated
by the non-wavelike tide in the case of zero steps, but increase in amplitude
as we increase the number of steps and the outer interface moves outwards. The
case with one step behaves visually like a solid core at the location of the
interface, with inertial waves in the envelope and little activity inside the
interface within the “core”. We will explore this similarity further in the
next section. As we increase the number of steps, the solution appears to be
dominated by inertial waves, presumably excited at the critical latitudes on
each interface (where rays are tangent to it), that subsequently bounce
between the interfaces to form these modes. The solution becomes increasingly
dominated by wavelike tides, rather than the larger scale non-wavelike tide,
as we increase the number of steps in these examples.

Figure 15: Comparison of total dissipation between profiles with interfaces,
uniformly stratified layers and a solid core, with $\Omega=0.4$ and
$\nu=\kappa=10^{-6}$ in all cases. The black solid lines are the single and
five step cases in panels 15 and 15 respectively. The dashed coloured lines
describe a stratified layer with $\alpha=0.1$ and $\bar{N}=1$, and solid
coloured lines a solid core with $\bar{N}=0$. The blue aligns the
core/stratification with the first (or only) interface and the red the last
interface.

We are beginning to see that the size of the convective envelope is key to the
behaviour of the inertial waves and the corresponding dissipation, which is
consistent with expectation from incompressible models of the unstratified
case (e.g. Ogilvie, 2009; Goodman & Lackner, 2009; Rieutord & Valdettaro,
2010; Ogilvie, 2013). Therefore, we now compare cases with a solid or
stratified core that extends to the same radii as the staircase interfaces.
The solid black lines on Figures 15 and 15 show the total dissipation as a
function of frequency for a single step and five steps, respectively. In
Figure 15, we have compared this case to both a solid core, and a stratified
layer extending to that radius, i.e. $\alpha=0.55$ and $\bar{N}=0$ and
$\alpha=0.1$, $\beta=0.55$ and $\bar{N}=1$. We can see that at low frequencies
there is very good agreement between these three different profiles,
suggesting that the forced wave response in the outer envelope only depends
weakly on what is below the envelope, provided the buoyancy frequency of the
stably stratified layers is sufficiently strong. We will explore this issue
further in the next section.

In Figure 15 we have similarly shown a solid core and a stably stratified
layer, this time extending to the first interface (blue, $\alpha/\beta=0.25$)
and last interface (red, $\alpha/\beta=0.85$). We notice that although the
agreement is not as good as in the single step case, there is still closer
agreement between the cases where the core corresponds to the last interface
of the staircase. This suggests that this is the key interface in dictating
the dissipation due to inertial waves. We note that dissipation is larger in
the staircase model due to the additional interior interfaces where inertial
waves can be excited from critical latitudes. In both figures there is
significantly different behaviour between the various models in the mid to
high frequencies where it is expected that the buoyancy effects dominate. This
comparison suggests that the importance of buoyancy forces compared with
Coriolis forces will strongly depend on the forcing frequency.

##  5 Comparison of dilute core models

Motivated by the relevance of giant planet interior models that consist of an
extended “dilute” core, and that the properties of this core are highly
uncertain, in this section we explore further the consequences of different
buoyancy frequency profiles describing stable stratification. We consider four
cases that have different buoyancy profiles that each represent a stratified
outer core extending to half of the planetary radius surrounded by a
convective envelope. We compare:

  * •

Case 1: a large solid core with a convective envelope, $\alpha=0.5$,
$\bar{N}=0$.

  * •

Case 2: a uniformly stably stratified layer extending from an inner core
boundary to an outer core boundary, $\alpha=0.1$, $\beta=0.5$, $\bar{N}=1$.

  * •

Case 3: a single stable interface at the outer core boundary $\beta r_{0}$,
$\alpha=0.1$, $\beta=0.5$, $\bar{N}=1$.

  * •

Case 4: a staircase extending from an inner core boundary to an outer core
boundary, $\alpha=0.1$, $\beta=0.5$, $\bar{N}=1$, steps $=3$.

All other parameters are kept constant with $\Omega=0.4$ and
$\nu=\kappa=10^{-6}$. We plot the total dissipation in each case in Figure 16.

We find that when considering the low frequency inertial range (i.e. observed
to be $|\omega|\lesssim\Omega=0.4$), the frequency-dependent dissipation is
remarkably similar in all four cases. The inertial wave behaviour in the
convective envelope then appears to dominate the behaviour and is little
affected by the form of the stratified layer (or solid core) beneath it. We
see that the stable layer, single interface, and staircase, each act like a
solid boundary for the propagation of inertial waves in the convective
envelope, and they enhance the dissipation for low frequencies in a very
similar manner.

Figure 16: Total dissipation for different dilute core modes, where in all
cases $\Omega=0.4$, $\nu=\kappa=10^{-6}$.

Case 1

Case 2

Case 3

Case 4

(a) $\omega=0.21$

(b) $\omega=0.73$

(c) $\omega=0.90$

(d) $\omega=0.21$

(e) $\omega=0.73$

(f) $\omega=0.73$

(g) $\omega=0.21$

(h) $\omega=0.73$

(i) $\omega=0.90$

(j) $\omega=0.21$

(k) $\omega=0.73$

(l) $\omega=0.90$

Figure 17: Examples of the spatial structure for different dilute cores and
forcing frequencies, with $\nu=\kappa=10^{-6}$ and $\Omega=0.4$.

Within the mid-frequency range between $0.8\lesssim\omega\lesssim 1$, the
behaviour varies significantly, as it is in this frequency range that the
gravito-inertial modes within the stratified region are dominant when
$\bar{N}=1$. These modes are sensitive to the form of the stratification
adopted, and we observe peaks corresponding to gravito-inertial modes in the
case of a uniformly stratified layer, as well as the interfacial modes that
are characteristic of a staircase structure in cases with interfaces.

In Figure 17, we compare the spatial structure of the response in all four
cases at three different forcing frequencies, using colour-scales that differ
in each panel to most clearly illustrate the variation in space. Considering
the first column for which the forcing frequency is low ($\omega=0.21$),
within the inertial wave range, we observe the solution in the convective
envelope to be similar in each case. Their structures are very similar, and
their amplitudes are mostly similar but do differ to some extent. In all
examples the stable stratification is acting effectively as a solid boundary
for the propagation of inertial waves, showing that an extended (sufficiently
stably) stratified core acts like a large solid core, enhancing the
dissipation over cases with a small core. Inertial waves are presumably
excited in a similar way at the critical latitudes on the interface at the
outer boundary of the core.

In the second column, at a higher forcing frequency of $\omega=0.73$, the
spatial structure observed in the convective region again remains consistent.
However, we now see the different modes that form within the stratified
region, which vary significantly, contributing to the differences between the
resultant dissipation. At both forcing frequencies, for the staircase and
interface cases we faintly observe additional inertial modes in the deeper
convective layers, as well as gravito-inertial modes in the uniformly
stratified case. These additional modes can explain the increase in
dissipation observed for Cases 2–4 over that of a solid core (Case 1).

Finally, in the last column we show the solution with the highest forcing
frequency $\omega=0.9$, at which we are outside the inertial wave range
($|\omega|<2\Omega=0.8$) but within the gravito-inertial wave ranges. We
clearly observe different responses in each case that depend on the properties
of the stratified region, with wavelike behaviour only observable in the cases
with a stratified layer and a staircase structure. Purely non-wavelike
behaviour is observed in the unstratified case, as expected since inertial
waves do not propagate with this frequency.

To summarise this section and § 4.3, we find that a sufficiently strongly
stratified core behaves very similarly to a rigid core for the inertial wave
response, and for the corresponding dissipation at low frequencies
($|\omega|\ll\bar{N}$), and it is largely insensitive to the properties of the
stratified layer beneath. For larger frequencies
$\Omega\lesssim|\omega|\sim\bar{N}\sim\omega_{d}$, we find more substantial
differences between these cases due to the decreasing importance of rotation
relative to internal and surface buoyancy forces, and the transmission of wave
energy between the convective and stably stratified layers. For the forcing
due to Jupiter’s and Saturn’s moons, we typically expect $|\omega|\sim\Omega$,
but it is not certain how large $\bar{N}$ should be. Weaker values of
$\bar{N}\sim\Omega$, such as those considered by Lin (2023); Dewberry (2023)
for Jupiter would permit substantial connection between convective and stable
layers for these tidal frequencies, such that the response would be expected
to differ somewhat from the case of a solid core. On the other hand the larger
values inferred for Saturn (Mankovich & Fuller, 2021) for which
$\bar{N}\gg\Omega$, would typically predict that the inertial wave response in
the envelope may be better represented by a solid core at the outer interface
of the stratified “dilute core”. Similarly, extrapolating our result to slowly
rotating solar-type stars, for which $\bar{N}\gg\Omega$, the radiative zone
would be expected to behave quite similarly to a rigid core for the inertial
waves in the envelope, thereby motivating use of rigid boundary conditions in
studies of their properties (e.g. Astoul & Barker, 2022).

##  6 Application to Saturn’s tidal dissipation

Figure 18: Example of total dissipation and modified tidal quality factor for
Saturn-like parameters, with the tidal forcing frequency of six of Saturn’s
moons over-plotted: Mimas, Enceladus, Tethys, Dione, Rhea and Titan.

Finally, we explore cases with parameter values that are as consistent as
possible with the latest observational constraints for Saturn. Our model is
the simplest global one able to capture the dynamics of stratified layers,
since we adopt the Boussinesq approximation, and hence neglect realistic
variations in density throughout the planet. It is still instructive however
to compute the tidal response in this model and to compare with observational
constraints. To do this, we adopt values similar to those suggested by
Mankovich & Fuller (2021) and consider an example where $\alpha=0.1$ (tiny
solid core for which the precise value of $\alpha$ is unlikely to be
important), $\beta=0.6$ (extended stably stratified “dilute core”),
$\bar{N}=2$ and $\Omega=0.4$. We consider both a uniformly stratified case as
well as a staircase structure with one, five and nine steps. Figure 18 shows
the dissipation rate using these four profiles and Figure 18 shows the
modified tidal quality factor $Q^{\prime}$ for the case of a uniformly
stratified layer, where $Q^{\prime}=15|\omega|/(16\pi D_{total})$ – see
section 4 of paper 1 for further details. The black vertical lines show the
tidal frequency of six of Saturn’s major moons: Mimas, Enceladus, Tethys,
Dione, Rhea and Titan, as points of reference for the relevant frequency
regimes. The tidal forcing frequencies are also reported in Table
LABEL:tab:Saturn_moons.

We see that all four cases show qualitatively and quantitatively similar
dissipation profiles, and the typical level of dissipation (e.g. quantified by
the frequency average) is almost unchanged as the step number is varied. In
these examples all three wave frequency ranges – gravito-inertial, inertial,
and surface gravity modes – overlap, making it difficult to separate the
behaviour of each. However, given the sensitivity to the tidal frequency due
to the moons exhibited, stably stratified layers could have important
implications for Saturn’s tidal dissipation rates.

The tidal quality factor $Q^{\prime}$ obtained in Figure 18 ranges from
approximately $10^{2}$ to $10^{4}$ at the frequencies relevant for Saturn’s
moons. This is comparable to observational constraints from the migration
rates of Saturn’s moons (e.g. Lainey et al., 2017), which provide
$Q^{\prime}\approx(0.94\pm 0.44)\times 10^{4}$. Our idealised calculations
therefore highlight the importance of considering stably-stratified layers on
the excitation and dissipation of inertial and internal waves in planets. This
figure demonstrates that efficient tidal dissipation rates – sufficient to
explain the observed migration rates of Saturn’s moons – are predicted at the
frequencies of the orbiting moons due to the excitation of inertial and
gravito-inertial waves in our models with stable layers. Note that we do not
require resonance locking to operate (e.g. Fuller et al., 2016), though the
possibility of resonance locking and whether or not it is a viable mechanism
in giant planets should be explored further. The presence of a stably
stratified dilute core, and its effects in enhancing inertial wave excitation
in the overlying convective envelope, as well as – to a lesser extent – the
additional gravito-inertial mode excitation in the fluid core itself, may be
the key mechanisms of tidal dissipation in Saturn, and which could explain
observations. We envisage similar results may apply to Jupiter also (see also
Lin, 2023; Dewberry, 2023), though see the caveats in § 5.

##  7 Conclusions

We have presented new theoretical models of giant planets containing stable
layers similar to those constrained observationally for Saturn, and
hypothesised for Jupiter, to explore dissipation of tidal flows inside these
planets. We have studied the role of stably-stratified and semi-convective
layers on tidal dissipation in rotating giant planets, extending our prior
work without rotation (in paper 1, Pontin et al., 2023) to account for
Coriolis forces. Rotation permits the propagation of inertial waves that can
significantly enhance tidal dissipation in neutrally-stratified convective
regions. These can be tidally forced for frequencies $|\omega|\leq 2|\Omega|$,
which is typically the most relevant range for solar and extrasolar giant
planets. Rotation also modifies the properties of internal, surface and
interfacial gravity modes (the former are then commonly referred to as
gravito-inertial waves). With our idealised (Boussinesq) model of a rotating
and tidally-forced planet, we analysed the dissipative fluid response in a
spherical shell using both linear theoretical analysis and numerical
calculations (with high-resolution spectral methods). Our parameter study
analysed the properties of both the dissipative forced response and free
oscillation modes as we varied the properties of any stably-stratified layers
in the planet (including their sizes, strengths and compositions – layered or
smooth), the sizes of any solid core, the viscosity and thermal diffusivity
and the rotation rate (relative to the dynamical frequency), in addition to
scanning the full range of relevant tidal frequencies.

We found the presence of an extended stably-stratified fluid core in a giant
planet significantly enhanced tidal wave excitation of both inertial waves in
the convective envelope and gravito-inertial waves in the core. We have
demonstrated that efficient tidal dissipation rates – sufficient to explain
the observed migration rates of Saturn’s moons – are predicted at the
frequencies of the orbiting moons due to the excitation of inertial waves in
convective envelopes in our models with interior stable layers, and to a
lesser extent gravito-inertial waves in the fluid core itself. Stable layers
could also be important for tidal evolution of hot and warm Jupiters, and hot
Neptunes, providing efficient tidal circularisation rates (see discussion in
Pontin et al., 2023).

We analysed both the frequency-dependent and frequency-averaged response to
establish some overarching trends in this problem, building upon our non-
rotating study in Pontin et al. (2023). We established that increasing the
rotation rate typically enhances the inertial wave response, in turn
increasing the total dissipation. Increasing the size of the core, whether it
is a solid core or a stably stratified layer (layered or smooth),
significantly increases the dissipation rate. We find that gravito-inertial
waves excited in a stably stratified layer can enhance the dissipation
compared to that of a solid core with the same radius, depending on the tidal
frequency that is relevant. As in non-rotating cases, we established that
provided a sufficient number of steps are present in a staircase-like density
structure, the region will behave like a uniformly stably stratified layer
when considering any frequency-integrated quantities. However, important
differences in the frequencies of the free modes, and hence the enhancement of
dissipation at the frequencies of these resonances, are found, which can lead
to significant differences in the response to different stratified models at a
given tidal frequency.

We found that a key parameter in the excitation of inertial waves is the size
of the outer convective envelope (i.e. the radius to which the dilute core
extends). The dissipative properties of the envelope were shown to be
approximately independent of the properties of the stratified layer beneath
it, whether it is a stably stratified layer (layered or smooth), single
interface or solid core for Saturn-like parameter values. The buoyancy
frequency profile beneath this layer can alter the frequencies of the free
modes however. This conclusion is valid for low frequencies relative to the
buoyancy frequency of the stable layer, as expected for Saturn. Different
results might be found for more weakly stratified layers.

Future work should study more sophisticated planetary models that also account
for magnetism and differential rotation (Baruteau & Rieutord, 2013; Guenel et
al., 2016; Wei, 2018; Lin & Ogilvie, 2018; Astoul et al., 2021; Astoul &
Barker, 2022), as well as the interaction of inertial waves with turbulent
convection (as opposed to large-scale tidal flows e.g., Duguid et al., 2020).
It would also be of interest to separate the thermal and compositional
contributions to the buoyancy and study whether double-diffusive effects could
be important in this problem.

CMP was supported by STFC PhD studentship 2024753. AJB was supported by STFC
grants ST/R00059X/1, ST/S000275/1 and ST/W000873/1. RH was supported by STFC
grants ST/S000275/1 and ST/W000873/1. We thank Stéphane Mathis and Quentin
André for discussions at an early stage in this project, and Aurélie Astoul,
Chris Jones and Gordon Ogilvie for helpful feedback. We would like to thank
the two reviewers for their constructive and positive referee reports.

##  Appendix A Analytical derivation of f-mode frequencies with rotation

Here we outline a derivation of the f-mode frequencies of a uniformly rotating
homogeneous spherical fluid body i.e. neglecting centrifugal deformations,
consistently with the model adopted throughout our paper. We follow Barker et
al. (2016) and use the elegant Lagrangian perturbation theory of Lebovitz
(1989a, b). For this section only we define
$\boldsymbol{\xi}\in\mathbb{R}^{3}$ to be the Lagrangian displacement vector,
satisfying (in the rotating frame)

$\displaystyle\partial_{t}^{2}\boldsymbol{\xi}+2\boldsymbol{\Omega}\times\partial_{t}\boldsymbol{\xi}-\nabla(\boldsymbol{\xi}\cdot\nabla
p)+\nabla\Delta p=0,$ (A1) $\displaystyle\nabla\cdot\boldsymbol{\xi}=0,$ (A2)

where $\Delta p$ is the Lagrangian pressure perturbation (which vanishes on
the free surface) and we neglect perturbations to the gravitational potential.
For our basic state, we assume a fixed spherically-symmetric gravitational
potential such that $\nabla p=-\omega_{d}^{2}\boldsymbol{r}$. We seek
solutions in the form of solenoidal vector fields whose components are
polynomials in the Cartesian coordinates up to a specified harmonic degree
$\ell=2$ (for the pressure perturbation), such that

$\boldsymbol{\xi}(\boldsymbol{x},t)=\sum_{i=1}^{i_{\mathrm{max}}}\alpha_{i}(t)\boldsymbol{\xi}_{i}(\boldsymbol{x}),$
(A3)

where $\alpha_{i}(t)\in\mathbb{C}$ is an amplitude and
$i_{\mathrm{max}}\in\mathbb{Z}^{+}$ is the number of basis elements
considered. To analyse surface gravity modes, we consider irrotational motions
that perturb the boundaries of the body (i.e. basis elements belonging to the
subspace $\mathrm{U}_{2}$ in Barker et al. 2016). To exactly represent all
surface gravity modes with $\ell=1$, we require $i_{\mathrm{max}}=3$ linearly-
independent vectors,

$\displaystyle\left(1,0,0\right)^{T},\quad\left(0,1,0\right)^{T},\quad\left(0,0,1\right)^{T}.$
(A4)

To exactly represent all such modes up to $\ell=2$, we require
$i_{\mathrm{max}}=8$, of which the first 3 are those in (A4), with 5
additional linearly independent vector fields

$\displaystyle\left(y,x,0\right)^{T},\quad\left(0,z,y\right)^{T},\quad\left(z,0,x\right)^{T},\quad\left(x,-y,0\right)^{T},\quad\left(x,y,-2z\right)^{T}.$
(A5)

We can project (A1) onto the basis given in (A3) using the inner product
(involving integration over volume $V$)

$\langle\boldsymbol{\xi}_{i},\boldsymbol{\xi}_{j}\rangle=\int_{V}\boldsymbol{\xi}^{*}_{i}\cdot\boldsymbol{\xi}_{j}\,\mathrm{d}V,$
(A6)

to obtain

$\langle\boldsymbol{\xi}_{i},\boldsymbol{\xi}_{j}\rangle\ddot{\alpha}_{j}+\langle\boldsymbol{\xi}_{i},2\boldsymbol{\Omega}\times\boldsymbol{\xi}_{j}\rangle\dot{\alpha}_{j}+\langle\boldsymbol{\xi}_{i},-\nabla(\boldsymbol{\xi}_{j}\cdot\nabla
p)\rangle\alpha_{j}=0,$ (A7)

where a sum over $j$ is implied. Seeking solutions
$\alpha_{j}\propto\mathrm{e}^{-\mathrm{i}\omega t}$, this problem is converted
to a quadratic eigenvalue problem of the form

$(-\omega^{2}\mathrm{M}-\mathrm{i}\mathrm{A}\omega+\mathrm{B})\boldsymbol{\alpha}=\boldsymbol{0},$
(A8)

for appropriate matrices M, A and B and column vector $\boldsymbol{\alpha}$
with elements $\alpha_{j}$. To obtain all $\ell=1$ surface gravity modes, we
restrict $i_{\mathrm{max}}=3$ to obtain the 6 eigenvalues

$\omega=\underbrace{\pm\Omega\pm\sqrt{\Omega^{2}+\omega_{d}^{2}}}_{m=1},\quad\mathrm{or}\quad\underbrace{\pm\omega_{d}}_{m=0},$
(A9)

where we have identified the corresponding azimuthal wavenumber magnitude $m$.
These are unimportant but strictly unphysical modes that involve the body
oscillating about a fixed position in space, and arise because we have fixed
the gravitational potential. On the other hand for $i_{\mathrm{max}}=8$ we
obtain the additional eigenvalues

$\omega=\underbrace{\pm\Omega\pm\sqrt{\Omega^{2}+2\omega_{d}^{2}}}_{m=2},\quad\mathrm{or}\quad\underbrace{\pm\frac{\Omega}{2}\pm\frac{1}{2}\sqrt{\Omega^{2}+8\omega_{d}^{2}}}_{m=1},\quad\mathrm{or}\quad\underbrace{\pm\sqrt{2}\omega_{d}}_{m=0}.$
(A10)

These are the frequencies of all surface gravity modes with $\ell=2$. When
$\Omega^{2}\ll\omega_{d}^{2}$, i.e. for slow rotation, we find the latter can
be written

$\omega=\pm\sqrt{2}\omega_{d}\pm\frac{m}{2}\Omega,$ (A11)

which agree with the results of Lebovitz (1961) for a (self-gravitating)
Maclaurin spheroid in the same limit. In this paper we have adopted a fixed
gravitational potential for simplicity, and because the effects of surface
perturbations on the gravitational potential are unrealistically enhanced in a
Maclaurin spheroid over a more realistic planet model that is denser at its
centre than near its surface (see Appendix A of Barker et al., 2016).

## References

  * André et al. (2017) André, Q., Barker, A., & Mathis, S. 2017, A&A, 605, A117, doi: https://doi.org/10.1051/0004-6361/201730765
  * André et al. (2019) André, Q., Mathis, S., & Barker, A. J. 2019, A&A, 626, A82, doi: 10.1051/0004-6361/201833674
  * Astoul & Barker (2022) Astoul, A., & Barker, A. J. 2022, MNRAS, 516, 2913, doi: 10.1093/mnras/stac2117
  * Astoul et al. (2021) Astoul, A., Park, J., Mathis, S., Baruteau, C., & Gallet, F. 2021, A&A, 647, A144, doi: 10.1051/0004-6361/202039148
  * Barker & Astoul (2021) Barker, A. J., & Astoul, A. A. V. 2021, MNRAS, 506, L69, doi: 10.1093/mnrasl/slab077
  * Barker et al. (2016) Barker, A. J., Braviner, H. J., & Ogilvie, G. I. 2016, Monthly Notices of the Royal Astronomical Society, 459, 924, doi: 10.1093/mnras/stw701
  * Baruteau & Rieutord (2013) Baruteau, C., & Rieutord, M. 2013, Journal of Fluid Mechanics, 719, 47, doi: 10.1017/jfm.2012.605
  * Bolton et al. (2017) Bolton, S. J., Adriani, A., Adumitroaie, V., et al. 2017, Science, 356, 821, doi: 10.1126/science.aal2108
  * Boyd (2001) Boyd, J. P. 2001, Chebyshev and Fourier Spectral Methods, 2nd edn. (Mineola, New York: Dover)
  * Braviner & Ogilvie (2014) Braviner, H. J., & Ogilvie, G. I. 2014, MNRAS, 441, 2321, doi: 10.1093/mnras/stu704
  * Braviner & Ogilvie (2015) —. 2015, MNRAS, 447, 1141, doi: 10.1093/mnras/stu2521
  * Connerney et al. (2022) Connerney, J. E. P., Timmins, S., Oliversen, R. J., et al. 2022, Journal of Geophysical Research (Planets), 127, e07055, doi: 10.1029/2021JE007055
  * de Vries et al. (2023) de Vries, N. B., Barker, A. J., & Hollerbach, R. 2023, MNRAS, 524, 2661, doi: 10.1093/mnras/stad1990
  * Dewberry (2023) Dewberry, J. W. 2023, arXiv e-prints, arXiv:2301.07097, doi: 10.48550/arXiv.2301.07097
  * Dewberry & Lai (2022) Dewberry, J. W., & Lai, D. 2022, The Astrophysical Journal, 925, 124, doi: 10.3847/1538-4357/ac3ede
  * Dintrans et al. (1999) Dintrans, B., Rieutord, M., & Valdettaro, L. 1999, Journal of Fluid Mechanics, 398, 271, doi: 10.1017/S0022112099006308
  * Duguid et al. (2020) Duguid, C. D., Barker, A. J., & Jones, C. A. 2020, Monthly Notices of the Royal Astronomical Society, 497, 3400, doi: 10.1093/mnras/staa2216
  * Durante et al. (2020) Durante, D., Parisi, M., Serra, D., et al. 2020, Geophys. Res. Lett., 47, e86572, doi: 10.1029/2019GL086572
  * French et al. (2012) French, M., Becker, A., Lorenzen, W., et al. 2012, ApJS, 202, 5, doi: 10.1088/0067-0049/202/1/5
  * Fuentes et al. (2022) Fuentes, J. R., Cumming, A., & Anders, E. H. 2022, Physical Review Fluids, 7, 124501, doi: 10.1103/PhysRevFluids.7.124501
  * Fuller (2014) Fuller, J. 2014, Icarus, 242, 283, doi: 10.1016/j.icarus.2014.08.006
  * Fuller et al. (2016) Fuller, J., Luan, J., & Quataert, E. 2016, Monthly Notices of the Royal Astronomical Society, 458, 3867, doi: 10.1093/mnras/stw609
  * Garaud (2018) Garaud, P. 2018, Annual Review of Fluid Mechanics, 50, 275, doi: 10.1146/annurev-fluid-122316-045234
  * Goldreich & Nicholson (1977) Goldreich, P., & Nicholson, P. D. 1977, Icarus, 30, 301, doi: 10.1016/0019-1035(77)90163-4
  * Goldreich & Soter (1966) Goldreich, P., & Soter, S. 1966, Icarus, 5, 375, doi: 10.1016/0019-1035(66)90051-0
  * Goodman & Lackner (2009) Goodman, J., & Lackner, C. 2009, ApJ, 696, 2054, doi: 10.1088/0004-637X/696/2/2054
  * Guenel et al. (2016) Guenel, M., Baruteau, C., Mathis, S., & Rieutord, M. 2016, A&A, 589, A22, doi: 10.1051/0004-6361/201527621
  * Guillot et al. (2022) Guillot, T., Fletcher, L. N., Helled, R., et al. 2022, arXiv e-prints, arXiv:2205.04100, doi: 10.48550/arXiv.2205.04100
  * Guillot et al. (2004) Guillot, T., Stevenson, D. J., Hubbard, W. B., & Saumon, D. 2004, Jupiter: The Planet, Satellites and Magnetosphere, 35, 57
  * Guillot et al. (2018) Guillot, T., Miguel, Y., Militzer, B., et al. 2018, Nature, 555, 227
  * Hedman & Nicholson (2013) Hedman, M. M., & Nicholson, P. D. 2013, The Astronomical Journal, 146, 12, doi: 10.1088/0004-6256/146/1/12
  * Ingersoll (2020) Ingersoll, A. P. 2020, Space Sci. Rev., 216, 122, doi: 10.1007/s11214-020-00751-1
  * Ivanov et al. (2013) Ivanov, P. B., Papaloizou, J. C. B., & Chernov, S. V. 2013, MNRAS, 432, 2339, doi: 10.1093/mnras/stt595
  * Jet Propulsion Laboratory (2022) Jet Propulsion Laboratory. 2022.  https://ssd.jpl.nasa.gov/sats/elem/
  * Lainey et al. (2009) Lainey, V., Arlot, J.-E., Karatekin, Ö., & Van Hoolst, T. 2009, Nature, 459, 957
  * Lainey et al. (2012) Lainey, V., Özgür Karatekin, Desmars, J., et al. 2012, The Astrophysical Journal, 752, 14, doi: 10.1088/0004-637x/752/1/14
  * Lainey et al. (2017) Lainey, V., Jacobson, R. A., Tajeddine, R., et al. 2017, Icarus, 281, 286, doi: https://doi.org/10.1016/j.icarus.2016.07.014
  * Lainey et al. (2020) Lainey, V., Casajus, L. G., Fuller, J., et al. 2020, Nature Astronomy, 4, 1053, doi: https://doi.org/10.1038/s41550-020-1120-5
  * Lebovitz (1961) Lebovitz, N. R. 1961, ApJ, 134, 500, doi: 10.1086/147175
  * Lebovitz (1989a) —. 1989a, Geophysical and Astrophysical Fluid Dynamics, 46, 221, doi: 10.1080/03091928908208913
  * Lebovitz (1989b) —. 1989b, Geophysical and Astrophysical Fluid Dynamics, 47, 225, doi: 10.1080/03091928908221822
  * Lin (2023) Lin, Y. 2023, arXiv e-prints, arXiv:2301.02418, doi: 10.48550/arXiv.2301.02418
  * Lin & Ogilvie (2018) Lin, Y., & Ogilvie, G. I. 2018, MNRAS, 474, 1644, doi: 10.1093/mnras/stx2764
  * Mankovich & Fuller (2021) Mankovich, C. R., & Fuller, J. 2021, Nature Astronomy, 5, 1103 , doi: https://doi.org/10.1038/s41550-021-01448-3
  * Marley & Porco (1993) Marley, M. S., & Porco, C. C. 1993, Icarus, 106, 508, doi: https://doi.org/10.1006/icar.1993.1189
  * Mazevet et al. (2015) Mazevet, S., Tsuchiya, T., Taniuchi, T., Benuzzi-Mounaix, A., & Guyot, F. 2015, Phys. Rev. B, 92, 014105, doi: 10.1103/PhysRevB.92.014105
  * Ogilvie (2009) Ogilvie, G. I. 2009, MNRAS, 396, 794, doi: 10.1111/j.1365-2966.2009.14814.x
  * Ogilvie (2013) —. 2013, MNRAS, 429, 613, doi: 10.1093/mnras/sts362
  * Ogilvie (2014) —. 2014, ARA&A, 52, 171, doi: 10.1146/annurev-astro-081913-035941
  * Ogilvie & Lin (2004) Ogilvie, G. I., & Lin, D. N. C. 2004, The Astrophysical Journal, 610, 477 , doi: 10.1086/421454
  * Papaloizou & Pringle (1978) Papaloizou, J., & Pringle, J. E. 1978, MNRAS, 182, 423, doi: 10.1093/mnras/182.3.423
  * Papaloizou & Savonije (1997) Papaloizou, J. C. B., & Savonije, G. J. 1997, MNRAS, 291, 651, doi: 10.1093/mnras/291.4.651
  * Pontin et al. (2023) Pontin, C. M., Barker, A. J., & Hollerbach, R. 2023, The Astrophysical Journal, 950, 176, doi: 10.3847/1538-4357/accd67
  * Pontin et al. (2020) Pontin, C. M., Barker, A. J., Hollerbach, R., André, Q., & Mathis, S. 2020, Monthly Notices of the Royal Astronomical Society, 493, 5788, doi: 10.1093/mnras/staa664
  * Remus et al. (2012) Remus, F., Mathis, S., Zahn, J.-P., & Lainey, V. 2012, A&A, 541, A165, doi: 10.1051/0004-6361/201118595
  * Remus et al. (2015) Remus, F., Mathis, S., Zahn, J. P., & Lainey, V. 2015, A&A, 573, A23, doi: 10.1051/0004-6361/201424472
  * Rieutord (2009) Rieutord, M. 2009, Approaching the low-frequency spectrum of rotating stars (Berlin, Heidelberg: Springer Berlin Heidelberg), 101 – 121, doi: 10.1007/978-3-540-87831-5_4
  * Rieutord & Valdettaro (1997) Rieutord, M., & Valdettaro, L. 1997, Journal of Fluid Mechanics, 341, 77, doi: 10.1017/S0022112097005491
  * Rieutord & Valdettaro (2010) —. 2010, Journal of Fluid Mechanics, 643, 363, doi: 10.1017/S002211200999214X
  * Storch & Lai (2014) Storch, N. I., & Lai, D. 2014, MNRAS, 438, 1526, doi: 10.1093/mnras/stt2292
  * Storch & Lai (2015) —. 2015, MNRAS, 450, 3952, doi: 10.1093/mnras/stv904
  * Terquem (2021) Terquem, C. 2021, MNRAS, 503, 5789, doi: 10.1093/mnras/stab224
  * Vazan et al. (2018) Vazan, A., Helled, R., & Guillot, T. 2018, A&A, 610, L14, doi: https://doi.org/10.1051/0004-6361/201732522
  * Wahl et al. (2017) Wahl, S. M., Hubbard, W. B., Militzer, B., et al. 2017, Geophysical Research Letters, 44, 4649, doi: 10.1002/2017GL073160
  * Wei (2018) Wei, X. 2018, ApJ, 854, 34, doi: 10.3847/1538-4357/aaa54d
  * Wood et al. (2013) Wood, T. S., Garaud, P., & Stellmach, S. 2013, ApJ, 768, 157, doi: 10.1088/0004-637X/768/2/157
  * Wu (2005) Wu, Y. 2005, The Astrophysical Journal, 635, 688, doi: 10.1086/497355
  * Zaqarashvili et al. (2021) Zaqarashvili, T. V., Albekioni, M., Ballester, J. L., et al. 2021, Space Sci. Rev., 217, 15, doi: 10.1007/s11214-021-00790-2