2.1 Orbital forcing

We are interested in the orbital problem of a companion body of mass $m$ (e.g. a moon, a gaseous planet or a low mass star), which moves on an eccentric Kepler orbit of eccentricity $e$ around an attractor of mass $M\gg m$ . The centre of mass of the attractor coincides with the centre of mass of the two-body system, which is also a focus of the eccentric Kepler orbit described by the companion body.

Figure 1.

Eccentric Kepler orbit for our two-body problem. The companion body rotates around an attractor at an orbital angular velocity orthogonal to the orbital plane of amplitude $\unicode[STIX]{x1D6FA}_{orb}(t)>0$ . The eccentric Kepler orbit, of eccentricity $e$ (thick black line) and of geometric centre $O$ , has semi-axes $a_{orb}$ and $b_{orb}$ . The perihelion (respectively aphelion) point of the orbit is $\unicode[STIX]{x1D6F1}$ (respectively $A$ ). The dashed black line is the circumscribed circle of radius $a_{orb}$ . The orbital position vector of the orbiting companion body, relating the centre of mass of the attractor to the one of the companion body, is $r(t)$ . We denote the true anomaly $\unicode[STIX]{x1D703}(t)$ and the mean anomaly $E(t)$ . Using Cartesian coordinates centred on the attractor, the position of the orbiting body is $x_{orb}=a_{orb}(\cos E-e)$ and $y_{orb}=a_{orb}\sqrt{1-e^{2}}\sin E$ . (a) A fluid ellipsoidal companion body orbiting around a point-source mass attractor. (b) A companion point-source mass orbiting around a fluid ellipsoidal attractor.

Depending on the astrophysical configuration, we consider that either the companion or the attractor is a tidally deformed homogeneous fluid body. The latter has rotating internal motions of time-dependent uniform vorticity $\unicode[STIX]{x1D74E}(t)$ . The fluid is incompressible, of uniform density $\unicode[STIX]{x1D70C}$ and kinematic viscosity $\unicode[STIX]{x1D708}$ . The other celestial body is then modelled as a rigid point-source mass. The two situations are sketched in figure 1. Under these circumstances the fluid body experiences its self-gravitating acceleration, the disturbing tidal acceleration and the centrifugal, Coriolis and Poincaré accelerations. Following Aizenman (Reference Aizenman1968), Chandrasekhar (Reference Chandrasekhar1969), Nduka (Reference Nduka1971), we keep only the quadratic terms in the expansion of the tidal potential generated by a point-source mass. Then a mathematically exact description of the fluid boundary is achieved by considering a triaxial ellipsoid, denoting the principal semi-axes $(a(t),b(t),c(t))$ . They depend on time $t$ because of the time-dependent gravitational force exerted along the eccentric orbit. The fluid ellipsoid is characterised by its equatorial and polar ellipticities

The limit $\unicode[STIX]{x1D6FD}_{ac}\rightarrow 1$ corresponds to the limit case of a disk ( $c=0$ ), whereas $\unicode[STIX]{x1D6FD}_{ac}\rightarrow -1$ corresponds to an infinite cylinder ( $c\rightarrow \infty$ ).

To describe the orbital and fluid motions, we introduce two reference frames. We define an inertial frame with fixed axes and whose origin is the centre of mass of the attractor. The horizontal plane defines the orbital plane and the vertical axis $\widehat{\boldsymbol{z}}$ is parallel to the direction of the orbital angular velocity of the companion body (of amplitude $\unicode[STIX]{x1D6FA}_{orb}(t)$ ). It is the natural frame for describing the orbital motions. A tractable frame to describe the fluid motions is the rotating body frame, whose origin is the centre of mass of the fluid ellipsoid (either the attractor or the companion). Its main axes coincide with the directions of its principal ellipsoidal axes $(a(t),b(t),c(t))$ . The body frame is rotating at the angular velocity $\unicode[STIX]{x1D734}^{{\mathcal{B}}}(t)$ with respect to the inertial frame. The general problem with an arbitrary orientation of $\unicode[STIX]{x1D74E}(t)$ with respect to $\unicode[STIX]{x1D6FA}_{orb}(t)\widehat{\boldsymbol{z}}$ is of great mathematical complexity (see the general equations (44) to (52) of Nduka Reference Nduka1971). It could be solved but it is beyond the scope of the present paper. Instead following Aizenman (Reference Aizenman1968), Chandrasekhar (Reference Chandrasekhar1969), Nduka (Reference Nduka1971), we assume a null obliquity such that $\unicode[STIX]{x1D74E}(t)=\unicode[STIX]{x1D714}(t)\widehat{\boldsymbol{z}}$ and $\unicode[STIX]{x1D734}^{{\mathcal{B}}}(t)$ is along $\widehat{\boldsymbol{z}}$ .

To make the problem dimensionless, we choose $L=\sqrt{(a_{0}^{2}+b_{0}^{2})/2}$ as the length scale, where $(a_{0},b_{0})$ are characteristic hydrostatic semi-axes of the fluid ellipsoid, and $\unicode[STIX]{x1D6FA}_{s}^{-1}$ as the time scale, where $\unicode[STIX]{x1D6FA}_{s}$ is the steady internal fluid spin rate in the inertial frame. For clarity, the dimensionless variables will be also noted as their dimensional counterparts.

The time dependencies of $\unicode[STIX]{x1D6FA}_{orb}(t)$ and $\unicode[STIX]{x1D6FD}_{ab}(t)$ are given by the orbital dynamics (figure 1). We introduce the dimensionless mean orbiting angular velocity $\unicode[STIX]{x1D6FA}_{0}$ of the body along the elliptical orbit. The orbit has main orbital semi-axes ( $a_{orb},b_{orb}$ ). Following Murray & Dermott (Reference Murray and Dermott1999), an elliptical orbit is described by Kepler’s equation at a given time $t$ ,

with $E(t)$ the eccentric anomaly. The orbital rotation rate on the elliptical orbit is

where $\unicode[STIX]{x1D703}(t)$ is the true anomaly defined by

The orbital position $r(t)$ , describing the position of the centre of mass of the companion with respect to the attractor, is

The fluid ellipsoid may have a relative orientation with respect to the orbital position vector (2.5). However the relative orientation is extremely small in the null obliquity case (Nduka Reference Nduka1971). So we assume that the tidal bulge is always aligned with the orbital position vector (instantaneous bulge response), i.e.  $\unicode[STIX]{x1D734}^{{\mathcal{B}}}(t)=\unicode[STIX]{x1D6FA}_{orb}(t)\widehat{\boldsymbol{z}}$ . We estimate at first order the equatorial ellipticity (2.1) of the fluid ellipsoid with an hydrostatic balance. Following Cébron et al. (Reference Cébron, Le Bars, Moutou and Le Gal2012a ) it reads

with $\unicode[STIX]{x1D6FD}_{0}$ a characteristic equatorial ellipticity, ${\mathcal{M}}$ the ratio between the mass of the celestial body responsible for the disturbing tidal potential and the mass of the fluid ellipsoid, $D$ the mean spherical radius of the fluid ellipsoid and $k_{2}$ the potential Love number. The latter can be computed with the Clairaut–Radau theory (e.g. Van Hoolst et al. Reference Van Hoolst, Rambaux, Karatekin, Dehant and Rivoldini2008). A typical value is $k_{2}=3/2$ for an incompressible homogeneous body in hydrostatic equilibrium (Greff-Lefftz, Métivier & Legros Reference Greff-Lefftz, Métivier and Legros2005). To take into account all the possible triaxial geometries, the polar ellipticity $\unicode[STIX]{x1D6FD}_{ac}(t)$ is a free parameter.

Note that $\unicode[STIX]{x1D6FD}_{0}$ is the ellipticity of a body with the same mass $m$ but moving on a circular orbit of radius $a_{orb}$ (dashed circle in figure 1). It is not the time averaged value of $\unicode[STIX]{x1D6FD}_{ab}(t)$ . The ellipticity $\unicode[STIX]{x1D6FD}_{0}$ refers to the static (tidal) bulge or equilibrium tide (Zahn Reference Zahn1966). The fluctuations in time superimposed on this equilibrium tide are the dynamical tides (Zahn Reference Zahn1975), which are excited by the periodic terms of the disturbing tidal potential. From formula (2.6), minimum and maximum values $(\unicode[STIX]{x1D6FD}_{min},\unicode[STIX]{x1D6FD}_{max})$ of the ellipticity at the aphelion point $A$ ( $\unicode[STIX]{x1D703}=\unicode[STIX]{x03C0}$ ) and the perihelion point $\unicode[STIX]{x1D6F1}$ ( $\unicode[STIX]{x1D703}=0$ ) are

Because the ellipticity is bounded ( $\unicode[STIX]{x1D6FD}_{ab}(t)<1$ ), the upper bound of the maximum allowable eccentricity for a given ellipticity $\unicode[STIX]{x1D6FD}_{0}$ , denoted $e_{max}$ , is

However from a physical ground, the maximum allowable ellipticity is governed by a balance between the internal cohesion model of the fluid body (e.g. self-gravitation, molecular… $\!$ ) and the disturbing tidal and centrifugal accelerations. For homogeneous self-gravitating ellipsoids, the maximum allowable ellipticity is given by the Roche limit (Aizenman Reference Aizenman1968). The classical Roche problem considers an homogeneous self-gravitating ellipsoidal body moving on a circular orbit. However as long as the eccentric orbit remains outside the Roche limiting circle, the variations of the ellipsoidal figure are small (Nduka Reference Nduka1971). Therefore in our framework we estimate a lower bound of the orbital eccentricity $e_{R}$ as the eccentricity of the first orbit crossing the Roche limiting circle. It reads

where $\unicode[STIX]{x1D6FD}^{\ast }=0.59$ is the lower bound of the equatorial ellipticity of unstable homogeneous ellipsoids moving on circular orbits in the Roche limit (estimated from point B in figure 3 of Aizenman Reference Aizenman1968). When $0\leqslant e\leqslant e_{R}$ the ellipsoidal configuration is assumed to be stable, whereas when $e_{R}\leqslant e<e_{max}$ some ellipsoidal configurations could be unstable (and hence physically unrealistic) for self-gravitating bodies.

For a circular orbit ( $e=0$ ), the orbital rotation rate is steady $\unicode[STIX]{x1D6FA}_{orb}(t)=\unicode[STIX]{x1D6FA}_{0}$ and $\unicode[STIX]{x1D6FD}_{ab}(t)=\unicode[STIX]{x1D6FD}_{0}$ . For an eccentric orbit ( $e\neq 0$ ), we determine $\unicode[STIX]{x1D6FA}_{orb}(t)$ by solving Kepler’s equation (2.2) numerically using an iterative Newton’s algorithm (starting with $E=0$ as initial guess at $t=0$ ). We show in figure 2 the normalised ratio $\unicode[STIX]{x1D6FA}_{orb}(t)/\unicode[STIX]{x1D6FA}_{0}$ and its associated power spectral density for different eccentricities. $\unicode[STIX]{x1D6FA}_{orb}(t)$ has a fundamental angular frequency $\unicode[STIX]{x1D714}=\unicode[STIX]{x1D6FA}_{0}$ but, as $e$ is increased, more and more harmonics are required to properly describe the time dependence of $\unicode[STIX]{x1D6FA}_{orb}(t)$ .

Finally it is worth noting that the case $\unicode[STIX]{x1D6FA}_{0}=1$ corresponds to a synchronised fluid body, since the dimensional averaged orbital rate and the averaged fluid spin rate are equal. When $\unicode[STIX]{x1D6FA}_{0}\neq 1$ the fluid body is not synchronised. A mean differential rotation exists between the elliptical deformation and the fluid spin rate over one spin period.

Figure 2.

(a) Normalised orbital spin rate $\unicode[STIX]{x1D6FA}_{orb}(t)/\unicode[STIX]{x1D6FA}_{0}$ given by formula (2.3) for various eccentricities $e$ . Kepler’s equation (2.2) is solved with an iterative Newton’s algorithm at each time step. (b) Associated power spectral density (PSD) as a function of the normalised angular frequency $\unicode[STIX]{x1D714}/\unicode[STIX]{x1D6FA}_{0}$ of the Fourier transform.

2.2 Fluid equations

In the body reference frame rotating at the dimensionless angular velocity $\unicode[STIX]{x1D734}^{{\mathcal{B}}}(t)=\unicode[STIX]{x1D6FA}_{orb}(t)\,\widehat{\boldsymbol{z}}$ , the time-dependent fluid boundary is ellipsoidal at any time. It is described by the equation

Because the fluid is incompressible, we restrict ourselves (without loss of generality) to the case $a(t)b(t)c(t)=1$ , such that fluid ellipsoid has a dimensionless constant volume of $4\unicode[STIX]{x03C0}/3$ . To take into account all possible triaxial ellipsoids, we define the semi-axes as

where $R$ is a free parameter governing the polar ellipticity $\unicode[STIX]{x1D6FD}_{ac}(t)$ . Note that many celestial bodies are flattened at their poles. The flattening condition, valid at each time, leads to the sufficient condition $R\geqslant R_{m}$ for a given orbit, with

and $\unicode[STIX]{x1D6FD}_{max}(e)$ defined by the expression (2.7).

For a fluid mechanics study, the knowledge of the axes $(a(t),b(t),c(t))$ , of the forcing $\unicode[STIX]{x1D734}^{{\mathcal{B}}}(t)$ and of the internal vorticity $\unicode[STIX]{x1D714}(t)$ is sufficient to fully determine the internal dynamics. The basic flow driven by the orbital motions in the fluid ellipsoid, expressed in the body frame, is

with $\boldsymbol{r}=(x,y,z)^{\text{T}}$ the position vector in the body frame, $(\widehat{\boldsymbol{x}},\widehat{\boldsymbol{y}},\widehat{\boldsymbol{z}})$ the unit Cartesian basis vectors in that frame and $\unicode[STIX]{x1D6FD}_{ab}(t)$ the time-dependent equatorial ellipticity defined by formula (2.6). It is an incompressible ( $\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{U}=0$ ) and laminar flow of uniform vorticity $\unicode[STIX]{x1D714}(t)=2[1-\unicode[STIX]{x1D6FA}_{orb}(t)]$ in the body frame. It is an exact solution of the dimensionless and nonlinear Navier–Stokes equation in the body frame

with $E_{k}=\unicode[STIX]{x1D708}/(\unicode[STIX]{x1D6FA}_{s}L^{2})$ the dimensionless Ekman number and $P$ the modified pressure. Equation (2.14) is supplemented with the impermeability condition $\boldsymbol{U}\boldsymbol{\cdot }\boldsymbol{n}=0$ , at the boundary with $\boldsymbol{n}$ the unit vector normal to the boundary (2.10). However the viscous boundary condition (either no slip or stress free) is violated. Moreover, the exact kinematic boundary condition is in fact $\boldsymbol{U}\boldsymbol{\cdot }\unicode[STIX]{x1D735}F+\unicode[STIX]{x2202}F/\unicode[STIX]{x2202}t=0$ , with $F(\boldsymbol{r},t)=1-(x/a(t))^{2}-(y/b(t))^{2}-(z/c(t))^{2}$ . Neglecting the heterogeneous term $\unicode[STIX]{x2202}F/\unicode[STIX]{x2202}t$ in the boundary condition is relevant when typically $\Vert \boldsymbol{U}\Vert \gg \Vert \unicode[STIX]{x2202}F/\unicode[STIX]{x2202}t\Vert /\Vert \unicode[STIX]{x1D735}F\Vert$ , i.e. when $\Vert \boldsymbol{U}\Vert \gg e\unicode[STIX]{x1D6FD}_{0}\unicode[STIX]{x1D6FA}_{0}$ . The latter condition is the astrophysically relevant limit. Indeed, celestial fluid bodies are typically characterised by $e\ll 1$ and $\unicode[STIX]{x1D6FD}_{0}\ll 1$ . At leading order, the basic flow is made of a solid body rotation of order $O(1)$ and a tidal correction of order $O(\unicode[STIX]{x1D6FD}_{0})$ . Corrections of the basic flow only appear at the next order $O(e\unicode[STIX]{x1D6FD}_{0})\ll 1$ and are neglected in the following.

In the literature, the stability of basic flows (2.13) has only been studied for steady values of $\unicode[STIX]{x1D6FD}_{ab}$ . We relax here this assumption. To investigate whether the basic flow (2.13) is stable against small perturbations, we perform a linear stability analysis. We expand the total velocity field into the sum of the basic flow $\boldsymbol{U}(\boldsymbol{r},t)$ (2.13) and a perturbation $\boldsymbol{u}(\boldsymbol{r},t)$ . The inviscid ( $E_{k}\ll 1$ ), linearised governing equations for the perturbation in the body frame are

with $\unicode[STIX]{x03C0}$ the modified pressure perturbation. By consistency with the basic flow, we also neglect $\unicode[STIX]{x2202}F/\unicode[STIX]{x2202}t$ in the boundary condition for the velocity perturbation. However, because the stability problem is linear, our results are not affected. Indeed the heterogeneous term $\unicode[STIX]{x2202}F/\unicode[STIX]{x2202}t$ in the boundary condition could only generate additional instabilities. Consequently we emphasise that our stability study gives sufficient conditions for instability. The basic flow $\boldsymbol{U}(\boldsymbol{r},t)$ is linearly unstable if the amplitude $\Vert \boldsymbol{u}(\boldsymbol{r},t)\Vert$ grows without bound with time.

3.1 Global method in triaxial ellipsoids

To remove the pressure term in (2.15a ), we take the curl of (2.15a ) and obtain the governing equation for the vorticity of the perturbation $\unicode[STIX]{x1D701}=\unicode[STIX]{x1D735}\times \boldsymbol{u}$ , i.e.

As originally devised by Gledzer & Ponomarev (Reference Gledzer and Ponomarev1978, Reference Gledzer and Ponomarev1992), if $\boldsymbol{u}$ is a Cartesian polynomial of maximum degree $n$ in the Cartesian coordinates $(x,y,z)$ , then each term in (3.1) is a polynomial in the Cartesian coordinates of maximum degree $n-1$ . This suggests to look for perturbations $\boldsymbol{u}$ which belong to a finite-dimensional vector space of Cartesian polynomials.

We consider the finite-dimensional vector space $\boldsymbol{{\mathcal{V}}}_{n}$ , such that an element $\boldsymbol{v}\in \boldsymbol{{\mathcal{V}}}_{n}$ is of maximum degree $n$ and satisfies (at any time) $\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{n}=0$ at the ellipsoidal boundary (2.10) and $\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{v}=0$ . Elements of $\boldsymbol{{\mathcal{V}}}_{n}$ represent vortical perturbations that are tangential to the ellipsoidal boundary at any time. The dimension of $\boldsymbol{{\mathcal{V}}}_{n}$ is (Lebovitz Reference Lebovitz1989; Backus & Rieutord Reference Backus and Rieutord2017a ; Ivers Reference Ivers2017)

Finding an appropriate basis is a difficult task. Any polynomial basis of $\boldsymbol{{\mathcal{V}}}_{n}$ is a complete basis for velocity fields defined over triaxial ellipsoids and meeting the impermeable boundary condition (Lebovitz Reference Lebovitz1989; Backus & Rieutord Reference Backus and Rieutord2017b ; Ivers Reference Ivers2017), i.e. any velocity field can be projected in theory onto $\boldsymbol{{\mathcal{V}}}_{n}$ in the limit $n\rightarrow \infty$ . Ellipsoidal harmonics, which are the eigenfunctions of the Laplace operator in ellipsoidal coordinates, are known to form a complete basis (Dassios Reference Dassios2012). Unfortunately, ellipsoidal harmonics have neither explicit expressions nor known recurrence relationships to generate them.

For a given degree $n$ , alternative bases have been proposed. Vantieghem (Reference Vantieghem2014), Backus & Rieutord (Reference Backus and Rieutord2017a ), Ivers (Reference Ivers2017) established that inertial modes, i.e. the eigenmodes of rotating fluids restored by the Coriolis force, form a basis of $\boldsymbol{{\mathcal{V}}}_{n}$ in ellipsoids rotating at steady angular velocities. Hence, the global method describes the dynamics of the perturbation $\boldsymbol{u}$ in terms of a superposition of inertial (and geostrophic) modes. This approach has been considered by Kerswell (Reference Kerswell1993b ), Kerswell & Malkus (Reference Kerswell and Malkus1998), Zhang, Chan & Liao (Reference Zhang, Chan and Liao2010), Zhang et al. (Reference Zhang, Chan and Liao2012, Reference Zhang, Chan, Liao and Aurnou2013), Zhang, Chan & Liao (Reference Zhang, Chan and Liao2014). The latter studies look at the effects of various mechanical forcings in steady spheroidal containers ( $a=b$ ), using the explicit formula of spheroidal inertial modes (Greenspan Reference Greenspan1968; Zhang, Liao & Earnshaw Reference Zhang, Liao and Earnshaw2004). On the other hand, Vantieghem et al. (Reference Vantieghem, Cébron and Noir2015) study global instabilities driven by latitudinal libration in triaxial ellipsoids, but only considering the inertial modes of degree $n\leqslant 3$ (Vantieghem Reference Vantieghem2014). The explicit spatial dependence of inertial modes is not available in triaxial ellipsoids for higher degrees (Vantieghem Reference Vantieghem2014; Backus & Rieutord Reference Backus and Rieutord2017a ), and also in spheroidal containers as soon as the time dependence of the figure axes (2.10) is taken into account.

Instead we build an algebraic polynomial basis of $\boldsymbol{{\mathcal{V}}}_{n}$ in time-dependent ellipsoids, denoted $\{\boldsymbol{v}_{i}(\boldsymbol{r},t)\}$ . The basis does not required to satisfy any dynamical equation, such that it holds at any time provided that the boundary (2.10) is ellipsoidal in the body frame. It is thus an alternative to ellipsoidal harmonics to perform spectral computations in ellipsoids. This basis has two main advantages over other ellipsoidal harmonics: (i) the Cartesian coordinate system is easier to tackle than the ellipsoidal one; (ii) the basis is explicit and can be generated for any polynomial degree $n$ .

We consider linearly independent Cartesian monomials $x^{i}y^{j}z^{k}$ of degree $i+i+k\leqslant n-1$ . The number of such independent monomials is $N_{2}=n(n+1)(n+2)/6$ . Among them, there are $N_{1}=n(n+1)/2$ monomials independent of $z$ , denoted $g_{i}(\boldsymbol{r})$ . The other monomials, denoted $h_{i}(\boldsymbol{r})$ , contain $z$ as factor. We index the set of these polynomials as (Lebovitz Reference Lebovitz1989; Barker et al. Reference Barker, Braviner and Ogilvie2016)

Next we consider the linearly independent basis elements at each time with $F(\boldsymbol{r},t)=1-(x/a(t))^{2}-(y/b(t))^{2}-(z/c(t))^{2}$ . Note that the total number of basis elements (3.2) satisfies $N_{{\mathcal{V}}}=N_{1}+2N_{2}$ . The polynomial set (3.3) ensures that basis elements (3.4) are linearly independent. It is worth noting that the basis (3.4) is neither orthogonal nor normalised, which is not necessary to build the stability equations. They can be a posteriori orthonormalised with the modified Gram–Schmidt algorithm. The basis (3.4) is explicit and can be built analytically for any degree $n$ .

We have also implemented another algorithm to build the basis of $\boldsymbol{{\mathcal{V}}}_{n}$ for arbitrary $n$ . It relies on spherical harmonics, after transforming the triaxial ellipsoid into a sphere with the Poincaré transform (Poincaré Reference Poincaré1910; Wu & Roberts Reference Wu and Roberts2011). The method is described in appendix A. However the algorithm is less efficient than the above procedure, since the basis is built numerically.

We consider perturbations $\boldsymbol{u}(\boldsymbol{r},t)\in \boldsymbol{{\mathcal{V}}}_{n}$ at any time and expand them as linear combinations of the $N_{{\mathcal{V}}}$ basis elements

where $\{\unicode[STIX]{x1D6FC}_{i}(t)\}$ is a set of arbitrary time-dependent coefficients. In the expansion (3.5) we emphasise that the basis polynomial elements are also time dependent, since the ellipsoid in the body frame has time-dependent axes. We substitute the expansion (3.5) into the stability equation (2.15) and project the resulting equation on the polynomial basis (3.4), using the (real) inner product defined by the integral over the ellipsoidal volume

Then the stability problem (2.15) yields a finite number of ordinary differential equations

where $\unicode[STIX]{x1D615}_{ij}$ , $\unicode[STIX]{x1D613}_{ij}$ and $\unicode[STIX]{x1D614}_{ij}$ are the time-dependent elements of squares matrices $\unicode[STIX]{x1D649}$ , $\unicode[STIX]{x1D647}$ and $\unicode[STIX]{x1D648}$ of size $N_{{\mathcal{V}}}\times N_{{\mathcal{V}}}$ . Explicitly these elements are given by

Note that the pressure term does not contribute to (3.8b ). We compute the elements (3.8) explicitly using the formula (misprint corrected from Lebovitz Reference Lebovitz1989)

with $2\unicode[STIX]{x1D6FE}_{1}=i$ , $2\unicode[STIX]{x1D6FE}_{2}=j$ , $2\unicode[STIX]{x1D6FE}_{3}=k$ , $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D6FE}_{1}+\unicode[STIX]{x1D6FE}_{2}+\unicode[STIX]{x1D6FE}_{3}$ and $\unicode[STIX]{x1D6FE}!=\unicode[STIX]{x1D6FE}_{1}!\unicode[STIX]{x1D6FE}_{2}!\unicode[STIX]{x1D6FE}_{3}!$ . Stability equations (3.7) are written in canonical matrix form

with the unknown vector $\unicode[STIX]{x1D736}(t)=(\unicode[STIX]{x1D6FC}_{1}(t),\unicode[STIX]{x1D6FC}_{2}(t),\ldots )^{\text{T}}$ and $\unicode[STIX]{x1D645}$ is the time-dependent Jacobian matrix of the system.

The basic flow $\boldsymbol{U}(\boldsymbol{r},t)$ is unstable if the perturbation $\boldsymbol{u}(\boldsymbol{r},t)$ has at least one modal coefficient $\unicode[STIX]{x1D6FC}_{i}(t)$ , governed by the stability equation (3.10), which grows without bound in time. The most unstable perturbation is associated with the fastest growth rate denoted $\unicode[STIX]{x1D70E}$ . Since the basic flow (2.13) is periodic of period $T=2\unicode[STIX]{x03C0}/\unicode[STIX]{x1D6FA}_{0}$ , its stability can be determined using the Floquet theory. We compute the eigenvalues (Floquet exponents) $\{\unicode[STIX]{x1D707}_{i}\}$ of the fundamental matrix $\unicode[STIX]{x1D731}(t)$ evaluated at time $T$ . The fundamental matrix is solution of

with $\unicode[STIX]{x1D644}$ the identity matrix. Then the growth rates $\{\unicode[STIX]{x1D70E}_{i}\}$ and the frequencies $\{\unicode[STIX]{x1D714}_{i}\}$ associated with the flow perturbations are then

The fastest growth rate, associated with the most dangerous unstable flow, is $\unicode[STIX]{x1D70E}=\max _{i}\unicode[STIX]{x1D70E}_{i}$ and its associated angular frequency is $\unicode[STIX]{x1D714}$ . Note that in the special case of a circular orbit ( $e=0$ ) the orbital forcing (2.3) is steady, such that the above Floquet analysis reduces to a classical eigenvalue stability analysis.

A key point of the global method is that the vector space $\boldsymbol{{\mathcal{V}}}_{n}$ is invariant under the action of the perturbation stability equation (2.15) (Kerswell Reference Kerswell1993b ; Lebovitz Reference Lebovitz1989; Vantieghem Reference Vantieghem2014; Backus & Rieutord Reference Backus and Rieutord2017b ; Ivers Reference Ivers2017). Thus the stability equation (2.15) exactly reduces to (3.10) for any finite value of the degree $n$ in the expansion (3.5). It is not an approximation and it is not necessary to replace $n$ by $\infty$ in the expansion (3.5), and then to truncate at finite $n$ , to get the stability equation (3.10). It is a main difference with classical spectral Galerkin expansions in various geometries (e.g. spherical harmonics), even in the linear framework. For our purposes, it is another advantage of our properly chosen polynomial basis over the spectral basis of ellipsoidal harmonics.

Finally, because the expansion (3.5) is exact, the global method at a given degree $n$ gives only exact sufficient conditions for inviscid instability (to the numerical precision of the time step solver). New tongues of instability of the basic flow generally appear in the parameter space when $n$ increases. The largest growth rate also generally increases by considering larger and larger $n$ , reaching progressively its asymptotic value. This phenomenon has already been noticed in previous global analyses performed at lower degrees $n\leqslant 7$ (e.g. Kerswell Reference Kerswell1993b ; Wu & Roberts Reference Wu and Roberts2011, Reference Wu and Roberts2013; Vantieghem et al. Reference Vantieghem, Cébron and Noir2015; Barker et al. Reference Barker, Braviner and Ogilvie2016). Indeed, more resonances are expected when $n$ increases. So we conclude that the global method at finite values of $n$ cannot prove the stability of the basic flow, but it gives sufficient conditions for instability. Finally, it is usually expected that the upper bounds of global growth rates, in the asymptotic limit $n\rightarrow \infty$ , coincide with the growth rates of local perturbations of short wavelength (see § 3.2). However no general mathematical proof is available to justify it.

3.2 Local method in unbounded fluids

To get a complementary physical understanding of fluid instabilities growing upon the basic flow (2.13), we also perform a local Wentzel–Kramers–Brillouin (WKB) stability analysis. It probes the stability of any inviscid, three-dimensional and time-dependent basic flow in an unbounded fluid, considering localised plane wave perturbations of small wavelength which are advected along the basic flow (Friedlander & Vishik Reference Friedlander and Vishik1991; Lifschitz & Hameiri Reference Lifschitz and Hameiri1991; Friedlander & Lipton-Lifschitz Reference Friedlander, Lipton-Lifschitz, Friedlander and Serre2003). Because the orbitally driven basic flow (2.13) is linear in space coordinates, the short-wavelength perturbations exactly reduce to Kelvin waves (Bayly Reference Bayly1986; Craik & Criminale Reference Craik and Criminale1986; Craik Reference Craik1989; Waleffe Reference Waleffe1990)

with $\boldsymbol{k}(t)$ the time-dependent wavenumber and $\boldsymbol{a}(t)$ the time-dependent amplitude of the velocity perturbation. Kelvin waves (3.13) are exact inviscid, nonlinear and incompressible solutions upon the basic flow (2.13) in the body frame if

and the incompressibility condition $\boldsymbol{k}(t)\boldsymbol{\cdot }\boldsymbol{a}(t)=0$ hold. The latter is always satisfied if it holds for the initial condition $(\boldsymbol{k}_{0},\boldsymbol{a}_{0})$ . The existence of an unbounded solution for $\boldsymbol{a}(t)$ is a sufficient condition for instability (Lifschitz & Hameiri Reference Lifschitz and Hameiri1991; Friedlander & Lipton-Lifschitz Reference Friedlander, Lipton-Lifschitz, Friedlander and Serre2003).

Equations (3.14) are independent of the magnitude of $\boldsymbol{k}_{0}$ . So we restrict the initial wave vector to the spherical surface of unit radius $\boldsymbol{k}_{0}=(\sin (\unicode[STIX]{x1D703}_{0})\cos (\unicode[STIX]{x1D719}),\sin (\unicode[STIX]{x1D703}_{0})\sin (\unicode[STIX]{x1D719}),\cos (\unicode[STIX]{x1D703}_{0}))^{\text{T}}$ , where $\unicode[STIX]{x1D719}\in [0,2\unicode[STIX]{x03C0}]$ is the longitude and $\unicode[STIX]{x1D703}_{0}\in [0,\unicode[STIX]{x03C0}]$ is the colatitude between the spin axis $\widehat{\boldsymbol{z}}$ and the initial wave vector $\boldsymbol{k}_{0}$ . In practice, the equation (3.14a ) is time stepped with an ordinary differential equation solver from a range of initial wave vectors. Then we compute the maximum growth rate $\unicode[STIX]{x1D70E}$ of (3.14b ) as the fastest growing solution from all possible initial wave vectors.

3.3 Numerical implementation

We have developed for the global stability analysis the SIREN code (Stability with IneRtial EigeNmodes), freely available at https://bitbucket.org/vidalje/siren. It handles any mechanical forcing and basic flow of uniform vorticity. The matrices $\unicode[STIX]{x1D649}$ , $\unicode[STIX]{x1D647}$ and $\unicode[STIX]{x1D648}$ are first computed symbolically with Sympy (http://www.sympy.org/), a computer algebra system (CAS) for Python, which is used to manipulate the Cartesian polynomials $x^{i}y^{j}z^{k}$ in a symbolic way. Then they are converted to Fortran subroutines with the Sympy fcode function and finally wrapped with f2py (Peterson Reference Peterson2009) for fast numerical evaluation inside Python using Numpy (Van Der Walt, Colbert & Varoquaux Reference Van Der Walt, Colbert and Varoquaux2011). The Jacobian matrix $\unicode[STIX]{x1D645}$ is computed numerically, because we cannot compute the symbolic inverse $\unicode[STIX]{x1D649}^{-1}$ for arbitrary $n$ . Because of the difficulty to build the Jacobian matrix for an arbitrary forcing, previous global studies have only considered less than 200 basis elements ( $n\leqslant 7$ ) (Kerswell Reference Kerswell1993b ; Lebovitz & Lifschitz Reference Lebovitz and Lifschitz1996a ; Wu & Roberts Reference Wu and Roberts2011; Vantieghem et al. Reference Vantieghem, Cébron and Noir2015; Barker Reference Barker2016a ; Barker et al. Reference Barker, Braviner and Ogilvie2016). In practice, we have built and solved numerically the stability system (3.10) for degrees as large as $n=25$ , yielding more than 6000 basis elements.

Figure 3.

Characteristic CPU time to compute a growth rate for tidally driven and libration-driven instabilities presented below in § 4 with the different numerical solvers. To smooth out the variability of computation time between different parameters, we compute a stability map of 100 points in the plane $(\unicode[STIX]{x1D6FD}_{ac},\unicode[STIX]{x1D6FA}_{0})$ for tides and in the plane $(\unicode[STIX]{x1D6FD}_{0},2e)$ for longitudinal libration to extract an average time for one iteration. Circle symbols stand for tidal forcing and square ones for longitudinal libration. The magenta solid line shows the power law $\propto n^{6}$ , in good agreement with the numerical scaling.

For the local stability analysis we have also developed the SWAN code (Short-Wavelength stability Analysis), freely available at https://bitbucket.org/vidalje/swan. It gives sufficient conditions for inviscid instability of any basic flow (not necessarily of uniform vorticity) expressed in Cartesian coordinates. The stability equations (3.14) are built using Sympy, then converted to a Fortran subroutine with the Sympy fcode function and finally wrapped with f2py for fast numerical evaluations with Numpy.

Both numerical codes use an explicit Runge–Kutta time step solver with adaptive step size (available in the Python library Scipy) to integrate the stability differential equations. Performing a survey in parameter space is an embarrassingly parallel problem, and our implementation takes full advantage of this situation using mpi4py (http://mpi4py.scipy.org/).

To validate our codes, we have first considered the precession of a steady spheroid ( $a=b\neq c$ ). This benchmark is described in appendix B. We perfectly recover previous studies (Kerswell Reference Kerswell1993b ; Wu & Roberts Reference Wu and Roberts2011). We also get a very good agreement between local and global analyses, because we can reach large enough degrees $n$ with the SIREN code. Then we assess the performance of our SIREN code in figure 3, which shows the evolution of CPU time with $n$ for the tidally driven and libration-driven flows considered below (see §§ 4.1 and 4.2). We observe that the CPU time scales as $n^{6}$ , in agreement with formula (3.2). Indeed the number of basis elements scales as $n^{3}$ and so the number of elements in matrices $\unicode[STIX]{x1D649}$ , $\unicode[STIX]{x1D647}$ and $\unicode[STIX]{x1D648}$ is of order $n^{6}$ . As expected the eigenvalue solver is faster than the Floquet solver.

4.1 Tidally driven elliptical instability on circular orbits

We focus here on the effect of the equilibrium tide on a circular orbit ( $e=0$ ). The fluid ellipsoid has steady semi-axes ( $a,b,c$ ) and rotates at the steady orbital rate $\unicode[STIX]{x1D6FA}_{0}$ . Basic flow (2.13) thus reduces to the tidally driven basic flow

This flow can be unstable if $\unicode[STIX]{x1D6FA}_{0}\neq 1$ , leading to the classical tidally driven elliptical instability (TDEI).

On one hand, the TDEI has been widely studied with a local analysis in unbounded domains (Bayly Reference Bayly1986; Craik Reference Craik1989; Waleffe Reference Waleffe1990; Cébron et al. Reference Cébron, Le Bars, Moutou and Le Gal2012a ). Note that in the asymptotic limit $\unicode[STIX]{x1D6FD}_{0}\rightarrow 0$ , equations (3.14) can be solved analytically using a multiple-scale analysis in $\unicode[STIX]{x1D6FD}_{0}$ to get a theoretical growth rate (see § C.1). Le Dizès (Reference Le Dizès2000) shows that the TDEI exists in the range $(\unicode[STIX]{x1D6FD}_{0}~+~1)/(\unicode[STIX]{x1D6FD}_{0}~-~1)<\unicode[STIX]{x1D6FA}_{0}<3$ . Outside this range, the flow is stable and lies in the classical forbidden zone for $\unicode[STIX]{x1D6FD}_{0}\ll 1$ , hereafter denoted $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ .

Figure 4.

Areas of instability of the tidally driven flow in the $(\unicode[STIX]{x1D6FD}_{0},\unicode[STIX]{x1D6FA}_{0})$ plane. Colour map shows $\log _{10}(\unicode[STIX]{x1D70E}/[\unicode[STIX]{x1D6FD}_{0}|1-\unicode[STIX]{x1D6FA}_{0}|])$ . Triaxial geometry $a=\sqrt{1+\unicode[STIX]{x1D6FD}_{0}}$ , $b=\sqrt{1-\unicode[STIX]{x1D6FD}_{0}}$ and $c=1/(ab)$ such that the triaxial container has a constant dimensionless volume $4\unicode[STIX]{x03C0}/3$ . On the vertical white line $\unicode[STIX]{x1D6FA}_{0}=1$ the TDEI does not exist. Black dashed lines $\unicode[STIX]{x1D6FA}_{0}=(1~+~\unicode[STIX]{x1D6FD}_{0})/(1~-~\unicode[STIX]{x1D6FD}_{0})$ and $\unicode[STIX]{x1D6FA}_{0}=3$ are the bounds of the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ .

Figure 5.

(a,b) Areas of instability of the tidally driven flow in the $(\unicode[STIX]{x1D6FD}_{ac},\unicode[STIX]{x1D6FA}_{0})$ plane at degree $n=15$ . Colour map shows $\log _{10}(\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D6FD}_{0})$ . White areas correspond to marginally stable regions. Triaxial geometry $a=\sqrt{1+\unicode[STIX]{x1D6FD}_{0}}$ and $b=\sqrt{1-\unicode[STIX]{x1D6FD}_{0}}$ and $c=a\sqrt{(1-\unicode[STIX]{x1D6FD}_{ac})/(1+\unicode[STIX]{x1D6FD}_{ac})}$ . Vertical dashed black lines represent the lower and upper bounds of the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ . The solid black line indicates the synchronised case $\unicode[STIX]{x1D6FA}_{0}=1$ (no instability). White dashed lines correspond to the isoline $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D6FD}_{0}=0.01$ for the stability problem reduced to degree $n=1$ , such that the spin-over instability is excited in between (in this case, $\unicode[STIX]{x1D70E}$ is analytically known, see e.g. Roberts & Wu Reference Roberts and Wu2011).

On the other hand, the global stability analysis of tidal basic flow (4.1) has been mainly performed for weakly deformed spheroids (Lacaze, Le Gal & Le Dizès Reference Lacaze, Le Gal and Le Dizès2004) or cylinders (Malkus Reference Malkus1989; Eloy, Le Gal & Le Dizès Reference Eloy, Le Gal and Le Dizès2003). Triaxial ellipsoids have also been considered (Gledzer & Ponomarev Reference Gledzer and Ponomarev1978, Reference Gledzer and Ponomarev1992; Kerswell Reference Kerswell2002; Roberts & Wu Reference Roberts and Wu2011; Barker Reference Barker2016a ; Barker et al. Reference Barker, Braviner and Ogilvie2016), but only allowing perturbations of small polynomial degrees ( $n\leqslant 7$ ). Such large-wavelength instabilities do not compare well with the aforementioned local stability analyses. Using our framework, we can reach much larger polynomial degrees $n$ . A comparison between local and numerical analyses is given in figure 4. As expected, the local results are upper bounds of global results (with expected matches for large enough polynomial degrees). A more in-depth discussion is given in § C.1.

Global stability results at maximum degree $n=15$ are shown in figure 5, where the ratio $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D6FD}_{0}$ of the instability is computed for two equatorial ellipticities ( $\unicode[STIX]{x1D6FD}_{0}=0.15$ and $\unicode[STIX]{x1D6FD}_{0}=0.6$ ). We always find $\unicode[STIX]{x1D70E}=0$ when $\unicode[STIX]{x1D6FA}_{0}=1$ as expected, because $\boldsymbol{U}=\mathbf{0}$ from the expression (4.1). When $\unicode[STIX]{x1D6FD}_{0}$ increases, the zone of instability extends but it is still outside of the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ .

We observe that ellipsoids spinning in the retrograde direction ( $\unicode[STIX]{x1D6FA}_{0}<0$ ) are more unstable than the prograde ones ( $\unicode[STIX]{x1D6FA}_{0}>0$ ). By varying the polynomial degree, we note that the TDEI for prograde rotation ( $\unicode[STIX]{x1D6FA}_{0}>0$ ) appears at larger $n$ than the TDEI for retrograde rotation ( $\unicode[STIX]{x1D6FA}_{0}<0$ ). As an example the spin-over mode (Kerswell Reference Kerswell2002), associated with the linear basis ( $n=1$ ), appears only for $\unicode[STIX]{x1D6FA}_{0}<0$ . Similarly the largest $\unicode[STIX]{x1D70E}$ is reached at smaller $n$ for retrograde rotation than for prograde rotation (not shown). We also observe an effect of $\unicode[STIX]{x1D6FD}_{ac}$ at large values of $|\unicode[STIX]{x1D6FD}_{ac}|$ , not predicted by the local analysis (insensitive to $\unicode[STIX]{x1D6FD}_{ac}$ ). There, higher polynomial degrees may be needed to reach the asymptotic limit of the local analysis to completely fill in the map with new global unstable tongues.

Figure 6.

Three-dimensional renderings of the most unstable flows associated with figure 5. Degree $n=15$ and amplitude of equilibrium tide $\unicode[STIX]{x1D6FD}_{0}=0.15$ . Velocity magnitude $\Vert \boldsymbol{u}\Vert$ is shown in meridional/equatorial planes and at the ellipsoidal surface. The colour map is saturated for $\Vert \boldsymbol{u}\Vert \geqslant 3$ .

In figure 6 we show the most dangerous unstable flows of the TDEI, as varying $\unicode[STIX]{x1D6FA}_{0}$ when $\unicode[STIX]{x1D6FD}_{0}=0.15$ . When $\unicode[STIX]{x1D6FD}_{0}=0.6$ the spatial structures of the flows, at the same values of $\unicode[STIX]{x1D6FA}_{0}$ , are similar (not shown). In all these flows, the motions seem to be concentrated in conical layers tilted from the fluid rotation axis. Between these layers, the flow has low or zero amplitude. Some flows also exhibit one or several nodes in azimuth.

In figure 6(a) computed at $\unicode[STIX]{x1D6FA}_{0}=-1$ , the flow has the particular structure of a stack of pancakes (SoP). Note that the modal angular frequency of SoP is $\unicode[STIX]{x1D714}=0$ in the body frame, as predicted theoretically by Lebovitz & Lifschitz (Reference Lebovitz and Lifschitz1996b ), Barker et al. (Reference Barker, Braviner and Ogilvie2016). SoP structures have also been observed in experiments (Grannan et al. Reference Grannan, Le Bars, Cébron and Aurnou2014) and in nonlinear direct simulations (Favier et al. Reference Favier, Grannan, Le Bars and Aurnou2015; Barker Reference Barker2016a ). For this instability, each pancake moves horizontally in the direction opposite to its neighbours (horizontal epicyclic motions), independently of all other heights, in a plane at $45^{\circ }$ from the main equatorial axis where the stretching is maximum (Waleffe Reference Waleffe1990). SoP structures are illustrated by the streamlines in figure 7. High amplitudes are located near the poles. Note also that the number of pancakes increases as $n$ increases, as suggested by figure 7. However this number seems to be insensitive to the amplitude of the equilibrium tide $\unicode[STIX]{x1D6FD}_{0}$ and to $\unicode[STIX]{x1D6FD}_{ac}$ (not shown). Such a small-scale flow will undoubtedly lead to turbulence if it reaches high enough amplitudes.

Figure 7.

‘Stack of pancakes’-like instability for $\unicode[STIX]{x1D6FD}_{0}=0.15,\unicode[STIX]{x1D6FD}_{ac}=0.5$ and $\unicode[STIX]{x1D6FA}_{0}=-1$ . Magnitude $\Vert \boldsymbol{u}\Vert$ and velocity streamlines in the meridional plane at $45$ degrees from the long axis where the stretching is maximum.

When $\unicode[STIX]{x1D6FA}_{0}=-0.5$ the unstable flow shown in figure 6(b) for $n=15$ is identical for degrees $n=10$ and $n=6$ . So we expect the observed flow to be the most unstable flow when $\unicode[STIX]{x1D6FA}_{0}=-0.5$ . It is mostly an equatorially symmetric mode, with high intensities located on the rotation axis and within a surface band at mid-latitudes. Finally flows in figure 6(c,d) share a common structure which could evolve with $n$ .

4.2 Libration-driven elliptical instability

We investigate here the stability of a synchronised fluid body ( $\unicode[STIX]{x1D6FA}_{0}=1$ ) on an eccentric orbit ( $0<e\leqslant 1$ ). The associated forcing, called longitudinal librations, leads to the libration-driven elliptical instability (LDEI). We distinguish the following two limit cases of longitudinal librations.

If the fluid is enclosed in a solid container with a strong enough rigidity (e.g. a silicate mantle), the entire body rigidly rotates with a fixed shape. Dynamical tides can be neglected with respect to the equilibrium tide, such that $\unicode[STIX]{x1D6FD}_{ab}(t)=\unicode[STIX]{x1D6FD}_{0}$ . The forcing bears the name of physical librations. A differential rotation exists between the fluid spin rate and the equilibrium tide, rotating at leading order at the angular velocity (equal to the orbital angular velocity in our model)

with $\unicode[STIX]{x1D716}\leqslant 2e$ the libration amplitude. This amplitude depends on the rheology of the celestial body. The LDEI driven by physical librations has been studied amongst others by Cébron et al. (Reference Cébron, Le Bars, Noir and Aurnou2012b ), Noir et al. (Reference Noir, Cébron, Le Bars, Sauret and Aurnou2012), Grannan et al. (Reference Grannan, Le Bars, Cébron and Aurnou2014) and Favier et al. (Reference Favier, Grannan, Le Bars and Aurnou2015). Note that the physical libration forcing (4.2) could actually contain multiple frequencies due to the presence of several attracting bodies (Rambaux & Williams Reference Rambaux and Williams2011). In the limit $\unicode[STIX]{x1D6FD}_{0}\ll 1$ , the local growth rate of this physical LDEI is (Herreman, Le Bars & Le Gal Reference Herreman, Le Bars and Le Gal2009; Cébron et al. Reference Cébron, Le Bars, Noir and Aurnou2012b ; Cébron, Vantieghem & Herreman Reference Cébron, Vantieghem and Herreman2014)

According to the local formula (4.3) the LDEI is triggered for any non-zero $\unicode[STIX]{x1D716}$ and $\unicode[STIX]{x1D6FD}_{0}$ .

On the other hand, if the fluid spin rate is low enough for the shape of the fluid body to have time to adapt to the gravitational constraints or is enclosed within a solid container with a small enough rigidity (e.g. a thin icy shell), the ellipsoidal cavity always points toward the attractor. Then the container has a time-dependent equatorial ellipticity given by the expression (2.6). A differential rotation exists between the fluid spin rate and the dynamical tides (superimposed on the equilibrium tide). In the inviscid framework of this work, we call this forcing optical librations (because the amplitude of optical librations is $2e$ , see e.g. Murray & Dermott Reference Murray and Dermott1999). In the limit $e\rightarrow 0$ , this forcing simply associates a prescribed time evolution of $(a(t),b(t),c(t))$ to the forcing (4.2), rather than considering a constant ellipsoidal shape. At this first order in $e$ , the time dependence of the dynamical tides is monochromatic, in agreement with numerical results of figure 2 at small $e$ . Physical librations with maximum amplitude $\unicode[STIX]{x1D716}=2e$ are recovered if we neglect the dynamical tides, yielding $\unicode[STIX]{x1D6FD}_{ab}(t)=\unicode[STIX]{x1D6FD}_{0}$ .

Figure 8.

Libration-driven elliptical instability ( $\unicode[STIX]{x1D6FA}_{0}=1$ ) for eccentric Kepler orbits. (a) Ratio $\log _{10}(\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb})$ with $\unicode[STIX]{x1D70E}_{wkb}$ the local growth rate of LDEI given by the formula (4.3). Colour bar saturated at $\log _{10}(\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb})=1$ . Larger values ( $\log _{10}(\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb})\leqslant 2$ ) occur for small enough $\unicode[STIX]{x1D6FD}_{0}$ when $[\unicode[STIX]{x1D6FD}_{max}(e)-\unicode[STIX]{x1D6FD}_{min}(e)]/\unicode[STIX]{x1D6FD}_{0}\geqslant 3$ (not shown). Solid red line shows the Roche limiting circular orbit. Ellipsoids moving along orbits located above the Roche limit could be unstable against free-surface perturbations (not considered in this work). (b) Flow magnitude $\Vert \boldsymbol{u}\Vert$ for $\unicode[STIX]{x1D6FD}_{0}=0.05,e/e_{max}=0.4$ . Flow computed at $\unicode[STIX]{x1D703}(t)=\unicode[STIX]{x03C0}/2$ on the orbit (see figure 1).

Figure 9.

Survey of the orbitally driven elliptical instability. (a,b) Range of values of $\unicode[STIX]{x1D6FD}_{ab}$ (solid blue line) and $\unicode[STIX]{x1D6FD}_{ac}$ (dashed red line) for the various eccentric orbits considered. Blue thick (respectively red dashed) line shows the mean value of $\unicode[STIX]{x1D6FD}_{ab}$ (respectively $\unicode[STIX]{x1D6FD}_{ac}$ ) along the orbits. The second vertical axis shows the maximum of ratio $\unicode[STIX]{x1D6FA}_{orb}/\unicode[STIX]{x1D6FA}_{0}$ (green dotted line). (c,d) Growth rates $\unicode[STIX]{x1D70E}$ of the ODEI in the plane $(e/e_{max},\unicode[STIX]{x1D6FA}_{0})$ for degree $n=10$ . The colour bar is saturated at $\unicode[STIX]{x1D70E}\geqslant 0.2$ in (c) and $\unicode[STIX]{x1D70E}=0.6$ in (d). White areas correspond to marginally stable regions. The containers considered are oblate ellipsoids with $R=R_{m}+0.05$ . Vertical black line corresponds to the synchronised case ( $\unicode[STIX]{x1D6FA}_{0}=1$ ) driving the LDEI (see § 4.2). The horizontal line $e=0$ corresponds to the TDEI (see § 4.1). Vertical dashed black lines are the bounds of the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ of the classical TDEI, valid for $e=0$ and $\unicode[STIX]{x1D6FD}_{ab}=\unicode[STIX]{x1D6FD}_{0}$ . Horizontal magenta line is $e_{R}/e_{max}$ , with $e_{R}$ defined by formula (2.9). In (a) the line is outside of the plot ( $e_{R}/e_{max}=0.89$ ).

Figure 10.

Velocity magnitude $\Vert \boldsymbol{u}\Vert$ of several unstable flows of the figure 11(c). $n=10$ , $e/e_{max}=0.4$ and $\unicode[STIX]{x1D6FD}_{0}=0.05$ . $\Vert \boldsymbol{u}\Vert$ is shown in meridional/equatorial planes and at the outer ellipsoidal surface. The colour map is saturated for $\Vert \boldsymbol{u}\Vert \geqslant 3$ . Flows are computed at $\unicode[STIX]{x1D703}(t)=\unicode[STIX]{x03C0}/2$ on the orbit (see figure 1).

We consider the general optical LDEI, taking into account the exact orbital motion (2.3) and associated dynamical tides. We survey in figure 8(a) the optical LDEI on eccentric orbits, varying the equilibrium tide $\unicode[STIX]{x1D6FD}_{0}$ and the eccentricity $e$ from the circular case to $e/e_{max}=0.8$ . Two distinctive behaviours occur. The transition is associated with a physical change in the main tidal effect. To compare the effects of dynamical and equilibrium tides, we introduce the normalised ratio

with $\unicode[STIX]{x1D6FD}_{max}(e)$ (respectively $\unicode[STIX]{x1D6FD}_{min}(e)$ ) the maximum (respectively minimum) equatorial ellipticity for a given eccentricity as defined in expressions (2.7). Physically when $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\ll 1$ the equilibrium tide $\unicode[STIX]{x1D6FD}_{0}$ is of prime importance compared with dynamical tides. For weakly eccentric orbits ( $e\rightarrow 0$ ), the growth rates of the optical LDEI coincide with the ones of the physical LDEI predicted by formula (4.3), as shown by the unit ratio $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}=1$ . However we observe that new tongues of instability, with normalised growth rates $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}\geqslant 1$ , appear when $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\leqslant 1$ ( $e\leqslant 0.15$ ). These new tongues are not predicted by the local formula (4.3). This phenomenon is already visible at large $\unicode[STIX]{x1D6FD}_{0}$ in figure 8(a), computed at $n=10$ . For smaller values of $\unicode[STIX]{x1D6FD}_{0}$ , several tongues also appear but higher degrees $(n\geqslant 15)$ are required to catch them all (not shown). Thus the LDEI can be more vigorous than predicted before, even in the range of small eccentricities relevant in geo and astrophysics, with growth rates as high as $10\,\unicode[STIX]{x1D70E}_{wkb}$ .

On the other hand when $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\geqslant 1$ , the effects of dynamical tides overcome the ones of the equilibrium tide. The eccentricity of the orbit ( $e\geqslant 0.15$ ) plays now a fundamental role in the tidal effects and the fluid body tends to forget its equilibrium tide. In figure 8(b) we show the most dangerous unstable flow for an equilibrium tide of amplitude $\unicode[STIX]{x1D6FD}_{0}=0.05$ and on an orbit of eccentricity $e/e_{max}=0.4$ . Violent instabilities occur with growth rates $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}\geqslant 10$ figure 8(a), whatever the value of $\unicode[STIX]{x1D6FD}_{0}$ , and can even reach extreme values $10\leqslant \unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}\leqslant 100$ when $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\geqslant 3$ (not shown in a). This latter effect occurs for small enough $\unicode[STIX]{x1D6FD}_{0}$ and highly eccentric orbits ( $e\geqslant 0.34$ ), a situation relevant for planets.

4.3 Survey of the orbitally driven elliptical instability

The general case of a fluid ellipsoid orbiting on a Kepler orbit of eccentricity $0\leqslant e<1$ is now considered. Since the Kepler equation (2.2) is solved at any time step (as in § 4.2), the computational cost is more expensive than for computations done in § 4.1. As in § 4.2, we fix the polynomial degree to the value $n=10$ to survey the whole parameter space.

We survey the stability of the orbitally driven basic flow in figure 9 as each parameter in the set $(\unicode[STIX]{x1D6FA}_{0},e)$ varies. We arbitrary fix the equilibrium to a small value ( $\unicode[STIX]{x1D6FD}_{0}=0.05$ ) and to a larger one ( $\unicode[STIX]{x1D6FD}_{0}=0.3$ ). In figure 9(a,b) we show the average, maximum and minimum values of the equatorial and polar ellipticities ( $\unicode[STIX]{x1D6FD}_{ab}(t),\unicode[STIX]{x1D6FD}_{ac}(t)$ ) along the orbits. We consider here only oblate containers ( $b>c$ ), which typically describe the shapes of celestial bodies. The maximum value of the orbital angular velocity (normalised by $\unicode[STIX]{x1D6FA}_{0}$ ) is also shown.

Figure 11.

Growth rates of the ODEI for eccentric Kepler orbits at $e/e_{max}=0.4$ compared to the growth rates of the TDEI (red thick solid line). Shaded areas represent the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ for the TDEI ( $\unicode[STIX]{x1D6FA}_{0}\leqslant (1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)$ and $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ ).

Figure 12.

Comparison between local and global stability analyses of orbitally driven flows (2.13) for $\unicode[STIX]{x1D6FD}_{0}=0.3$ and $e/e_{max}=0.4$ . We compute the local maximum growth rate as the fastest growing solution from a range of initial wave vectors $\boldsymbol{k}_{0}=(\sin \unicode[STIX]{x1D703}\cos \unicode[STIX]{x1D719},\sin \unicode[STIX]{x1D703}\sin \unicode[STIX]{x1D719},\cos \unicode[STIX]{x1D703})$ , with 0.5 degree spacing in $\unicode[STIX]{x1D703}\in [0,90^{\circ }]$ . The initial azimuthal angle of the local wave vector $\unicode[STIX]{x1D719}=\unicode[STIX]{x03C0}/4$ . The latter maximises the growth rate of the classical TDEI (Le Dizès Reference Le Dizès2000). We have checked that the largest growth rates are insensitive to the value of $\unicode[STIX]{x1D719}$ . (a) Local (dashed blue line) and global $n=15$ (solid red line) growth rates $\unicode[STIX]{x1D70E}$ in function of $\unicode[STIX]{x1D6FA}_{0}$ . (b) Numerical angle $\unicode[STIX]{x1D703}_{0}$ of the initial local wave vector leading to the maximum growth rate as varying $\unicode[STIX]{x1D6FA}_{0}$ . Thin green line shows the classical destabilising angle leading to the TDEI and the LDEI.

Figure 9(c,d) shows the maximum growth rates of the most unstable modes. Some of the associated unstable flows are shown in figure 10. First, the maximum growth rate in each panel tends to increase when $\unicode[STIX]{x1D6FD}_{0}$ increases from (c) to (d). Then several aspects of figure 9(c,d) are worthy of comment.

We recover the TDEI considered in § 4.1, corresponding to the horizontal line $e=0$ . We also show the bounds of the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ , i.e.  $\unicode[STIX]{x1D6FA}_{0}=(1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)$ and $\unicode[STIX]{x1D6FA}_{0}=3$ (dashed grey lines). The instability with the largest wavelength ( $n=1$ basis) is the spin-over mode (not shown). It is similar to the ‘middle-moment-of-inertia’ instability of rigid bodies. The spin-over develops on any circular orbits ( $e=0$ ) for retrograde rotations $-1\leqslant \unicode[STIX]{x1D6FA}_{0}\leqslant 0$ , as expected from previous global analyses (e.g. Roberts & Wu Reference Roberts and Wu2011; Barker et al. Reference Barker, Braviner and Ogilvie2016), but also for any eccentric orbit ( $0\leqslant e$ ) with an almost constant growth rate (see figure 19 in appendix D). Then an increasing region of the parameter space becomes unstable as the polynomial degree $n$ is increased from $n=1$ to $n=10$ , within the expected allowable range $(1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}<3$ where the classical TDEI develops on circular orbits. We observe that the eccentricity has little effect on the growth rates of the TDEI for retrograde orbits (exceptions occur for large eccentricities) within the allowable range $(1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)\leqslant \unicode[STIX]{x1D6FA}_{0}\leqslant 0$ . For instance we recover the SoP unstable modes at $\unicode[STIX]{x1D6FA}_{0}=-1$ in figure 10(b). For other values of $\unicode[STIX]{x1D6FA}_{0}$ , the exact flow structure of the unstable modes depends on the triaxial shape. So it prevents from directly comparing with flows in figure 6 (also obtained at larger $n$ ), although showing broad common patterns. This first unexpected result justifies a posteriori the validity of TDEI mechanism in tidally disturbed planets or stars, since it can be extended to eccentric retrograde orbits within the classical allowable range of the TDEI ( $(1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}\leqslant 0$ ).

However, our survey also shows that dynamical tides strongly enhance two instability tongues for prograde eccentric orbits. The growth rates are indeed much larger than the ones predicted on circular orbits. The first one is associated with the LDEI at $\unicode[STIX]{x1D6FA}_{0}=1$ as previously discussed in § 4.2. The second tongue seems to be centred on $\unicode[STIX]{x1D6FA}_{0}=2+\unicode[STIX]{x1D6FD}_{0}$ . The most unstable flow at $n=10$ is shown in figure 10(f). It exhibits intense motion localised in patches around the equator. It is very different from the TDEI flow at $\unicode[STIX]{x1D6FA}_{0}=2$ in figure 6(d). We expect this localisation to increase as $n$ is increased further. The enhancing of the growth rate first appears at degree $n=2$ , for large enough $e$ (not shown). So it is not associated with the spin-over mode ( $n=1$ ). Then the instability band moves towards smaller eccentricities as $n$ increases (not shown), even when the effects of dynamical tides are a priori small ( $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\leqslant 1$ ).

Another striking result is observed in figure 9(c,d). We uncover new violent instabilities within the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ for both retrograde ( $\unicode[STIX]{x1D6FA}_{0}\leqslant (1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)$ ) and prograde eccentric orbits ( $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ ). The unstable tongue for retrograde orbits appears first at large enough eccentricities and it is not initially associated with the spin-over mode (see figure 19 in appendix D). Then increasing the degree $n$ shows that this new tongue is replaced by more unstable tongues which merge with the tongue of the TDEI near $\unicode[STIX]{x1D6FA}_{0}=(1+\unicode[STIX]{x1D6FD}_{0})/(\unicode[STIX]{x1D6FD}_{0}-1)$ . The latter tongue also spreads towards more retrograde orbits for large enough eccentricities (when $n$ is large enough). An example of an unstable flow in this tongue at $\unicode[STIX]{x1D6FA}_{0}=-2$ and $\unicode[STIX]{x1D6FD}_{0}=0.05$ is shown in figure 10(a). It displays vertical stripes that seem similar to the SoP observed at $\unicode[STIX]{x1D6FA}_{0}=-1$ (but here stacked along an equatorial axis). For $\unicode[STIX]{x1D6FA}_{0}=-3$ and $\unicode[STIX]{x1D6FD}_{0}=0.3$ , the unstable flow is instead a SoP flow (see the discussion of figure 12 in § 5 below). For prograde eccentric orbits ( $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ ), these new instabilities are initially associated with an unstable tongue of degree $n=3$ . Then, there is also a merging between this tongue and the one appearing near $\unicode[STIX]{x1D6FA}_{0}=2+\unicode[STIX]{x1D6FD}_{0}$ at $n=2$ , which spreads out towards more and more prograde orbits ( $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ ) for large eccentricities when $n$ increases. Note that these new unstable tongues exist for orbits of eccentricities $e\leqslant e_{R}$ in figure 9. So these new instabilities may physically exist in astrophysical fluid bodies.

We compare now more quantitatively the strength of these instabilities in figure 11 by pushing the degree to $n=15$ . We show the growth rates of the TDEI on circular orbits (red solid curve) and the ones of orbitally driven instabilities on eccentric Kepler orbits for a finite value of the eccentricity ( $e/e_{max}=0.4$ ). When $(\unicode[STIX]{x1D6FD}_{0}+1)/(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}\leqslant 1$ , the growth rate of the ODEI has almost the same value as for the classical TDEI at $e=0$ . We note that increasing the degree $n$ yields small variations in $\unicode[STIX]{x1D70E}$ . The maximum $\unicode[STIX]{x1D70E}$ is well predicted by the local WKB analysis of the TDEI (see formula (C 1) in § C.1), showing a scaling in $\unicode[STIX]{x1D6FD}_{0}$ . Around $\unicode[STIX]{x1D6FA}_{0}=1$ we see the peak corresponding to the LDEI, previously computed for various $\unicode[STIX]{x1D6FD}_{0}$ in figure 8(a) at $n=10$ . The largest growth rates are obtained for the instability located at $\unicode[STIX]{x1D6FA}_{0}=2+\unicode[STIX]{x1D6FD}_{0}$ (in the limit $e\ll 1$ ). Its growth rates are approximately ten times larger than the ones predicted by the classical TDEI circular orbits at the same values of $\unicode[STIX]{x1D6FA}_{0}$ . This instability becomes stronger as $n$ is increased from $n=6$ to $n=15$ , suggesting a rather small-scale instability. We expect its growth rate $\unicode[STIX]{x1D70E}$ to reach an upper bound for large enough $n$ . However we recall that the global method gives only sufficient conditions for instability. So even though the growth rate has not reached yet its asymptotic value, it does not physically rule out the enhanced strength of this instability. The new instabilities driven by the non-circular Kepler orbits within the classical forbidden zone of the TDEI are also clearly visible. The unstable tongues extend deeply inside the forbidden zone, even at low $\unicode[STIX]{x1D6FD}_{0}$ . The growth rate is almost insensitive to the chosen $n$ from $n=10$ to $n=15$ . It suggests that the asymptotic growth rates have already been (at least at $e/e_{max}=0.4$ ).

To sum up, we have found new sufficient conditions for inviscid instability (for the values of $n$ considered here). They show that the orbitally driven basic flow (figure 9) can be unstable in the allowable region of the classical TDEI, but also inside the classical forbidden zone for both retrograde and prograde eccentric orbits. We also show in appendix D that these instabilities are recovered in three-dimensional viscous numerical simulations. Therefore we expect this phenomenon to hold in astrophysical bodies.

5.1 Local approach

In this section, we discuss the physical mechanism responsible for the orbitally driven instabilities. First we perform a local stability analysis by solving equations (3.14) in unbounded fluids. Indeed the nature of an unstable tongue is be related to the colatitude $\unicode[STIX]{x1D703}_{0}$ of the initial wave vector $\boldsymbol{k}_{0}$ (see § 3.2) leading to the largest growth rate.

We aim at solving analytically the stability equations in the limit of weakly eccentric orbits ( $e\ll 1$ ). We expand at first order in $e$ the orbital forcing (2.3) and (2.6) to get

At this order, the time dependence of $\unicode[STIX]{x1D6FD}_{ab}(t)$ is monochromatic, in agreement with numerical results of figure 2 at small $e$ . When $e=0$ , the basic state is a pure solid-body rotation and it admits plane inertial wave perturbations (e.g. Greenspan Reference Greenspan1968). They have periodic wave vectors which are orthogonal to the local velocity vector. Plane inertial waves exist when $-1<\unicode[STIX]{x1D6FA}_{0}<3$ , which is the allowable range of the classical TDEI on circular orbits when $e\rightarrow 0$ .

The basic mechanism of the elliptical instability is a parametric resonance between a pair of inertial waves and the basic flow, provided that certain resonance conditions are met (e.g. Le Dizès Reference Le Dizès2000; Kerswell Reference Kerswell2002). This mechanism also applies here. Indeed when $e\rightarrow 0$ , two inertial waves can resonate with the orbitally driven basic flow (2.13) to drive an instability. The latter is governed by an Hill–Schrödinger equation, which can be readily obtained following Naing & Fukumoto (Reference Naing and Fukumoto2009). It is not written here for the sake of clarity. However the forcing term in the Hill–Schrödinger equation is not strictly periodic, as for tidally driven and libration-driven basic flows. Instead it is quasi-periodic with multiple forcing frequencies, such that many resonances are possible. For a given forcing frequency $f$ , the condition of perfect resonance yields

where the left hand side is actually the inertial waves pulsation (Doppler shifted in the rotating frame) and $\widetilde{\unicode[STIX]{x1D6FA}}_{0}=\unicode[STIX]{x1D6FA}_{0}/(1-\unicode[STIX]{x1D6FA}_{0})$ . The forcing frequencies are $f=1,2,\widetilde{\unicode[STIX]{x1D6FA}}_{0}\ldots$ and possible combinations of them (through cosine and sine products). The nature of the unstable tongues, determined by $\unicode[STIX]{x1D703}_{0}$ , depends on the considered forcing frequency $f$ . The frequency $f=2$ is associated with the classical TDEI on circular orbits (Waleffe Reference Waleffe1990; Le Dizès Reference Le Dizès2000) and $f=1$ with the LDEI on weakly eccentric orbits (Herreman et al. Reference Herreman, Le Bars and Le Gal2009; Cébron et al. Reference Cébron, Le Bars, Noir and Aurnou2012b , Reference Cébron, Vantieghem and Herreman2014). For the TDEI, SoP instabilities have for instance wave vectors aligned with the spin axis $\unicode[STIX]{x1D703}_{0}=0$ (Lebovitz & Lifschitz Reference Lebovitz and Lifschitz1996a ; Barker et al. Reference Barker, Braviner and Ogilvie2016).

Equation (5.2) shows that resonances associated with a given frequency $f$ only occur for values of $\unicode[STIX]{x1D6FA}_{0}$ located outside of the forbidden band $|f/(4(1+\widetilde{\unicode[STIX]{x1D6FA}}_{0}))|\geqslant 1$ . We recover the existence of the classical TDEI ( $f=2$ ) inside the allowed region $-1<\unicode[STIX]{x1D6FA}_{0}<3$ . When $f=\widetilde{\unicode[STIX]{x1D6FA}}_{0}$ , the forbidden band is $|\unicode[STIX]{x1D6FA}_{0}|\geqslant 4$ . However for finite values of $\unicode[STIX]{x1D6FD}_{0}$ , the unstable tongues have finite widths of order $O(\unicode[STIX]{x1D6FD}_{0})$ in the limit $0\leqslant e\ll 1$ . Hence it is possible to excite imperfect resonances by geometric detuning, even though the condition (5.2) is not strictly satisfied. A wider range of the parameter space is thus unstable when $\unicode[STIX]{x1D6FD}_{0}$ increases, as observed in §§ 4.1–4.3. For instance the classical TDEI is excited inside the allowable range $(\unicode[STIX]{x1D6FD}_{0}+1)(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}<3$ for finite values of $\unicode[STIX]{x1D6FD}_{0}$ (see figure 4). Therefore considering all the possible frequencies, the resonance condition (5.2) shows that the orbitally driven instabilities may a priori be triggered well outside the allowed region of the TDEI.

We further simplify the orbital forcing (5.1) to consider two limiting simplified configurations. Firstly we neglect the dynamical tides (i.e.  $\unicode[STIX]{x1D6FD}(t)=\unicode[STIX]{x1D6FD}_{0}$ ) to isolate the modulation of the background rotation (i.e.  $\unicode[STIX]{x1D6FA}_{orb}(t)=\unicode[STIX]{x1D6FA}_{0}[1+2e\cos (\unicode[STIX]{x1D6FA}_{0}t)]$ ). Such a forcing refers to the tidal forcing of a telluric (i.e. rigid) planet moving on an eccentric Kepler orbit. At leading order in $e$ , the associated forcing frequency is $f=2$ and the resonance gives the classical growth rate of the TDEI (see § C.1). It shows that the time modulation of the background rotation does not destabilise further the tidal basic flow (at leading order in $e$ ). Higher-order terms in $e$ may be necessary to handle possible new effects. Secondly we disable the Coriolis force ( $\unicode[STIX]{x1D6FA}_{orb}(t)=0$ ), but retains the time dependence of the ellipticity along the orbit. In this case, we obtain that the angle $\unicode[STIX]{x1D703}_{0}=\unicode[STIX]{x03C0}/3$ is the most destabilising one for rapidly oscillating tides ( $|\unicode[STIX]{x1D6FA}_{0}|\gg 1$ ), leading to $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D6FD}_{0}=9/16$ (as confirmed by solving numerically equations (3.14) in this configuration). The latter formula is identical to the growth rate of the classical TDEI without background rotation (see § C.1). It shows that dynamical tides are the key physical mechanism responsible for the instabilities located outside of the allowable range of the classical TDEI, as observed in figures 9 and 11. However we note that the associated growth rates are overestimated with respect to the growth rates obtained numerically for the full problem, suggesting that the Coriolis force has a stabilising effect.

To get quantitative local predictions in agreement with the global results, we consider the combined effect of rotation and dynamical tides. However obtaining a growth rate formula for any possible resonance, predicted by the equation (5.2), is not of practical interest. Indeed resonances may appear or be modified when the full orbital forcing is considered (even at small $e$ ), possibly leading to more unstable tongues. Consequently we solve numerically the local stability equation (3.14) with the SWAN code, taking into account the full orbital forcing (2.3)–(2.6). We show the comparison between local and global analyses in figure 12 for $\unicode[STIX]{x1D6FD}_{0}=0.3$ and $e/e_{max}=0.4$ . Results obtained at smaller $\unicode[STIX]{x1D6FD}_{0}$ are similar and do not change the interpretation. We first note in figure 12(a) a good agreement between local and global growth rates for retrograde orbits ( $-4\leqslant \unicode[STIX]{x1D6FA}_{0}\leqslant 0$ ) and for prograde orbits within the forbidden zone ( $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ ). When $(\unicode[STIX]{x1D6FD}_{0}+1)/(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}\leqslant 0$ , the angle $\cos \unicode[STIX]{x1D703}_{0}=[2(1+\widetilde{\unicode[STIX]{x1D6FA}}_{0})]^{-1}$ in figure 12(b) shows that the ODEI reduces to the classical TDEI, which is not modified by the orbital eccentricity. We also recover the SoP instabilities ( $\unicode[STIX]{x1D703}_{0}=0$ ) when $-2\leqslant \unicode[STIX]{x1D6FA}_{0}\leqslant -1.5$ for $\unicode[STIX]{x1D6FD}_{0}=0.3$ . Within the forbidden zone of the classical TDEI for retrograde orbits, we find two new tongues of instability not predicted by the TDEI resonance. The modulation of the global rotation is responsible for these instabilities when $-3\leqslant \unicode[STIX]{x1D6FA}_{0}\leqslant -2$ (not shown). These instabilities were not obtained analytically in the asymptotic limit $e\rightarrow 0$ , because they are due to higher-order terms. When $\unicode[STIX]{x1D6FA}_{0}=-3$ we find a SoP instability, which is also obtained in the global analysis and direct numerical simulations (see appendix D).

Then we find a new tongue of instability within the classical forbidden zone for both rapidly oscillating prograde orbits $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ and retrograde orbits $\unicode[STIX]{x1D6FA}_{0}\leqslant -4$ . Dynamical tides are responsible for these instabilities. The numerical local growth rates are in much closer agreement to the global ones than the analytical growth rate $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D6FD}_{0}=9/16$ obtained without background rotation. We conclude that the Coriolis force has a stabilising effect on these instabilities. Angle $\unicode[STIX]{x1D703}_{0}$ , initially $\unicode[STIX]{x1D703}_{0}=\unicode[STIX]{x03C0}/3$ without background rotation (i.e. for $\unicode[STIX]{x1D6FA}_{orb}(t)=0$ ), is modulated by the rotation. It explains the observed linear trend with $\unicode[STIX]{x1D6FA}_{0}$ in figure 12(b).

Figure 13.

Survey of the libration-driven elliptical instability for physical librations (4.2) and first-order optical librations in the plane $(\unicode[STIX]{x1D6FD}_{0},2e)$ . Polynomial degree $n=10$ . Colour map shows the ratio $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}$ with $\unicode[STIX]{x1D70E}_{wkb}$ given by formula (4.3). White areas correspond to marginally stable regions. Triaxial geometry $a(t)=\sqrt{1+\unicode[STIX]{x1D6FD}_{ab}(t)}$ , $b(t)=\sqrt{1-\unicode[STIX]{x1D6FD}_{ab}(t)}$ and $c(t)=[a(t)b(t)]^{-1}$ . (a) $\unicode[STIX]{x1D6FD}_{ab}(t)=\unicode[STIX]{x1D6FD}_{0}$ and $\unicode[STIX]{x1D716}=2e$ . (b) $\unicode[STIX]{x1D6FD}_{ab}(t)=\unicode[STIX]{x1D6FD}_{0}(1+3e\cos t)$ .

The last striking result in figure 12 is that the local analysis does not predict the enhancing of growth rates at $\unicode[STIX]{x1D6FA}_{0}=1$ (LDEI) and $\unicode[STIX]{x1D6FA}_{0}=2+\unicode[STIX]{x1D6FD}_{0}$ . For any finite value of $e$ , no local instability is found. A possible explanation is that we find numerical wavenumbers $\boldsymbol{k}(t)$ with secular growths, which challenge the validity of the local analysis. Indeed, local instabilities are generally obtained under the assumption of bounded and asymptotically non-decaying periodic or quasi-periodic wavenumbers (e.g. Eckhardt & Yao Reference Eckhardt and Yao1995). We also recall that the local analysis gives only sufficient conditions for instability. Finally it has already been observed that a local analysis can be in disagreement with a global analysis. For instance, global radiative instabilities in compressible Rankine vortex (Broadbent & Moore Reference Broadbent and Moore1979) are not predicted by a local WKB analysis (Le Duc Reference Le Duc2001).

5.2 Global approach

The enhancing of the growth rates when $\unicode[STIX]{x1D6FA}_{0}=1,2+\unicode[STIX]{x1D6FD}_{0}$ is left unexplained by the local analysis. Nevertheless, the global analysis provides an explanation for this phenomenon. In ellipsoids, the elliptical instability is also a parametric resonance between a pair of inertial modes and the basic flow, provided certain resonance conditions are met (e.g. Kerswell Reference Kerswell2002; Lacaze et al. Reference Lacaze, Le Gal and Le Dizès2004; Le Bars et al. Reference Le Bars, Lacaze, Le Dizès, Le Gal and Rieutord2010). However, the present orbital forcing challenges this classical instability mechanism. Indeed, we cannot define properly inertial modes in our time-dependent fluid ellipsoids. So identifying the possible resonant couplings is a difficult task. We have to isolate the most unstable modes from the computations and to try to relate them to some inertial modes of a well-chosen ellipsoidal shape (for instance the one associated with the equilibrium tide $\unicode[STIX]{x1D6FD}_{0}$ ). Such an approach relies on an empirical modal decomposition (e.g. Schmid Reference Schmid2010; Sieber, Paschereit & Oberleithner Reference Sieber, Paschereit and Oberleithner2016), which is beyond the scope of the present study.

The key phenomenon is the time-dependent ellipticity $\unicode[STIX]{x1D6FD}(t)$ , even in the limit $e\ll 1$ . To illustrate this point we focus on the first-order optical libration forcing, making use of the forcing (5.1) with $\unicode[STIX]{x1D6FA}_{0}=1$ . In the asymptotic local (WKB) analysis of the LDEI on weakly elliptical orbits ( $e\rightarrow 0$ ), $\unicode[STIX]{x1D6FD}_{0}$ and $e$ are supposed to be of the same order of magnitude in the asymptotic expansion (Herreman et al. Reference Herreman, Le Bars and Le Gal2009; Cébron et al. Reference Cébron, Le Bars, Noir and Aurnou2012b , Reference Cébron, Vantieghem and Herreman2014). So the leading-order effect is the physical libration forcing (4.2), whereas the time-variable tidal effect (of order $e\unicode[STIX]{x1D6FD}_{0}$ ) is a priori of second order. However, the latter effect can become of primary importance if it is large enough, which cannot be probed by the local analysis for the coupled forcing (as explained before).

We compare in figure 13 physical librations (a) with optical librations (b), assuming an amplitude $\unicode[STIX]{x1D716}=2e$ for physical librations. We consider only perturbations of maximum degree $n=10$ . Perturbations of higher degrees are not essential for this comparison. We show the ratio $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}$ to compare the global growth rates $\unicode[STIX]{x1D70E}$ with the local ones $\unicode[STIX]{x1D70E}_{wkb}$ predicted by the formula (4.3). The global growth rates in figure 13(a) do not reach yet the asymptotic local growth rates (4.3) for all the values of $\unicode[STIX]{x1D6FD}_{0}$ , as expected with a global analysis at $n=10$ . In figure 13(b) the unstable tongues generated by physical librations coincide with the unstable tongues generated by optical librations in the limit $e\rightarrow 0$ . However for finite values of $e$ (even small), we observe that the tongues in figure 13(b) are much wider because of the dynamical tides. Moreover new violent instabilities are triggered, with growth rates much larger than those predicted by formula (4.3). It clearly illustrates the enhancing of the LDEI. Note that a similar behaviour is obtained for the unstable tongue appearing at $\unicode[STIX]{x1D6FA}_{0}=2+\unicode[STIX]{x1D6FD}_{0}$ (not shown).

Figure 14.

Growth rates of instabilities obtained for synchronised bodies $\unicode[STIX]{x1D6FA}_{0}=1$ and pseudo-synchronised bodies with $\unicode[STIX]{x1D6FA}_{0}(e)$ given by the formula (5.3). Ratio $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}$ with $\unicode[STIX]{x1D70E}_{wkb}$ given by formula (4.3). Triaxial geometry $a(t)=R\sqrt{1+\unicode[STIX]{x1D6FD}_{ab}(t)}$ , $b(t)=R\sqrt{1-\unicode[STIX]{x1D6FD}_{ab}(t)}$ and $c(t)=[a(t)b(t)]^{-1}$ with $R=R_{m}+0.05$ .

Finally for the libration forcing it is possible to partially remove the effect of dynamical tides. The tidal torque averaged over weakly eccentric orbits ( $e\rightarrow 0$ ) vanishes for $\unicode[STIX]{x1D6FA}_{0}=1$ . However it no longer holds for eccentric orbits. Following Hut (Reference Hut1981), the tidal torque averaged over an eccentric orbit of eccentricity $e$ vanishes when

When $\unicode[STIX]{x1D6FA}_{0}$ is given by the expression (5.3), the body has reached a pseudo-synchronised state. Pseudo-synchronisation is an important process in the dynamics of binary systems (e.g. Hut Reference Hut1981, Reference Hut1982). Indeed pseudo-synchronisation proceeds much faster than circularisation of the orbits (Zahn Reference Zahn2008). So the fluid spin rate of a celestial body in an eccentric orbit would tidally evolve towards pseudo-synchronisation (5.3), while the orbit remains eccentric of orbital period $2\unicode[STIX]{x03C0}/\unicode[STIX]{x1D6FA}_{0}$ .

We show in figure 14 the growth rates $\unicode[STIX]{x1D70E}$ , normalised by the growth rates of the classical LDEI given by the local formula (4.3). We consider both the synchronised case $\unicode[STIX]{x1D6FA}_{0}=1$ and the pseudo-synchronised case $\unicode[STIX]{x1D6FA}_{0}(e)$ given by formula (5.3). When the pseudo-synchronisation is reached, the growth rates remain bounded within the range $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}\leqslant 3$ for $\unicode[STIX]{x1D6FD}_{0}\leqslant 0.3$ . A normalisation with respect to the TDEI (see the formula (C 2) in § C.1) leads to erroneous growth rates for pseudo-synchronised bodies. So we conclude that the LDEI formula (4.3) gives a good estimate of the growth rates of orbitally driven instabilities in pseudo-synchronised bodies. Instead for synchronised states ( $\unicode[STIX]{x1D6FA}_{0}=1$ ) the instabilities can be much more vigorous than those predicted by the local formula (4.3), with $\unicode[STIX]{x1D70E}/\unicode[STIX]{x1D70E}_{wkb}\geqslant 10$ if $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\geqslant 1$ (as previously discussed in § 4.2).

Therefore the net non-zero tidal torque operating along eccentric Kepler orbits is responsible for the enhancement of the instabilities at $\unicode[STIX]{x1D6FA}_{0}=1$ . A similar effect also exists at $\unicode[STIX]{x1D6FA}_{0}=2+\unicode[STIX]{x1D6FD}_{0}$ . So it proves that the dynamical tides are essential to explain the observed enhancing of the growth rates. For the particular case of pseudo-synchronised orbits (5.3), the effects of dynamical tides do not overcome the ones of equilibrium tide on the fluid instabilities, as shown in figure 14. The onset of instabilities in pseudo-synchronised bodies is well predicted by the local analysis of the LDEI (4.3). So the quantitative predictions of Cébron et al. (Reference Cébron, Le Bars, Le Gal, Moutou, Leconte and Sauret2013) for the onset of the elliptical instability in bloated pseudo-synchronised hot Jupiters, based upon the formula (C 2)), may have to be reassessed.

6.1 Physical implications

Rotating fluid ellipsoids have been the subject of many works, going back to Riemann (Reference Riemann1860). Their stability is affected by free-surface perturbations, associated with surface gravity modes, and internal hydrodynamic perturbations. Surprisingly, free-surface perturbations weakly affect the stability of fluid ellipsoids (Lebovitz & Lifschitz Reference Lebovitz and Lifschitz1996a ; Barker Reference Barker2016a ; Barker et al. Reference Barker, Braviner and Ogilvie2016). Thus their stability is mainly governed by flow instabilities. Because the viscosity is extremely small in astrophysical bodies, an inviscid analysis is physically relevant. Previous hydrodynamic studies have to be completed, because Lebovitz & Lifschitz (Reference Lebovitz and Lifschitz1996a ,Reference Lebovitz and Lifschitz b ) consider isolated ellipsoidal fluid masses and Barker (Reference Barker2016a ), Barker et al. (Reference Barker, Braviner and Ogilvie2016) ellipsoids moving on circular orbits.

To simplify the problem, we have considered only the case where the mass of the attractor is much larger than the mass of the companion body. It is the simplest framework to model two-body systems such as synchronised moons around planets, gaseous extrasolar planets (Hot Jupiters) around stars or low massive stars orbiting around massive attractors. The general two-body problem could also be tackled (e.g. a binary stellar system), solely by changing the hydrostatic estimation of the equatorial ellipticity in formula (2.6). The radius $r(t)$ has to be replaced by the time-dependent distance between the centres of mass of the two bodies, which are both moving on eccentric orbits. We are confident that our main findings will not change qualitatively in that configuration.

We have revisited the hydrodynamic instabilities of homogeneous, incompressible and rotating ellipsoidal fluid masses subjected to a disturbing tidal potential. Several studies are devoted to the stability of fluid ellipsoids subjected to a tidal potential generated by orbital motions on circular orbits (Aizenman Reference Aizenman1968; Cébron et al. Reference Cébron, Le Bars, Moutou and Le Gal2012a , Reference Cébron, Le Bars, Le Gal, Moutou, Leconte and Sauret2013; Barker et al. Reference Barker, Braviner and Ogilvie2016). Thus our primary purpose was to study how the hydrodynamic stability of fluid ellipsoids is modified by considering a tidal potential generated by orbital motions on eccentric Kepler orbits. Our study is complementary to the hydrodynamic stability analysis of Barker (Reference Barker2016a ) of Roche–Riemann ellipsoids on circular orbits. We recover all the limiting cases of ellipsoidal flow instability (Lebovitz & Lifschitz Reference Lebovitz and Lifschitz1996a ,Reference Lebovitz and Lifschitz b ; Cébron et al. Reference Cébron, Le Bars, Le Gal, Moutou, Leconte and Sauret2013; Barker Reference Barker2016a ) and unify them into a global framework.

We may summarise our rather unexpected results in the following way. First, the classical TDEI is unaffected by the dynamical tides for retrograde eccentric orbits for $(1+\unicode[STIX]{x1D6FD}_{0})(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}\leqslant 0$ (outside of the forbidden zone). Second, instabilities excited on moderately eccentric orbits can have larger growth rates than those on nearly circular orbits. Indeed dynamical tides are responsible for the enhancing of the vigour of the TDEI near the 2:1 spin–orbit resonance ( $\unicode[STIX]{x1D6FA}_{0}\simeq 2$ ) and of the LDEI ( $\unicode[STIX]{x1D6FA}_{0}=1$ ). Finally, fluid ellipsoids exhibit new fluid instabilities which are triggered within the forbidden zone of the classical TDEI for retrograde ( $\unicode[STIX]{x1D6FA}_{0}\leqslant (1+\unicode[STIX]{x1D6FD}_{0})(\unicode[STIX]{x1D6FD}_{0}-1)$ ) and prograde ( $\unicode[STIX]{x1D6FA}_{0}\geqslant 3$ ) orbits. All these findings show that dynamical tides can drive new instabilities in fluid bodies moving on eccentric Kepler orbits.

We have updated the picture of the linear stability of tidally disturbed fluid ellipsoids. A complete view emerges now. They are prone to various local and global inviscid instabilities. On one hand, spheroids are only unstable against free-surface perturbations, associated with surface gravity modes (Lebovitz & Lifschitz Reference Lebovitz and Lifschitz1996a ; Barker et al. Reference Barker, Braviner and Ogilvie2016). Considering tidal effects generated by orbital motions on circular orbits, all ellipsoids (the Roche–Riemann ellipsoids) are unstable against the elliptical instability when $(1+\unicode[STIX]{x1D6FD}_{0})(\unicode[STIX]{x1D6FD}_{0}-1)<\unicode[STIX]{x1D6FA}_{0}<3$ , as predicted by previous analyses (Cébron et al. Reference Cébron, Le Bars, Moutou and Le Gal2012a , Reference Cébron, Le Bars, Le Gal, Moutou, Leconte and Sauret2013; Barker Reference Barker2016a ). Taking into account instabilities of all possible spatial complexity handled by the local and global theories, the parameter space of fluid ellipsoids subjected to a varying tidal torque (eccentric orbits) is unstable against orbitally driven instabilities for both retrograde and prograde eccentric orbits within the range $-10\leqslant \unicode[STIX]{x1D6FA}_{0}\leqslant 10$ (see figure 12). Although not considered in our computations, we also expect them to be intrinsically unstable on a wider range of $|\unicode[STIX]{x1D6FA}_{0}|$ (at least for large enough eccentricities).

Table 1.

Orbital parameters of some stars with companions orbiting on eccentric Kepler orbits. $P_{s}=2\unicode[STIX]{x03C0}/\unicode[STIX]{x1D6FA}_{s}$ is the spin period (in days) and $P_{orb}=2\unicode[STIX]{x03C0}/\unicode[STIX]{x1D6FA}_{orb}$ the orbital period (in days). The last column $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}$ is defined by formula (4.4). The given stars are located within the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ of the classical TDEI. Adapted from Cébron et al. (Reference Cébron, Le Bars, Le Gal, Moutou, Leconte and Sauret2013). Data have been updated from http://exoplanet.eu/.

Our findings may have important consequences for the tidal dissipation responsible for the circularisation and synchronisation of two-body systems. Some stars, located within the forbidden zone $\text{FZ}_{\unicode[STIX]{x1D6FD}_{0}}$ of the classical TDEI and with companions orbiting on eccentric Kepler orbits, are reported in table 1 as example. The effects of dynamical tides are not negligible ( $\unicode[STIX]{x0394}\unicode[STIX]{x1D6FD}_{ab}/\unicode[STIX]{x1D6FD}_{0}\sim 0.54$ for WASP-14). Consequently we expect these stars to be unstable for the orbitally driven instability (in spite of their presence in the classical forbidden zone). Two tidal dissipation processes have received most attention, namely tidal friction of the equilibrium tide (Zahn Reference Zahn1966) and tidal friction of eigenmodes forced by dynamical tides (Zahn Reference Zahn1975; Ogilvie & Lin Reference Ogilvie and Lin2004; Wu Reference Wu2005a ,Reference Wu b ; Ogilvie & Lin Reference Ogilvie and Lin2007; Goodman & Lackner Reference Goodman and Lackner2009; Ogilvie Reference Ogilvie2009; Rieutord & Valdettaro Reference Rieutord and Valdettaro2010). The elliptical instability is thus an alternative and promising mechanism (Rieutord Reference Rieutord2004, Reference Rieutord2008).

6.2 Perspectives

The nonlinear outcome of these fluid instabilities remains elusive in astrophysical bodies. Indeed it is not clear whether turbulent flows can develop and sustain an effective mixing in fluid interiors. A parameter survey of their nonlinear behaviours, using efficient numerical simulations, is necessary. For instance Le Reun et al. (Reference Le Reun, Favier, Barker and Le Bars2017) show that the saturation of the elliptical instability generates turbulence exhibiting many signatures of inertial wave turbulence, a regime possibly expected in planetary interiors. Barker (Reference Barker2016a ) also suggests that the elliptical instability may explain the spin synchronisation and circularisation of the shortest period hot Jupiters.

Future work is also required to adopt more realistic interior models. In particular the behaviour of these instabilities when a stable stratification (like a stellar radiative zone) is present is almost unknown, as well as the role of compressibility. For instance stars of mass larger than $1.8$ solar mass are stably stratified in their outer layers. Circularisation and synchronisation are also effective for these stars (e.g. Giuricin, Mardirossian & Mezzetti Reference Giuricin, Mardirossian and Mezzetti1984). Depending on the considered density profile, the stratification does favour or not the elliptical instability (Miyazaki Reference Miyazaki1993; Kerswell Reference Kerswell1993a ; Le Bars & Le Dizés Reference Le Bars and Le Dizés2006; Cébron, Maubert & Le Bars Reference Cébron, Maubert and Le Bars2010c ; Clausen & Tilgner Reference Clausen and Tilgner2014). However a unifying theory is still missing and has to be addressed. In addition, the capability of tidal effects to drive self-sustained magnetic field is still controversial. Tidal dynamos were first addressed by Barker & Lithwick (Reference Barker and Lithwick2013), Cébron & Hollerbach (Reference Cébron and Hollerbach2014), neglecting density effects.

Finally, we have developed two open-source numerical codes that may be useful for future linear stability studies of incompressible fluids. They are quite general and can be applied to several other situations. They are freely available for the community at https://bitbucket.org/vidalje/. The SWAN code allows to study the local stability of any time-dependent basic flow in unbounded fluids. The SIREN code allows to study the global stability of any mechanically driven flows of uniform vorticity in ellipsoids. Unlike previous studies, our code also handles arbitrary (time-dependent) ellipsoidal shapes (not limited to small departures from the sphere). Indeed, it handles ellipsoidal perturbations of unprecedented small wavelengths. We have considered in the present study polynomial degrees as large as $n=25$ . It corresponds to more than 6000 basis elements.

A fork of the SIREN code has been used to compute (i) the tilted hydromagnetic eigenmodes of a fluid in a corotating ellipsoid (Vidal, Cébron & Schaeffer Reference Vidal, Cébron and Schaeffer2016) and (ii) the viscous decay factors of inertial modes (Lemasquerier et al. Reference Lemasquerier, Grannan, Vidal, Cébron, Favier, Le Bars and Aurnou2017). Computing the inertial modes is the first step towards a complete and self-contained viscous stability analysis of inertial instabilities in arbitrary rotating ellipsoids. Indeed, viscous effects can be introduced as a correction of inviscid inertial modes. It will allow direct comparisons between the theory and simulations or experiments, all performed at finite values of viscosity and deformation.