Open Access
Issue
A&A
Volume 647, March 2021
Article Number A144
Number of page(s) 25
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/202039148
Published online 24 March 2021

© A. Astoul et al. 2021

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

Tidal interactions are known to drive the late evolution of short-period planetary systems, such as hot Jupiters orbiting around their host star and, in our Solar System, the satellites around Jupiter and Saturn (e.g. Ogilvie 2014; Mathis 2019). In particular, the dissipation of tides in the convective envelope of low-mass host stars and giant planets can modify the spin of the tidally perturbed body, the orbital period, and the spin-orbit angle of the perturber (e.g. Hut 1980; Ford & Rasio 2006; Lai 2012; Bolmont & Mathis 2016; Damiani & Mathis 2018). Inertial waves, which are driven by tidal forcing and restored by the Coriolis acceleration, are an important source of tidal dissipation in stellar (Ogilvie & Lin 2007; Barker & Ogilvie 2009; Bolmont & Mathis 2016) and planetary convective zones (Ogilvie & Lin 2004), where the action of turbulent motions on tidal flows is most often modelled as an effective frictional force or a viscous force with an effective viscosity that is much larger than the molecular viscosity (e.g. Zahn 1966, 1977; Duguid et al. 2020). For coplanar and circular systems, inertial waves are excited so long as the companion orbits beyond half its co-rotation radius (the orbit where the host’s rotation frequency is equal to the mean motion). Low-mass stars from K to F spectral type and giant gaseous planets both harbour a convective envelope surrounding a radiative and a solid (or diluted) core, respectively (e.g. Kippenhahn et al. 2012; Debras & Chabrier 2019). In these objects, inertial waves then propagate in a spherical shell and do not form regular normal modes of oscillation as they do in spherical and ellipsoidal geometries (Greenspan 1969; Bryan 1889, respectively). In contrast, they can focus on limit cycles, also called attractors of characteristics (Maas & Lam 1995), that are confined within the convective envelope (see also Rieutord & Valdettaro 1997). With a non-zero viscosity, attractors take the form of shear layers where the tidal wave’s energy and angular momentum can be deposited via viscous dissipation (Rieutord et al. 2001). Furthermore, viscous dissipation across shear layers can be more important as viscosity is weaker, as demonstrated notably by Ogilvie & Lin (2004) and Auclair Desrotour et al. (2015). In that respect, the tidal dissipation of inertial waves can compete with the dissipation of gravito-inertial waves in the radiative core, or it can be greater by several orders of magnitude than the dissipation of equilibrium tidal flows in the convective zone (i.e. the non-wave-like fluid’s response; see, e.g. Ogilvie & Lin 2007). The dissipation of tidally forced waves can have a great impact on the orbital and rotational evolution of the system (Auclair-Desrotour et al. 2014; Bolmont & Mathis 2016; Gallet et al. 2018; Benbakoura et al. 2019). Moreover, the dissipation of the stellar dynamical and equilibrium tides varies significantly along the evolution of the star and is highly dependent on stellar parameters, such as the mass, the angular velocity, and the metallicity of stars (Mathis 2015; Gallet et al. 2017; Bolmont et al. 2017). This makes it desirable to include all stellar processes on tidal interaction, in particular differential rotation.

The frequency-averaged tidal dissipation is often used to quantify the response of a body subject to tidal perturbations (Ogilvie & Lin 2004; Jackson et al. 2008). Yet, the dissipation of a tidally forced inertial wave is strongly correlated with the presence of an attractor at a specific eigenfrequency of the spherical shell (see Ogilvie 2009; Rieutord & Valdettaro 2010). Tidal dissipation at a given frequency may then alter each orbital and spin element of the two-body systems differently, as postulated, for instance, by Lai (2012) to explain the survival of hot Jupiters with completely damped spin-orbit angles; this idea was revisited by Damiani & Mathis (2018) with an improved treatment of dynamical tides in the convective region. Additionally, in the context of the Jupiter and Saturn moon systems, Fuller et al. (2016) and Luan et al. (2018) also investigated the dependence of tidal dissipation on frequency to explain the rapid outward migration of the moons through the resonant locking of tidally forced internal modes in the giant gaseous planets. This concept could, for example, explain the high dissipation observed in Saturn as derived from astrometric measurements at the frequency of Rhea (Lainey et al. 2017) and at the frequency of Titan (Lainey et al. 2020).

Furthermore, the fact that all layers in a star or a planet do not rotate at the same speed (i.e. differential rotation) is rarely taken into account in the determination of tidal dissipation. Yet, differential rotation seems ubiquitous in low-mass stars and giant gaseous planets. The Sun’s surface rotates in ~ 25 days at the equator versus ~35 days near the poles, and a latitude-dependent rotational gradient has also been observed in the Sun’s convective envelope thanks to helioseismology (Schou et al. 1998; Thompson et al. 2003). Through asteroseismology, latitudinal shears have been found to be comparable to that of the Sun for Sun analogues (Bazot et al. 2019), and they can be even larger for solar-like stars (Benomar et al. 2018). Essentially, differential rotation in low-mass stars depends on the effective temperature (Barnes et al. 2005, 2017) and seems to be more important when the convective envelope is thinner. Solar-like rotation profiles and anti-solar-like rotation profiles (with faster poles and a slower equator) are expected in G- and K-type stars depending on their rotation rates, based on three-dimensional numerical simulations (see in particular Brun et al. 2017; Beaudoin et al. 2018), while cylindrical rotation profiles are expected for fast rotators (Gastine et al. 2013). Regarding giant gaseous planets in our Solar System, the extent of zonal winds, which are visible on their surface as bands that run lengthwise, has been recently constrained by the probes Cassini and Juno. They extend to a depth of 3000 km for Jupiter (Kaspi et al. 2017) and penetrate down to 9000 km in Saturn (Galanti et al. 2019). Thus, the outermost molecular convective envelopes (Militzer et al. 2019; Debras & Chabrier 2019) are the seat of cylindrical differential rotation.

The study of the impact of differential rotation on the propagation and dissipation properties of inertial modes of oscillation began with the work of Baruteau & Rieutord (2013). They examined the impact of either a shellular (radial) or a cylindrical rotation profile on free inertial waves in an incompressible background by means of a Wentzel-Kramers-Brillouin-Jeffreys (WKBJ) linear analysis for an inviscid fluid and by solving the linearised hydrodynamics equations for a viscous fluid via a spectral code. Their linear analysis highlighted major differences relative to the case of solid-body rotation. Two regimes of propagation have been found in which inertial modes of oscillation can develop along curved paths of characteristics in the entire convective shell (which the authors named ‘D modes’) or in a restricted region of the convective shell, located between a turning surface and one of the shell’s boundaries (DT modes). Compared to solid-body rotation, the frequency range of the propagation of inertial modes is broader. Baruteau & Rieutord (2013) also pointed out the strong dissipation of wave energy at co-rotation resonances where the Doppler-shifted wave frequency vanishes within the fluid. All these new properties have been retrieved by Guenel et al. (2016a), who in turn examined a conical (latitudinal) rotation profile, which is typical of low-mass (F- to K-type) stars. They also confirmed the existence of unstable inertial modes (i.e. modes with positive growth rates) at co-rotation resonances, which were only found for shellular rotation in Baruteau & Rieutord (2013). Tidal forcing of inertial waves with conical rotation was introduced by Guenel et al. (2016b) within a linear numerical exploration, which also highlighted the strong dissipation of inertial waves at co-rotation resonances, particularly at low viscosities. Favier et al. (2014) also studied tidally forced inertial waves, but through non-linear numerical simulations. Differential rotation was triggered in their simulations by tidal waves depositing energy and angular momentum in an initially uniformly rotating spherical shell. In some cases, they observed hydrodynamical shear instabilities when the Ekman number (the ratio between the viscous and Coriolis accelerations) is sufficiently small.

Understanding how inertial waves interact with co-rotation resonances is thus a key issue in quantifying tidal dissipation, especially since waves may deeply interact with the background flow at this particular location, which in turn may alter the background flow (as was proposed first by Eliassen & Palm 1961, for terrestrial mountain waves). In binary systems and for late-type stars, Goldreich & Nicholson (1989) showed that the angular momentum transported by gravity waves and exchanged at co-rotation can lead to the successive synchronisation of the layers, from the base to the top of the radiative envelope. More generally, a body of work in various domains, from astrophysical disks (e.g. Goldreich & Tremaine 1979; Baruteau & Masset 2008; Latter & Balbus 2009; Tsang & Lai 2009) to geophysical fluid dynamics (e.g. Bretherton 1966; Yamanaka & Tanaka 1984), has tried to understand the properties of wave propagation and dissipation around co-rotation and, more generally, at all special locations in fluids that correspond to singularities in the linear wave propagation equation. We will refer to them as ‘critical levels’ in the following (Maslowe 1986), or as ‘critical layers’ in the case of a viscous medium. This distinction is analogous to the distinction between shear layers and attractors of characteristics that are kinds of singularities for the governing equation of inertial waves in a spherical shell. The aforementioned singularities can act very differently, with either severe absorption at the critical level (as in Booker & Bretherton 1967, for stratified vertical shear flows) or no attenuation if the wave propagates in a peculiar direction (Jones 1967; Acheson 1972; Grimshaw 1975a, for stratifiedvertical shear flows with rotation and magnetism). In other cases, a critical level may even give rise to wave amplification under certain conditions related to the first and second derivatives of the mean flow velocity (Lindzen & Tung 1978; Lindzen & Barker 1985, for barotropic and stratified shear flows, respectively). These studies all used an invariant quantity (the Reynolds stress or thewave action for rotating or magnetic flows) as a diagnostic tool to interpret the role of the critical level in terms of energy transmission and to quantify exchanges between the wave and the mean flow (Eliassen & Palm 1961; Bretherton 1966).

In light of these various studies, it is necessary to carefully consider co-rotation in differentially rotating convective zones. A local model can notably provide us a detailed understanding of physical processes at critical levels. While the propagation through a critical level of gravito-inertial waves in stratified shear flows and of Rossby waves in baroclinic and barotropic flows has been widely studied in the past decades, the behaviour of inertial waves in a latitudinal sheared flow with critical levels has so far been poorly investigated (e.g. Lindzen 1988, for a review). This is why we develop in this work a local Cartesian shearing box model to understand the complex interplay between tidal waves and zonal flows near critical levels. The concept of a shearing box for tidal flows was introduced by Ogilvie & Lesur (2012) to investigate the interactions between large-scale tidal perturbations and convective motions. In our model, we focus on the latitudinal differential rotation of the meanflow, varying the box orientation to model either cylindrical or conical rotation. The behaviour of free inertial waves in this framework is then examined near critical levels using both analytical and numerical approaches.

This paper is organised as follows. In Sect. 2, we describe the local shear model with its main assumptions and the system of governing equations. In Sect. 3, we establish a second-order ordinary differential equation (ODE) for the latitudinal perturbed velocity, and we derive the propagation properties of inertial waves for an inviscid fluid. This ODE is solved near each critical level for both conical and cylindrical rotation profiles, and we interpretenergy flux exchanges between the waves and the mean flow. We use, in Sect. 4, a three-layer numerical model to test our analytical predictions at critical levels. Frictional damping is included, and non-linear mean flow profiles are also used. Astrophysical applications with implications for low-mass stars hosting close exoplanets and for giant gaseous planets in our Solar System are discussed in Sect. 5. In Sect. 6, we summarise the main results of the paper and discuss some perspectives and caveats.

thumbnail Fig. 1

Sketch of the local Cartesian box in the convective region of a low-mass star or giant planet. Global spherical coordinates, including the depth r0, the inclination of the box θ0, and the co-latitude θ of a point of interest M inside the box, are shown to facilitate the analogy between the spherical and the Cartesian geometries.

Open with DEXTER
Table 1

Correspondence between local and global coordinate systems.

2 Local Cartesian model including differential rotation

2.1 Presentation of the model

The local model takes the form of an inclined sheared box, centred at a point C of a convective shell, as illustrated in Fig. 1. The inclined box model has already been used by Auclair Desrotour et al. (2015) to analytically characterise the properties of tidal gravito-inertial waves in the presence of viscous and thermal diffusion in stably stratified or convective regions, as well as by André et al. (2017) in layered semi-convective regions in giant planet interiors (see also Jouve & Ogilvie 2014, for two-dimensional numerical simulations of inertial wave attractors). The local coordinate system (x, y, z) corresponds to the local azimuthal, latitudinal, and radial directions of global spherical coordinates, respectively, as presented in Table 1. The mean flow velocity U is directed along the local azimuthal axis ex (we neglect possible meridional flows), and differential rotation is embodied by a latitudinal shear y U. As the box is tilted by an angle θ0 relative to the rotation axis, the rotation vector in the local coordinate system is (1)

where Ω0 is the rotation frequency of the star at the pole and and f are the normalised horizontal and vertical Coriolis components, respectively. It should be noted that the inclusion of both of these components means that we have gone beyond the traditional f-plane approximation (see also, Gerkema et al. 2008). Furthermore, we made several hypotheses to model wave propagation in a latitudinal shear flow. The buoyancy acceleration was kept in the fluid equations for the background flow. The effective gravity acceleration g also includes the centrifugal acceleration, and the fluid’s angular velocity is assumed to be small compared to the critical angular velocity , where G, M, and R are the gravitational constant, the mass, and the radius of the body, respectively. Thus, the geometry of the body is close to spherical. Furthermore, the vector g is supposed to be uniform and constant in the whole box. This requires that the typical length of the box L satisfies LHp, where Hp = −p(dz∕dp) is the vertical pressure scale height, with p being the pressure. We can assume this because tidally excited waves are expected to have small-scale structures (Ogilvie & Lin 2004; Rieutord & Valdettaro 2010; André et al. 2017). Moreover, the dimensions of the box were chosen to be small compared to the depth of the convective envelope so as to remove curvature effects.

2.2 Mean flow profile

In global spherical geometry, the mean flow based on a conical rotation profile Ω(θ) is written (e.g.in Guenel et al. 2016a) as: (2)

where eφ is the azimuthal unit vector and r and θ are the radius and co-latitude, respectively. We introduce u0 = rsinθΩ0eφ, the mean flow at a point M inside the box (see Fig. 1), without differential rotation. We also use the shear contrast δΩ = Ω(θ) − Ω0 (i.e. the difference between the angular frequency at co-latitude θ and at the pole). The shear contrast is positive for the Sun since the equator rotates faster than the pole, and negative for anti-solar-like rotating stars. Using the notations of Fig. 1, the centre C of the box is located at a distance r0 sinθ0 from the rotation axis. Accordingly, the latitudinal coordinate of the point M in the local frame is (3)

It should be noted that the radial coordinate r of the point M in spherical geometry can be written as r = r0 + z. Nevertheless, we neglected vertical displacements in the expression of the local shear because we are interested in how the (one-dimensional) horizontal shear affects the wave dynamics, contrary to many studies on differential rotation in stars that have focused on the vertical shear (e.g. Mathis et al. 2004, 2018; Decressin et al. 2009; Alvan et al. 2013). Since yr and thus θ0θ are small, we provide the correspondences in terms of mean flows and shears between the two geometries in Table 1.

As an example, the shear contrast from solid-body rotation used by Guenel et al. (2016a) was: (4)

where χ is the magnitude of the shear between the equator and the pole. Performing a second-order Taylor expansion around a fixed co-latitude θ0, such that θ = θ0yr0, and at a specified depth r0 inside the convective region, the local mean flow U can be recast as

We point out that the Taylor expansion must be pushed further at the pole θ0 = 0 (and at the pole θ0 = π with an opposite sign): (6)

Accordingly, we can approximate a conical shear as a linear mean flow at the first order when the box is tilted. We recall that conical shear has been observed in the solar convective zone and is expected in slowly and moderately rotating solar-like stars (we refer the reader to Sect. 5.1 for a detailed discussion; see also Brun et al. 2015; Beaudoin et al. 2018; Benomar et al. 2018; Bazot et al. 2019). When the box is at the pole, y becomes the distance from the rotation axis (hereafter the ‘axial distance’). Thus, the mean flow mimics a cylindrical differential rotation that can be modelled using a cubic y-profile given by Eq. (6). This rotation profile is found in Jupiter and Saturn, as well as in rapidly rotating stars, as demonstrated, for instance, by Gastine et al. (2013) and Brun et al. (2015).

2.3 System of equations

To derive the system of governing equations for tidal waves in the local reference frame, we made several hypotheses. Stratification terms, which usually drive the propagation of internal gravity waves, were kept for the sake of clarity and will be methodically kept or removed after applying the Boussinesq approximation and setting the equations for inertial waves. Moreover, we assume that the action of turbulence can be modelled as a Rayleigh friction term in the momentum equation with an effective frictional damping rate σf. This simplifies the analytical solution of the fluid equations compared to the usual modelling of turbulence as an effective viscous force (see in particular Ogilvie 2009). The momentum, continuity, and thermodynamic equations for tidal waves in a differentially rotating Cartesian framework are thus:

where u, p, ρ, and f denote the velocity, pressure, density, and volumetric tidal forcing, respectively. We have also introduced cs, the sound speed, and is the total derivative operator.

All variables are then linearised at the first order: Zero-order terms correspond to background equilibrium quantities, and first-order terms represent the leading perturbation. The local velocity, density, and pressure are therefore written as: (10)

where u = (u, v, w) in the local Cartesian basis. We have introduced the dimensionless parameter ϵ: (11)

where we use 1∕(2Ω0), a characteristic time scale, and L, a characteristic length scale of the mean flow. These notations are based on those of Grimshaw (1975a) and adapted to our model. In the following, we will work with dimensionless variables using the above scaling, including 2Ω0 L to scale velocity and ρTgL to scale pressure, with ρT the reference density. The dimensionless momentum equation of the mean flow is: (12)

with n the unit vector parallel to the rotation axis. Projecting Eq. (12) into Cartesian coordinates, one can derive: (13)

At the leading order in ϵ, one can recognise the hydrostatic balance and, at the first order, the geostrophic balance (the set is akin to the thermal-wind equilibrium assumption, see e.g. Grimshaw 1975a; Yamanaka & Tanaka 1984). We underline that tending ϵ to zero is similar to assuming the Boussinesq approximation. Indeed, all density variations are neglected, except the ones involved in the buoyancy force. The dimensionless Brunt-Väisälä frequency is (14)

where we introduce the dimensionless number , which is small when filtering acoustic waves. Consequently, the curl of Eq. (12) gives (15)

where we neglect the second-order terms in ϵ.

Now, we make several assumptions to treat the propagation of inertial waves. As the convective motions are essentially adiabatic, the convective zone can be assumed to be neutrally stratified to a first approximation. Hence, the Brunt-Väisälä frequency N is cancelled out in the third density relationship of Eq. (15). Moreover, we make the Boussinesq approximation, which means that we neglect terms in ϵ and F in the final set of perturbed equations. Thus, the dimensionless linearised momentum, continuity, and thermodynamic equations are ultimately:

where v is the latitudinal velocity perturbation. We emphasise that, although vertical stratification has been filtered in the limit where N goes to zero, a horizontal stratification term remains in Eq. (18). As a result, we consider the inertial waves propagating in the inclined shear box where the mean flow is maintained by the thermal-wind balance.

2.4 Equilibrium state of the background flow

It is worthwhile discussing the choice of keeping buoyancy forces in the zero-order momentum equation. Without gravitational forces, the momentum equation for mean dimensional variables is written as a geostrophic balance: (19)

This balance satisfies the Taylor-Proudman theorem (Rieutord 2015), namely the geostrophic flow is independent of the coordinate parallel to the rotation axis. When taking the x-axis (the only non-zero) projection of the curl of this equation, one obtains the following relationship: (20)

Without the vertical stratification embodied by the Brunt-Väisälä frequency, nor latitudinal stratification, the equilibrium of a y-dependent mean flow is thus not ensured for an incompressible fluid. An alternative to conserve the equilibrium without stratification would be to consider a z-dependence of the mean flow. Such a possibility is not considered in this paper since we are mainly interested in latitudinal mean flow profiles. Furthermore, in addition to maintaining differential rotation, the latitudinal stratification can allow for the construction of an invariant that is useful for studying energy transfer at critical levels: the wave action flux. This will be discussed further in Sect. 3. Lastly, since at the poles, the latitudinal stratification term will not appear in the perturbed fluid equations (as we can see from Eq. (18)).

3 Dynamics of inertial waves at critical levels: analytical predictions

In this section, we analytically investigate the behaviour of inertial waves at critical levels in a non-dissipative fluid at various co-latitudes. For this purpose, we consider perturbations q in the normal mode (21)

where ω is the complex inertial frequency, kx and kz are the real streamwise and vertical wavenumbers, respectively, and c.c. is the complex conjugate.

3.1 Wave propagation equation in the latitudinal direction

Using the modal form Eq. (21) for ρ, p, and u, we solved the set of hydrodynamic equations, Eqs. (16)–(18), for the latitudinal velocity v. Considering free inertial waves (i.e. without forcing terms), the set of perturbation equations can be recast into a single second-order ODE for v: (22)

where the prime now denotes the derivative according to y, and A, B, and C are the coefficients that can be simplified without friction as follows: (23)

where is the absolute wavenumber in the direction perpendicular to the y direction and σ = ωkxU is the (dimensionless) Doppler-shifted wave frequency. We refer the reader to Appendix A for the detailed ODE derivation with friction and tidal source terms. Equation (22) becomes singular when A = 0 or σ = 0, and these singular points are called critical levels (see e.g. Bretherton 1966; Grimshaw 1975a). The critical level where the Doppler-shifted frequency equals zero (i.e. σ = 0) can be met when the mean flow matches the local phase velocity, and this is also known as ‘co-rotation resonance’ (e.g. in Goldreich & Nicholson 1989; Goldreich & Tremaine 1979; Ogilvie & Lin 2004). When the Coriolis acceleration is not taken into account (i.e. when treating internal gravity waves), the co-rotation resonance is the unique critical level (see e.g. Booker & Bretherton 1967). At co-latitudes other than the poles, the critical levels come in three flavours, the co-rotation σ = 0 and two other critical levels that are defined, in our model, by (we recall that is the latitudinal component of the rotation vector). These critical levels were similarly reported for vertical shear flows, as in the studies of Jones (1967) and Grimshaw (1975a) for vertical and inclined rotation vectors, respectively. In these works, the Doppler-shifted frequency at critical levels other than the co-rotation resonance equals ± 2Ωv, where Ωv is the vertical component of the rotation vector.

3.2 Propagation properties

3.2.1 Dispersion relation, group, and phase velocities

The three-dimensional dispersion relation is a fourth-order equation in Doppler-shifted frequency when injecting wave-like solutions into the three directions, x,  y, and z, in Eq. (22). In order to understand the main properties of waves at the critical level, we made the short-wavelength approximation in the meridional plane, as in Baruteau & Rieutord (2013) and Guenel et al. (2016a). This involves keeping only the second-order derivatives in the y and z directions, and it reduces the relation dispersion to a second-order equation when injecting plane wave-like solutions. In the local meridional plane, the differential equation reduces to a Poincaré-like equation: (24)

where we recover the Poincaré equation (for the propagation of inertial waves in the inviscid limit; Cartan 1922) in the meridional plane when there is no shear (U′ = 0) and at the poles ( and f = 1). Moreover, we set v ∝ exp−i(kzzkxx) so as to write the wave dispersion relation for the Doppler-shifted frequency σ: (25)

where is the norm of the wave vector in the meridional plane (e.g. for fixed kx), as in Baruteau & Rieutord (2013). Compared to solid-body rotation (see e.g. Rieutord 2015), an additional term () is present, which accounts for the latitudinal shear. Assuming that σ2 takes positive values (as in Baruteau & Rieutord 2013; Guenel et al. 2016a), we therefore introduce (26)

We can then specify the phase velocity in the meridional plane: (27)

In the same way, we can derive the expression for the group velocity in the meridional plane: (28)

We note that without differential rotation, the group velocity reduces to its well-known expression for solid-body rotation (e.g. see Rieutord 2015): (29)

Moreover, as in solid-body rotation, the group velocity (Eq. (28)) and the phase velocity (Eq. (27)) lie in perpendicular planes: vgvϕ = 0.

When the box is located at the north pole (θ0 = 0 in Fig. 1), by setting κ = 1 − U′, we recover (30)

as in Latter & Balbus (2009) and Baruteau & Rieutord (2013), where κ can be identified with the epicyclic frequency of Baruteau & Rieutord (2013) and ky corresponds to the cylindrical component of the wavenumber (ks).

3.2.2 Phase and group velocity at singularities

In this section, we derive the conditions required to meet singularities in terms of wavenumbers and shear, and we examine the implications for the phase and group velocities. When the box is inclined, for σ → 0 we must have γ → 0, meaning that vϕ0 while  |vgey|→ and |vgez|→ (possibility 1).

Guenel et al. (2016a) found similar results by studying the propagation of free inertial waves in a global frame with conical shear, namely, when their parameter (which is homogeneous to a frequency and equivalent to our γ parameter) goes to zero, the group velocity goes to infinity while the phase velocity is cancelled out. According to their work, an inertial wave may propagate across the co-rotation.

Now to get , we either need |ky|→ at fixed kz, which implies vϕ0 and vg0 and means that inertial waves cannot get through the critical level (possibility 2), or |kz|→ 0 at fixed ky, which gives vϕez → 0 and vgey → 0, while and |vgez|→ fky; the wave may then cross the critical level with some preferential direction (possibility 3).

Again, these conditions share some similarities with those observed for co-rotation in a global spherical geometry. The first of the three possibilities is analogous to the global phase and group velocities tending to zero when ks, with ks the axial wavenumber in cylindrical coordinates (Baruteau & Rieutord 2013; Guenel et al. 2016a). This makes sense since here the axial distance is s = rsinθ and y ~ r0(θ0θ). However, the second possibility is slightly different from both these previous works, in that |kz | → 0 at fixed ks, where kz is the global vertical wavenumber along the rotation axis, unlike our local vertical wavenumber kz which is along the spherical radial coordinate.

We point out that the singularities at arise in our model because the rotation vector is inclined with respect to the local vertical axis of the box. In the global model of Guenel et al. (2016a), three conditions for a wave to meet the co-rotation exist, and these conditions are actually quite similar to the three above possibilities for waves in our model to interact either with the co-rotation σ = 0 or with the other critical levels at . Hence, the local critical levels at behave somewhat similarly to the co-rotation in the global framework, as if we partially broke the degeneracy in the local framework of the origin of the co-rotation found in the global framework.

When the box is at the north pole, the three possibilities to meet co-rotation are similar but lead to different relationships for the phase and group velocities: (i) κ → 0 (i.e. U′→ 1), meaning that vϕ0 and vg0: The wave is totally absorbed at co-rotation. (ii) |ky|→ at fixed kz, which implies vϕ0 and vg0: We can draw the same conclusion as in the previous case and this case is also analogous to the second possibility when the box is tilted. (iii) |kz | → 0 at fixed ky, which gives vϕ0, vgey → 0 while |vgez|→ κky: the wave energy does not cross the co-rotation in the latitudinal direction (equivalent to the vertical paths of characteristic in global cylindrical geometry, as in Baruteau & Rieutord 2013). At the north pole1, we actually have a perfect match with the conditions given by Baruteau & Rieutord (2013) when using a cylindrical rotation profile for the mean flow.

3.2.3 Energetical aspects

In this section, we examine the energetic balance associated with inertial waves in our inclined shear box model, without assuming the short-wavelength approximation. This energetic balance does not include potential energy because of the adiabaticity of the convective region, but two additional terms relative to the solid-body rotation case appear, which come from the differential rotation. We denote the displacements along the vertical and latitudinal directions with η and ζ, respectively. Considering that ω = kxc, where c is the longitudinal phase velocity (e.g. as in Booker & Bretherton 1967), we can use the first-order definition (31)

This first allows us to express the perturbed density from Eq. (18) as (as a reminder, the symbol ′ has been dropped for perturbed quantities). Then, by multiplying the momentum equation, Eq. (16), by ρ0 u, we get the energy balance equation: (32)

where ek = ρ0u2∕2 is the kinetic energy density and pu is the so-called acoustic flux. We now integrate the above energy balance equation over x and z, as well as over one wave period as the perturbed quantities have a wave-like form in these directions. Further assuming that the box is δ thick in the y direction, the energetic balance yields: (33)

where we introduce, from left to right, the power of the external pressures at the boundaries − δ∕2 and δ∕2 on the perturbed latitudinal flow, the work of the shear, the viscous dissipation, and the forcing power, which read, respectively: (34)

where the bar represents the average in the (z, x) plane over one period. We note that the energy density and the acoustic flux in the x and z directions drop out in Eq. (33) when integrating because of the wave periodicity in those directions. The quantity can also be seen as the power transferred from the mean flow to the perturbation (or conversely) by the Reynolds stress: (35)

where we use partial integration and the periodicity of perturbations in the x and z directions. At the pole, , so we recover the definition of the Reynolds stress in Miles (1961), who studied the stability of a two-dimensional stratified y-sheared flow (i.e. ). This quantity can also be called the latitudinal flux of horizontal momentum (in the (z, x) plane sense) in reference to the vertical flux of horizontal momentum in stratified z-sheared flows. Moreover, we emphasise that the latitudinal flux of energy is not conserved, even in the inviscid free-wave problem. This is due to the presence of the shear, as already stated for example by Eliassen & Palm (1961), who studied stratified vertically sheared flows. They underline that when the mean flow varies with height, the kinetic energy of the mean motion can be converted into wave energy. Without friction and forcing, the y-derivative of the latitudinal flux is: (36)

Using the same method as Broad (1995), we multiplied the x-projection of the inviscid force-free momentum equation by ζ: (37)

By multiplying by (Uc), the latitudinal flux of energy can thus be written as: (38)

By differentiating this relationship with respect to y, and by equalising with Eq. (36), one can obtain: (39)

that is with τ the Reynolds stress. Equation (39) is naturally satisfied at co-rotation, where Uc = 0, or if the Reynolds stress is uniform. Booker & Bretherton (1967) have shown that the Reynolds stress is discontinuous at a critical level, highlighting exchanges between wave energy and the mean flow. Compared to the analysis of Broad (1995) for three-dimensional stratified shear flows, Eq. (39) is not vectorial, because our base flow is unidirectional.

3.2.4 Polarisation relations

For the forthcoming analysis, it is useful to derive expressions of the perturbed projected velocities and the perturbed reduced pressure2 Π = pρ0, namely the polarisation relations (see Appendix B for more details). In the inviscid free-wave problem, these perturbed quantities can be written in terms of the latitudinal velocity, its derivative, and the shear: (40)

Without shear and at θ0 = 0, we recover the polarisation relations in the solid-body rotation case (see e.g. Rieutord 2015).

3.2.5 Conservation of the wave action flux

While the latitudinal flux of energy is not conserved in the whole domain, there is a conserved quantity, called the wave action flux, as introduced in Grimshaw (1975a): (41)

which is the latitudinal flux averaged over vertical and longitudinal wavelengths divided by the Doppler-shifted frequency. A general treatment for the derivation of the wave action as a conserved quantity can similarly be found in Andrews & McIntyre (1978). The wave action flux is related to the Reynolds stress τ as . By using the expression for the perturbed reduced pressure Π derived in the previous section, the wave action flux now reads: (42)

Unlike the latitudinal flux of energy, but similarly to the Reynolds stress, this wave action flux is conserved along the latitudinal direction. One can demonstrate that in the whole domain except at critical levels by using the expression for the reduced pressure in Eq. (40) and the ODE (Eq. (22)). Several works have shown that a properly defined (i.e. conserved) angular momentum transport parameter can be found in z-sheared mean flows without rotation (Booker & Bretherton 1967), with rotation under the traditional approximation (Jones 1967), and with rotation under the non-traditional approximation (Grimshaw 1975a). Verifying the conservation in the whole domain except at critical levels is really important because it brings to the fore energy transfers that are due to the critical levels. We specify that is a measure of wave energy through a surface (in the (z, x) plane) since is the energy density transported by the group velocity3 Vg in the latitudinal direction (e.g. Bretherton & Garrett 1968; Mathis & de Brye 2012). It should be underlined that the wave action flux has been defined in the inviscid limit and is not conserved when the friction is taken into account.

3.3 Inertial waves at critical levels when the box is tilted

In this section, we analytically investigate waves passing through the various critical levels in the tilted box. We examine the behaviour of the waves around the co-rotation σ = 0 and the critical levels when the box is tilted (for the co-rotation when the box is at the pole, see Sect. 3.4).

3.3.1 Critical levels at

In this subsection, we treat both singularities, , simultaneously. Although Eq. (22) does not have analytical solutions in general, it is still possible to study the behaviour of an inertial wave close to the critical levels defined by by approximating the ODE through its first-order Taylor expansion in the vicinity of these singularities, and then by applying the Frobenius method. We introduced y±, the location of the related critical level . For a linear mean flow profile U = Λy, with Λ a constant, y± are given by: (43)

Without any assumption on the mean flow profile, the first-order Taylor expansion of the ODE (Eq. (22)) near y± is: (44)

with (45)

where the symbol ± refers to the regular singularities4 y+ and y and is U′ evaluated at these singularities. The Frobenius method consists in injecting the power function into Eq. (44), with λ a constant to be determined (see e.g. Morse & Feshbach 1953). The corresponding indicial equation is then (46)

with solutions: (47)

Therefore, the two independent solutions of Eq. (44) can be written as follows: (48)

where an and bn are complex constants. Both solutions are valid in the vicinity of the critical level around which they are built in the complex plane, up to the next singularity if it exists. The coefficients a0 and b0 are unconstrained and depend on boundary conditions, unlike the other factors that can be determined by injecting the solutions from Eq. (48) into the linearised ODE (Eq. (22)) around y± at the right order for the desired coefficients.

Near the critical points y±, the total solution v± is well approximated by the lowest orders of v1,± and v2,±: (49)

Owing to the existence of a branch point at y± (since λ± is complex),reconnecting solutions on either side of the critical levels is not straightforward. This requires both physical and mathematical arguments (see in particular Miles 1961; Booker & Bretherton 1967; Ringot 1998). In order to remove the degeneracy of the path from positive to negative yy± (i.e. choosing either e+ or e), we made use of a complex inertial frequency ω = ωR + I, assuming the radiation condition ωI > 0. This condition ensures a non-growing wave towards infinity. The Taylor expansion of the base flow at the first order in yy± gives (50)

and, by definition, we have (51)

Consequently, the solution below the critical level is unambiguous in terms of the above solution coefficients and depends on (52)

In other words, when taking yy± to decrease from positive to negative values, its complex argument changes continuously from 0 to . Thus, the appropriate path for determining the branch of passes under(above) y± as long as () (the same reasoning can be found in Grimshaw 1975a). Therefore, the solution on both sides of the critical level y± is: (53)

The remaining issue is now to know in which direction the wave is propagating. The second part of the solution can be assimilated to a wave-like solution with the varying latitudinal wavenumber . Moreover, according to Eq. (42), the wave action flux on either side of y± is: (54)

The group velocity gives the direction towards which the energy is transported (we recall that , with Vg the group velocity and the local energy density). By consequence, for the solution featuring the coefficient b0. If kzf is positive, this wave transports energy downwards (upwards) across the critical level (). If kz f is negative, the wave transports energy upwards (downwards) across the critical level (). In all cases, the action flux of the wave with the amplitude |b0| will be transmitted (in the direction given by the sign of kzf and the critical level y+ or y) by a factor , where (55)

with αk = kzkx and , after passing through the critical level. Such a wave will always be attenuated since . The transmission factors and are displayed in Fig. 2 in terms of the absolute value of the shear Rossby number |Ro | and the ratio of wavenumbers αk. The lower the amplitude of the Rossby number and the lower the inclination, the more likely the wave is to be strongly attenuated at any αk. We reiterate that a low Rossby number refers either to fast-rotating stars or to low differential rotation. At the equator, one should note that f = 0, thus there is no transmission nor exchange of wave action flux near the critical levels y± in the inviscid limit (see Eq. (54)). Results are the same for θ0 + ∕2 with k ∈ {0, 1, 2, 3} and for negative Rossby numbers. However, it has to be emphasised that the cases where the inclination satisfies θ0 = with k ∈ {0, 1} are not well described by the attenuated factor and require a specific treatment, as discussed in Sect. 3.4.

It is important to note that, with fixed parameters {kz, θ0, Ro}, the attenuation of the wave action flux is specific to a single direction of wave propagation (i.e. the solution featuring the coefficient b0). The solution of coefficienta0 is not affected by the attenuation. This is the so-called valve effect introduced by Acheson (1972) in the context of hydromagnetic waves in a rotating fluid. It was also evidenced by Grimshaw (1975a) and further discussed in Grimshaw (1979) for magneto-gravito-inertial waves in an inviscid and compressible z-sheared fluid.

thumbnail Fig. 2

Transmission rate of a wave passing through any of the critical levels defined by as a functionof the absolute value of the shear Rossby number |Ro| = |U′| (where U′ is scaled by 2Ω0) and the ratioof the horizontal wavenumbers αk = kzkx for a co-latitude of the box θ0 = 10° (left panel) and θ0 = 80° (right panel).

Open with DEXTER

3.3.2 Inertial wave crossing co-rotation

We performed the same analysis as in the previous section to treat the co-rotation point y0 where σ = 0 (i.e. U(y0) = U0 = ωkx). The linearised ODE (Eq. (22)) near the co-rotation using the Taylor expansion of σ and U at the lowest orders is: (56)

where and are the first and second derivatives of the mean flow profile U evaluated at the critical level y0. The singularity at the co-rotation is a regular singularity, and we can again use the Frobenius method. The indicial equation has solutions λ = {2, 1}. Since the difference between the two values of the exponent λ is an integer, one expects a second independent solution v2 of Eq. (56) that includes a logarithmic part, such as (e.g. Schmid et al. 2002) (57)

with the first solution and an, bn, and L complex coefficients. However, when injecting v1 + v2 into Eq. (56), one finds L = 0, meaning that a sole polynomial solution in the form (58)

includes all the solutions of Eq. (56), with c0 = b0 and a1 = b1 + a0 determined by boundary conditions, and determined by recurrence via the expansion of Eq. (22) around y0. As a result, thewave action flux given by Eq. (42) becomes (59)

just below and above the co-rotation, and it is continuous there, similarly as in Grimshaw (1975a), but here without being restricted to a linear mean flow profile. Hence, no transfer of wave action flux is expected at co-rotation in the inviscid limit when the box is inclined relative to the rotation axis (i.e. for conical differential rotation), regardless of the mean flow profile. This result also holds true when the box is located at the equator.

As in the works of Grimshaw (1975a) and Jones (1967), it is tempting to investigate the asymptotic behaviour of a wave when y in order to better constrain the propagation of waves through one or multiple critical levels. Nevertheless, the term in the ODE (Eq. (22)), which cannot be overlooked as it was in the aforementioned studies since we do not have vertical stratification, makes the singularity y = an essential (or irregular) singularity, and the Frobenius method cannot be applied. This term also prevents us from applying an analysis such as the WKBJ approximation because, even far from critical levels, the coefficients C of the ODE (in Eq. (22)) still have a strong dependence on the latitudinal coordinate when the box is tilted.

3.4 Inertial waves when the box is at the poles

When the box is located at the north or south pole, and the ODE (Eq. (22)) is greatly simplified. For θ0 = 0, the dimensionless wave propagation equation becomes (60)

At the southpole (i.e. θ0 = π), the term 1 − U′ in Eq. (60) is replaced by 1 + U′. We note that this equation is reminiscent of the differential equation for Rossby waves in the β-plane, that is, and there is a constant df∕dy = β, with f the Coriolis parameter (e.g. Miles 1961; Grimshaw 1975b; Gliatto & Held 2020). However, we cannot make a direct comparison at co-rotation, because the singularity in the equations for Rossby waves and inertial waves is not of thesame order. We have a second-order pole around the co-rotation, while only first-order poles are found in the aforementioned studies. In fact, Eq. (60) is similar to the wave equation in stratified z-sheared flows (e.g. Jones 1968).

In our polar configuration, the y-coordinate is now the axial distance, and this means that the mean flow has a cylindrical profile. Sucha rotation profile is expected in giant planets such as Jupiter and Saturn (Kaspi et al. 2017; Galanti et al. 2019, respectively) as a natural outcome of the Proudman-Taylor theorem for fast-rotating bodies. The propagation and dissipation of inertial modes of oscillations in the presence of critical levels for this kind of mean flow have been investigated by Baruteau & Rieutord (2013) in a spherical shell.

3.4.1 Analytical solutions with constant shear

Analytical solutions of the ODE Eq. (60) are difficult to find for general profiles of the mean flow (e.g. a quadratic mean flow profile). A linear mean flow profile, on the other hand, has analytic solutions, which is why we use such a profile in this section (i.e. U = Roy, with Ro the shear Rossby number that is taken as constant here). Equation (60) then becomes: (61)

where y0 = ω∕(kxRo) and αk = kzkx is the vertical to longitudinal wavenumber ratio. When the box is located at the south pole, the left-hand term in the bracket is in the numerator. This equation takes the form of Whittaker’s equation (see Abramowitz & Stegun 1972), and solutions can be written in terms of the Whittaker functions M: (62)

where = 2k(yy0), (63)

and A and B are complex constants given by boundary conditions. The Whittaker function M0,−μ allows a quite straightforward analytic continuation: (64)

By consequence, the solution below the critical point y = y0 is: (65)

Although the Whittaker functions do not feature precisely as wave-like forms, we can already get a good idea of the attenuation factor thanks to analytic continuation, as will be shown in the following section.

It is important to point out that μ can be real or complex depending on the value of (66)

which will simply be denoted by R in the following. This can drastically change the behaviour of a wave passing through the co-rotation. A necessary, but not sufficient, condition to find an instability is that R < 1∕4, as we will demonstrate in Sect. 3.4.4. This condition is similar to the Miles-Howard theorem for stratified z-sheared flow (Miles & Howard 1964; Lindzen 1988). In these studies, the prerequisite for instability is that Ri < 1∕4, where Ri is the Richardson number, namely the squared ratio of the Brunt-Väisälä frequency to the vertical (or radial in global spherical geometry, Alvan et al. 2013) shear. In our model, unlike cases where the box is tilted, a WKBJ analysis can typically be performed for a linear mean flow, provided that |R| ≫ 1∕4, in line with the condition of stability derived in the coming sections and detailed in Appendix C.

These various situations regarding the value of at the north and south poles are illustrated in Fig. 3. We stress the particular case where Ro = 1 (Ro = −1) at the north (south) pole, and where the differential equation and its solutions take a quite simple form: (67)

Solutions are then fully evanescent for such shears. One can notice that Eq. (67) is the same far from co-rotation for any mean flow.

Finally, it is clear from Fig. 3 that wave propagation is the same at the north and south poles provided a Rossby number of opposite sign. As a result, only the equations at the north pole will be treated in the following, and the word ‘pole’ now refers to the north pole.

3.4.2 Frobenius method at the pole

Though analytic solutions are known, it is still useful to determine Frobenius solutions near co-rotation for two main reasons. First, these solutions can be derived for any mean flow profile near the co-rotation. Close to co-rotation, the mean flow is approximated by a Taylor expansion at the first-order . Secondly, Frobenius solutions may feature wave-like forms, which is helpful for physical interpretation. Therefore, Eq. (60) can be written near co-rotation as (68)

where and . The indicial equation gives: (69)

In the two next subsections, we examine the two cases where μ is imaginary or real.

thumbnail Fig. 3

Diagram showing how compares to 1∕4 for two positions of the box (θ0 = 0 and π), as well as for a range of Rossby numbers Ro and a range of the vertical to longitudinal wavenumber ratio αk. The solid and dashed black lines mark where at the north and south poles, respectively. The purple domain shows where , and the white region where , regardless of whether the box is at the north or the south pole. In the dark red region, and and vice versa in the orange region.

Open with DEXTER

3.4.3 Theoretical stable regime (R > 1/4)

We address here the case where R > 1∕4. The same analysis from Sect. 3.3 can be carried out to determine how a wave behaves upon crossing the co-rotation. The solutions of the indicial equation can be recast as (70)

The first-order solutions to Eq. (68) in the vicinity of y0 are: (71)

for . One can recover the same form of Frobenius solutions as in Alvan et al. (2013), who examined radially stratified mean flows in spherical geometry. As , the wave action flux Eq. (42) reduces to (72)

that is, injecting the solutions on both sides of the critical level: (73)

This formulation is quite similar to the expression of the Reynolds stress (τ) in vertically stratified mean flows, which can be found in Booker & Bretherton (1967) in Cartesian geometry. We recall that, indeed, . Moreover, given Eq. (35), the Reynolds stress in our model reads . Using the polarisation relations for u, we recover the wave action flux in Eq. (72).

The pre-factor i in the solutions (Eq. (71)) below the co-rotation does not affect the energy flow and simply indicates that the wave undergoes a phase shift of π∕2 through the critical level (see also Alvan et al. 2013). Above the critical level, the normalised Doppler-shifted frequency satisfies sign (σ) = − ς, as was the case for the co-rotation in the inclined box. The sign is reversed below the critical level. Thus, the first solution of main amplitude |a0 | carries its latitudinal flux of energy upwards (downwards) for ς = +1 (ς = −1), while the second solution transfers its energy in the opposite direction in the various cases. Therefore, the energy flux of an upward or downward wave is always attenuated by a factor (74)

This attenuation factor is shown versus Ro and αk in Fig. 4. We observe that the wave is largely absorbed at the critical level and thus deposits most, if not all, its energy for most couples (αk, Ro).

thumbnail Fig. 4

Transmission rate TF of the wave action flux across the co-rotation when the box is located at the pole for R > 1∕4. It is displayed against the Rossby number Ro and the ratioof the horizontal wavenumbers αk = kzkx. The forbidden region where R < 1∕4 is shown in white, and the dark red cone corresponds to values of TF that are lower than 10−9.

Open with DEXTER

3.4.4 Possible unstable regime (R < 1/4)

We now deal with the case where R < 1∕4 (i.e. μ is real). Contrary to the situation where R > 1∕4, we can no longer assimilate solutions to wave-like functions. The exponential form of solutions for R < 1∕4 near the critical level reads (75)

and makes this region fully evanescent. Furthermore, the associated wave action flux is (76)

Without knowing the direction of the wave or the energy flux, since , it is difficult to assess the impact of the critical level on wave propagation, that is, whether it will attenuate waves or, on the contrary, amplify them.

Lindzen & Barker (1985) found a way to investigate the behaviour of internal gravity waves in the presence of a vertical shear, passing through a critical level in a regime similar to ours (Ri < 1∕4) where solutions are fully evanescent. Their work, which was carried out in local Cartesian geometry, was taken up by Alvan et al. (2013) in global spherical geometry applied to the radiative zone of solar-like stars and evolved stars. The goal of the method is to determine the reflection and transmission coefficients in a three-zone model. The evanescent region where the Richardson number satisfies Ri < 1∕4 and where the critical level is located (zone II) is sandwiched between two propagating wave layers (zones I and III). Using a linear mean flow profile so as to establish solutions inside zone II, Lindzen & Barker (1985) and Alvan et al. (2013) both used continuity relations of the perturbed vertical (or radial) velocity and its derivatives at the interfaces between zones in order to get the transmission and reflection coefficients. The critical level is located in the middle of zone II which has a width δ. By consequence, the reflection and transmission coefficients depend, in their works, on the shear and, more precisely, the Richardson number, as well as on the width δ. They both found that, depending on the Richardson number and δ, the reflection and transmission coefficients can be greater than one, meaning that the wave can be over-reflected and/or over-transmitted and thus extract energy and angular momentum fluxes from the mean flow, which can, in turn, lead to potential shear instabilities after successive encounters of the wave with the critical layer. However, this result is conditioned by the geometry of the model. As shown by Lindzen (1988) in his review, and references therein, models with one or even two layers with evanescent waves, and potentially a wave-like region, do not allow such phenomena. A first region that allows the wave propagation is mandatory and is combined with a ‘sink’ that pulls the wave to cross the critical level. According to Lindzen & Barker (1985) and Lindzen (1988), the nature of the sink for wave flux can be either another propagative region or an evanescent region that is, as in zone II, subject to friction processes. Given this peculiar geometry, instabilities can occur under boundary conditions that allow the wave to successively return to the critical level. Many studies have tried to relate over-reflection and shear instability for a specific wave geometry (see in particular the reviews of Lindzen 1988; Harnik & Heifetz 2007, for internal gravity waves and Rossby waves).

In the present study, we do not further investigate shear instability by doing, for instance, a temporal analysis to estimate the waves’ growth rate (as in Lindzen & Barker 1985; Watts et al. 2004, who considered an initial value problem). On the contrary, we give arguments of necessary but not sufficient conditions to find instabilities, such as R < 1∕4. It is important to note that R is constant in the whole domain for a linear mean flow profile, and thus one is stuck with either a propagative (stable) or an evanescent regime. Therefore, finding an adequate geometry to allow over-reflection and over-transmission requires at minimum that the Rossby number is not the same in the whole domain by using, for instance, a non-linear mean flow profile. Furthermore, in the particular case where R = 0 (i.e. Ro = 1 or αk = 0 in Eq. (60) when the wave with kz = 0 propagates in the (x, y)-plane), a necessary condition for instability is given by the inflection point theorem (Schmid et al. 2002). This theorem is particularly used to study barotropic instabilities for Rossby waves (see e.g. Lindzen & Tung 1978). In other words, a necessary condition to have unstable modes for Ro = 1 is that U″ is cancelled out in the domain of wave propagation.

We summarise in Table 2 the main analytical results of Sects. 3.3 and 3.4 regarding wave and wave action flux transmission, when the box is tilted relative to the rotation axis with a random angle, at the north pole, and at the equator, in the inviscid limit. We note that when ‘no’ is given in both the attenuation and amplification columns, the wave is fully transmitted across the critical level, regardless of the wavenumbers and of the mean flow profile.

Table 2

Summary of analytical results at each critical level when the box is inclined or not with respect to the rotation axis, in the [0, π∕2] quadrant, and for positive wavenumbers.

4 A three-zone numerical model

In order to test the analytical predictions of the previous section, we built up a three-zone numerical model to simulate waves passing through critical levels. A similar model has been used, for instance, by Jones (1967) to explore the behaviour of internal gravity waves passing through critical levels in a fluid with rotation and vertical shear. In our model, we solve the two first-order ODEs satisfied by v and Π, the combinationof which led to the wave propagation equation, Eq. (22). By imposing boundary conditions such that waves satisfy the dispersion relations (see also Appendix D.1), we examine the dynamics of inertial waves propagating in the shear region. Also, whenever possible, we analytically calculate the wave transmission and reflection coefficients as the wave-like solution crosses the shear region.

thumbnail Fig. 5

Mean flow profiles used in the three-zone numerical model against y. In the no-shear regions I (y < 0) and III (y > 1), the mean flow is uniform and set respectively to U = 0 and U = Λ. The shearedregion II (grey-shaded) can have a linear, square, or cubic mean flow profile.

Open with DEXTER

4.1 Description of the model

The mean flow profile that is used in the three-zone model is illustrated in Fig. 5. The zone with shear (zone II) is surrounded by two no-shear regions, one with no mean flow (zone I) and one with a uniform mean flow (zone III). In the whole domain, the mean flow profile that we adopt is expressed as (77)

where U(y) is continuous at each interface and n is an integer: n = 1 for a linear shear flow, n = 2 for a square shearflow, or n = 3 for a cubic shear flow (see also Fig. 5). In zone I, we assume that there is an incident wave that enters the shear zone as well as a wave that is reflected at the interface between zones I and II or in zone II, that is: (78)

where AI and AR are the amplitudes and kI and kR are the wavenumbers of the incident and reflected waves, respectively. We further imposed as a boundary condition in zone III a transmitted wave that propagates towards positive y-values: (79)

where AT and kT are the amplitude and wavenumber of the transmitted wave, respectively. We imposed AT = 1 without loss of generality and computed the remaining amplitudes AI and AR. More detailson the solutions and the dispersion relations of the waves in zones I and III can be found in Appendix D.1. We ensured that the transmitted wave carries energy upwards by deriving the wave action flux in zones I and III (see Appendix D.2).

We imposed the continuity of the latitudinal velocity and reduced pressure at the interfaces (at y = 0 and y = 1). By doing so, the wave action flux is continuous at both interfaces. Thus, in the absence of critical points, the wave action flux is conserved in the whole domain, namely , where is the wave action flux of the transmitted wave, and , where and are the incident and reflected wave action fluxes, respectively.

To solve the ODE in the three zones, and in particular near singularities, we used MATLAB’s solver ode15s, which is suitable for solving stiff differential equations (Shampine & Reichelt 1997). To avoid strict singularities at , we added a small friction, σf = 10−8, to our set of units. Given the boundary conditions, the numerical solver deals with two first-order ODEs for v and Π, which take the form (80)

where (81)

and where we recall that s = σ + f is the modified Doppler-shifted frequency due to Rayleigh friction. While AT = 1 is imposed by the boundary condition, we computed AI and AR by comparing numerical solutions of the system Eq. (80) at y = 0 with the definition of velocity in zone I (Eq. (78)) and its associated reduced pressure (see Eq. (D.9)).

4.2 Numerical exploration at the pole for a constant shear

4.2.1 Reflection and transmission coefficients

In most cases, for any inclination of the box and any mean flow profile, there is no analytical solutions in zone II. Nevertheless, we have shown in Sect. 3.4.1 that, when the box is at the pole and for a linear mean shear flow, solutions can be found in terms of Whittaker functions. In this section, we will find reflection and transmission coefficients similarly to how it was done in Lindzen & Barker (1985) and Alvan et al. (2013), though there are a few differences. In particular, our present study differs from the latter due to the treatment of inertial waves in convective regions (instead of gravity waves in stably stratified radiative regions in their case) with a latitudinal shear (instead of a vertical or radial shear). Our study, however, uses a local Cartesian model as in Lindzen & Barker (1985). Moreover, our boundary conditions are different, as detailed in Sect. 3.4.4, and the thickness of our shear region is fixed to one in scaled units, while Lindzen & Barker (1985) and Alvan et al. (2013) leave the thickness δ as a control parameter. We also checked for the existence of a critical level in the shear zone and the frequency range that delineates the regimes with and without the critical level.

We considered the perturbed reduced pressure Π and velocity v to be continuous at the interfaces y = 0 and y = 1. In the presence of the critical level y0 in zone II, we have a set of four analytical solutions whose values at y = 0 and y = 1 allow us to determine the reflection and transmission coefficients. The solutions to the wave propagation equation in zones I, II (below and after the critical level), and III are: (82)

where A and B are complex coefficients that we will express below. We recall that AT = 1 in our numerical model. In the shear region (zone II), the reduced pressure perturbation is given by Eqs. (40), which, at the pole, can be recast as (83)

with Ro = Λ. In the regions with no shear (zones I and III), Π takes a simpler expression, with Ro = 0 and σ = ω in zone I, and σ = ωkxU(1) = ωkxRo in zone III (noted as σ3 in the following). We note that while the reduced pressure is kept continuous to conserve the wave action flux across the interfaces, the first derivative of the latitudinal velocity v′ is not necessarily continuous at the interfaces. To find the transmission and reflection coefficients, we solved the system of equations that consists of matching conditions at interfaces as follows:

  • 1.

    v is continuous at y = 0 : AI + AR = vW(0);

  • 2.

    Π is continuous at y = 0 :

  • 3.

    v is continuous at y = 1 : eikT = vW(1),

  • 4.

    Π is continuous at y = 1 :

with the Whittaker functions vW. At the interfaces below and above the critical level y0 (i.e. y = 0 and y = 1), we have (84)

It should be noted that the first derivative of the Whittaker functions can be computed either numerically or analytically via the relationships in Abramowitz & Stegun (1972).

The equations in the above continuity relationships 1 to 4 are independent two by two (1 & 2 and 3 & 4), and A and B can be found first: (85)

The amplitude of the incident and reflected waves can be written, in terms of A and B, as: (86)

where (87)

The transmission and reflection coefficients are then: (88)

We emphasise that these factors depend notably on the location of the critical level and on the inertial frequency, which was not the case in Lindzen & Barker (1985) and Alvan et al. (2013).

We display in Figs. 6 and 7 the transmission and reflection coefficients as a function of the Rossby number Ro and the normalised inertial frequency ω in the regimes R < 1∕4 and R > 1∕4, respectively. We chose kx = kz = 0.1 in these plots, and R = 1∕4 gives , which delineates the two regimes, as we can see in Fig. 3. Areas that are hatched do not possess the co-rotation point σ = 0. In addition, the wavenumberkT was chosen with a positive sign in these regions to maintain an upward propagating wave (see Appendices D.1 and D.2). Areas that are not hatched feature a critical level, according to the Table D.1. In the case where R < 1∕4, over-reflection and over-transmission are both possible (see Fig. 6). One should notice that for Ro > 1, we always have by definition of R (89)

regardless of y, which makes the solutions of Eq. (61) tend towards pure exponential functions (i.e. without any imaginary part). Also, we do not see any over-reflection, nor over-transmission, in the hatched areas where there is no co-rotation point. This highlights the essential role of the critical level in inducing the over-reflection or over-transmission of inertial waves crossing the shear region in this regime.

The analysis of the regime where R > 1∕4 (in Fig. 7) is trickier. According to our discussion in Sect. 3.4.3, we expect a strong attenuation of the wave and of the wave action flux, as shown in Fig. 4. From this figure and for αk = 0.1, the damping is very strong for low positive Rossby numbers. This tendency is also found for both transmission and reflection coefficients. Nevertheless, one can also observe an unexpected regime of over-transmission near Ro = 0.8 and low frequency ω. Still, we must not forget that solutions in this regime, even near the critical level (see Eq. (71)), are not rigorously equivalent to wave-like functions. In particular, the amplification term that can be found at the first order in the Frobenius solutions becomes more prominent as the thickness of the shear zone is larger. This is especially true for the transmission coefficient. Assuming that Eq. (71) holds throughout zone II and corresponds to upward and downward waves, the transmission coefficient is modulated by |1 − y0|1∕2∕|0 − y0|1∕2, the amplitude ratio between the transmitted and incident waves. This term can be greater than one in the shear region. In particular, it is always greater than one when y0 < 0 (i.e. no critical level in the regime R > 1∕4; hatched areas in Fig. 7). On the contrary, this ratio is not present for the reflection coefficient since |R| is a function of the incident and reflected waves evaluated at y = 0. Though this hand-waving explanation does not formally demonstrate the origin of this amplification, it stresses the important role of the shear-region thickness and, more generally, of the geometry of the model.

In order to clarify whether the amplification is due to the geometry or the critical level, we need to investigate how the wave action flux changes before and after the critical level. The wave action flux is indeed the relevant quantity to investigate energy flux exchanges at a critical level.

thumbnail Fig. 6

Transmission coefficient (|T|, left panel) and reflection coefficient (|R|, right panel) when the box is at the pole. The coefficients are plotted in the regime R < 1∕4 (a possibly unstable case, see Sect. 3.4.4) as a function of the Rossby number Ro and the inertial frequency ω. The hatched areas do not feature critical points and correspond to regions where ω > kx Ro in our peculiar geometry (see Appendix D.2 for this particular matter). Vertical and longitudinal wave numbers are fixed: kx = 0.1 and kz = 0.1. Moreover, the contours that correspond to coefficients |R| and |T| that are equal to one are indicated by solid black lines. Crosses mark the set of parameters used in Fig. 9for the analysis of the behaviour of the velocity in the three-layer model.

Open with DEXTER
thumbnail Fig. 7

Same as Fig. 6, but for R > 1∕4. When Ro < 0 in these panels, hatched areas mark instances where ω > 0 (see Appendix D.2).

Open with DEXTER
thumbnail Fig. 8

Ratio of the wave action flux above and below the critical level against the Rossby number and the inertial frequency for R >1∕4 (left panel) and R < 1∕4 (right panel) when the box is at the pole. As in Fig. 7, kx = kz = 0.1 and hatchedzones represent areas without critical levels. Again, crosses indicate the set of parameters chosen to analyse the behaviour of the velocity in the three-layer model.

Open with DEXTER

4.2.2 Wave action fluxes below and above the shear region

Since v and Π are continuous at the interfaces between the shear and no-shear regions, the wave action flux is preserved and continuous in all three zones in the absence of friction and critical levels. However, it is discontinuous at the co-rotation point, as demonstrated in Sects. 3.4.3 and 3.4.4. Given the amplitude of the incident and reflected waves (Eq. (86)), we can calculate the ratio of the wave action flux below and after the co-rotation (see Appendix D.2 for the detailed calculation): (90)

The signs + or − can be chosen in regards to the wave action flux of the transmitted wave that can be positive or negative, depending on the presence of the critical level, while the energy flux is always positive in order to have an upward propagating wave in zone III (see Appendix D.2 for a more detailed discussion). This wave action flux ratio is displayed in Fig. 8 in the two regimes, R ≶ 1∕4. As expected, this ratio is equal to one when no critical level is present (hatched areas). Unlike in the previous section, the regime where R > 1∕4 no longer has amplification areas, and everywhere. This supports the idea that the critical level has nothing to do with the amplification phenomenon observed in the left panel of Fig. 7. As already observed in Fig. 4, the damping due to the critical level is strong except for Rossby numbers close to the threshold between the two regimes. Moreover, means that since the minus sign is taken in Eq. (90). Therefore, no over-reflection due to the critical level is expected in this regime. The other regime (R < 1∕4, right panel) features areas where the wave is over-reflected for (i.e. when ) and areas where the wave is over-transmitted for . For the first inequality (), the threshold between under- and over-reflection (around Ro ≈ 0.9) is the same as for the reflection factor (in the right panel of Fig. 6). For the second one (), the comparison with the transmission factor (in the left panel of Fig. 6) is more questionable. Still, these two points suggest that the critical level can induce the over-reflection and over-transmission phenomena in the regime where R < 1∕4.

4.2.3 Numerical solutions

In the previous sections, we examined how the shear parameter and the inertial wave frequency impact the reflection and transmission coefficients as well as the wave action flux. We now study particular cases of wave propagation through the critical level for fixed sets of parameters in both regimes, R ≶ 1∕4. To do that, we numerically solved Eqs. (80) for the three-layer model described in Sect. 4.1, for θ0 = 0 and a linear shear flow (n = 1 in zone II). We selected three pairs of values for the inertial frequency and the shear, two in the regime R > 1∕4 and one in the regime R < 1∕4. These values are marked by crosses in Figs. 68. In each case, the latitudinal velocity and the wave action flux were successively calculated through the three zones and plotted in Fig. 9. The numerical solution, which was computed by imposing the boundary condition AT = 1 and the continuous interfacial conditions for v and Π at y = {0, 1}, is the sum of the incident and reflected waves in zone I and is equal to a transmitted wave in zone III. The expressions for the incident, reflected, and transmitted waves are given by Eqs. (78) and (79) (see also Appendix D.1). In the shear region (zone II) of the upper panel of Fig. 9, the Whittaker solution has been added and it matches perfectly with the numerical solution below and above the critical level in each case. Moreover, Frobenius approximations for the latitudinal velocity (Eq. (71) when R > 1∕4 and Eq. (75) when R < 1∕4) and for the wave action flux (Eq. (73) when R > 1∕4 and Eq. (76) otherwise) have also been included. The coefficients a0 and b0 were determined by matching the numerical solution and its derivative to the Frobenius approximation for the velocity close to the critical level. For both the latitudinal velocity and the wave action flux, this first-order approximation gives satisfactory agreement with numerical solutions, although a slight deviation (regarding velocity) from the numerical solution can be observed as one moves away from the critical level. In particular, it should be mentioned that for the far left panels, corresponding to an upward wave is sufficient to correctly fit the numerical solution, meaning that the counter-propagating wave in zone I is reflected at y = 0. However, for the middle and right panels, it is not clear whether the first-order Frobenius solutions can be reconnected to the incident and reflected waves at y = 0.

We now examine attenuation or amplification phenomena in each column of panels in Fig. 9. In the left panels, for which R > 1∕4, the latitudinal velocity is strongly attenuated at the critical levels and so is the wave flux action, as we can expect from the left panel of Fig. 8. While the transmitted wave is totally absorbed, the reflected wave remains, which is consistent with analytical values of the transmission and reflection coefficients in Fig. 7 (see white crosses). In the middle panels, where we also have R > 1∕4, the wave is over-transmitted but not over-reflected, which is also consistent with Fig. 7 (see black crosses). In view of the wave action flux, the amplification of the transmitted wave does not seem to be related to the critical level because this quantity is greatly reduced after the critical level (see also the white cross in the left panel of Fig. 8). The third column of panels now depicts the regime where R < 1∕4. The wave is over-reflected by a factor of ~1.5 and over-transmitted by a factor of ~2, in concordance with the reflection and transmission coefficients plotted in Fig. 6. The wave action flux is negative, and by a factor of three, as observed in Fig. 8. These three case studies reinforce the idea from Booker & Bretherton (1967) in the case of stratified z-sheared flows: that the wave energy can be lost to the mean flow or, on the contrary, that the wave can take energy from the mean flow.

thumbnail Fig. 9

Numerical outputs of the three-layer model. Top: real part of the latitudinal velocity v against y. The quantities vnum, vI, vR, vT, vF, and vW are the numerical, incident, reflected, transmitted, first-order Frobenius, and Whittaker velocities, respectively. Bottom: wave action flux against y. The quantities , , , and are the numerical, transmitted, incident and reflected, and Frobenius wave action fluxes, respectively. For all panels, θ0 = 0, the mean flow is linear in the grey-shaded shear regions, and the critical level is marked by dashed lines. The horizontal wave numbers are set to kx = kz = 0.1. From left to right: (i) ω = 0.02 and Ro = 0.3 (R > 1∕4), (ii) ω = 0.002 and Ro = 0.8 (R > 1∕4), and (iii) ω = 0.09 and Ro = 1.8 (R < 1∕4).

Open with DEXTER

4.3 Numerical exploration at constant shear when the box is inclined

We now investigate wave propagation through the different critical levels when the box is inclined with respect to the rotation axis. We still assume that the shear region (zone II) has a linear shear flow profile (n = 1). In contrast to the polar configuration, we do not have analytical solutions to the ordinary differential equation. Instead of going for an extensive numerical investigation of the parameter space, we instead focus on the dynamics of inertial waves in our three-layer model as they cross the critical levels and σ = 0. Our results are presented in Fig. 10 for a box inclination of 10° relative to the rotation axis, a shearfixed to Ro = 0.3, and wavenumbers set to kx = kz = 1. The value of the frequency ω determines the existence and the nature of the critical level as detailed in the Table D.1. As in Fig. 9, we plot in each column of Fig. 10 the latitudinal velocity and the wave action flux. From left to right, we illustrate our results for the critical levels , σ = 0, and (i.e. there are two critical levels in the rightmost panels). One can notice that the first-order Frobenius approximation is no longer in good agreement with the numerical solution in the entire shear region, though it remains a reasonable approximation in the vicinity of a critical level. Unlike the polar case, the discrepancy far outside the critical levels is due to the linear approximation that the governing ODE takes around critical levels (see Eqs. (44) and (56)).

In the left panels of Fig. 10, the reflected and transmitted waves are strongly attenuated at the critical level y = y+. Part of the wave energy is laid down to the mean flow, as corroborated by the drop in the wave action flux. In the middle panels, we do not see any discontinuity at the co-rotation y = y0, which is in line with the theoretical analysis in Sect. 3.3.2 for a constant shear. However, the wave is over-reflected and over-transmitted, possibly due to the polynomial form of the solutions in the Frobenius series around y = y0. In the right panels, the wave encounters successively critical levels at y = y0 and at y = y. Although the wave going through the shear region is not attenuated at the co-rotation y = y0, it is completely absorbed at the second critical level y = y where the wave action flux drops to zero. This is consistent with the transmission coefficient in the left panel of Fig. 2. The latitudinal velocities displayed in the top-left and top-right panels of Fig. 10 support the concept of a valve effect. Indeed, according to our analysis with the Frobenius method in Sect. 3.3.1, and given the shear and wavenumbers of Fig. 10, the attenuation is strong for a downward wave meeting the critical level y+ (first panel), while the attenuation is strong for an upward wave that meets the critical level y (third panel). Before the critical level y and after the critical level y+, we observe fast oscillations of shorter and shorter periods close to the critical level, as already evidenced by Booker & Bretherton (1967). The analysis to determine how the wave is reflected in the shear zone can hardly be taken any further because Frobenius solutions are not fully separable into upward and downward waves.

We emphasise that the behaviour of the wave at co-rotation y = y0 when the box is at the pole stands out as being clearly different from the case when the box is inclined for a linear mean flow profile. This is particularly true in terms of the absolute value of the wave action flux that is subject to rise and drop in the polar configuration, whereas it remains conserved when the box is inclined. In this inclined case, the only way for a wave to be attenuated without friction is for it to meet critical levels . Depending on the critical level encountered, an upward or downward wave will not be attenuated in the same way, which is a phenomenon also known as the valve effect. In addition, we no longer observe amplification due to the critical level, but there is still ‘geometric’ amplification (for instance, in the middle panels where we can observe over-transmission and over-reflection), which can be explained by the exponential form of the Frobenius series.

thumbnail Fig. 10

Same quantities as in Fig. 9 (when the box is at the pole), but for a box tilted by 10 degrees relative to the pole and for different values of the shear, wavenumbers, and inertial frequencies. We note that, unlike Fig. 9, we do not have analytical solutions in the shear region. In all panels, the horizontal wavenumbers are set to kx = kz = 1 and the shear is Ro = 0.3. From left to right: inertial frequencies set to ω = 0.31, ω = 0.16, and ω = 0.1. Third panels: we indicate y0 and y with dasheddotted and vertical dashed lines, respectively.

Open with DEXTER

4.4 Numerical results with a non-constant shear

The choice of a linear mean flow profile allows the resolution of the ordinary differential equation at the pole and a simpler implementation of the Frobenius method. Nevertheless, the correspondence between a global and a local mean flow as presented in Sect. 2.2 involves higher-order terms than a simple linear dependence. Therefore, it is important to examine the effect of different mean flow profiles with non-zero U″. Now assuming n > 1 for the shear flow profile in zone II (see Eq. (77)), the Frobenius method in the inclined and polar cases still holds provided that (91)

Moreover, we specify that the critical level in y < 0 for even values of n is not examined given our numerical set-up where the shear region is located in the range y ∈ [0, 1]. The conditions to have critical levels inside the shear region (zone II) are the same as those for a constant shear (see Appendix D.2).

We show in Table 3 numerical values of the reflection and transmission coefficients along with the ratio of the wave action flux below and after the critical level in the six parameter sets illustrated in Figs. 9 and 10 for linear, square, and cubic mean flow profiles. For all cases, R, T, and change quite significantly between n = 1, 2, and 3. At the pole, the (‘geometric’) over-transmission found for Ro = 0.8 disappears for n = 2 and n = 3, where the reflected and transmitted waves are strongly attenuated. Similarly, in the case where Ro = 1.8, the over-transmission and over-reflection disappear when n = 3, whereas , which entails that the transmitted wave has taken wave action flux from the mean flow just after the critical level. We find that the Frobenius method does not give consistent results for the case where Ro = 1.8 and n = {2, 3}. When the box is tilted and for the critical level y+, the reflected wave is more attenuated than the transmitted one, which is consistent with the stronger attenuation of the counter-propagating wave discussed in the previous section. At the co-rotation and when the box is inclined (fifth row of Table 3), the wave action flux remains the same in the whole domain for linear or non-linear mean flow profiles, in agreement with the wave action flux derived analytically in Eq. (59). We speculate that the over-transmission and over-reflection are due to the polynomial form of the solutions in the shear region. These wave amplifications may be related to shear instabilities in this particular three-layer configuration, and they are probably not linked with the presence of a critical level whose implications on the flow are well diagnosed by the wave action flux. In the last case (sixth row), where two critical levels co-exist, the transmitted wave is less attenuated at y = y from n = 1 to n = 3, but it remains more attenuated than the reflected wave, as discussed in the previous section. Again, no jump is found in the wave action flux at the co-rotation despite analytical predictions.

The same analysis was carried out for different amplitudes of the Rayleigh friction force, up to σf = 10−2. Of course, the wave action flux is no longer constant, but we observe that all three parameters of Table 3 change very little compared to the case where σf = 10−8. This result is consistent with Alvan et al. (2013) in the context of gravity waves and vertical shear. We comment that while a low friction is mandatory in the numerical code to solve the ODE at singularities y±, it is not the case at the co-rotation y = y0.

Table 3

Reflection R and transmission T coefficients and the ratio of the wave action flux above and below a critical level for a linear, square, and cubic mean flow (n = 1,  2, and 3, respectively).

5 Astrophysical discussion

5.1 Latitudinal differential rotation in stars

In stars, latitudinal differential rotation is often characterised by the difference in rotation frequency between the equator and the pole, that is the quantity ΔΩ = Ωeq − Ω0, where Ωeq is the rotation frequency at the equator (e.g. Barnes et al. 2017). We will now refer to ΔΩ as the shear contrast. Different regimes are distinguished according to the value of ΔΩ: anti-solar-like rotation for ΔΩ < 0, cylindrical rotation for |ΔΩ|≪ 1, and solar-like rotation rotation for (not-too-low) positive ΔΩ (e.g. ΔΩ∕Ω0 ≃ 0.3 for the Sun). Several works based on three-dimensional numerical simulations have explored the range of physical parameters leading to each aforementioned regime in stars and in giant planets (e.g. Gastine et al. 2013; Varela et al. 2016; Beaudoin et al. 2018) In particular, Brun et al. (2017) derived a criterion based on mixing length theory and calibrated with three-dimensional simulations that determines the rotation profile of a solar-like star. This criterion is based on the fluid Rossby number Rof, defined as (92)

where Rof,⊙ = 0.89 is the solar fluid Rossby number and Ω* and M* are the mean rotation and the mass of the star respectively, normalised with their solar values. Brun et al. (2017) highlighted the following three regimes: (i) Rof > 1 for anti-solar-like rotation, (ii) 0.3 < Rof < 0.9 for solar-like rotation, and (iii) Rof ≲ 0.3 for cylindrical rotation. Furthermore, they introduced the shear contrast ΔΩS at the co-latitude 30° as (93)

since the rotation frequency is often ill-defined at low co-latitudes in three-dimensional numerical simulations in spherical geometry. From their three-dimensional simulations, Brun et al. (2017) obtained the following scaling: (94)

where ΔΩS,⊙≃ 565 × 10−9 s−1 is the solar value of ΔΩS calculated from García et al. (2007). Using the mean flow profile (Eq. (4)), ΔΩS can be related to our shear contrast ΔΩ via (95)

We show in Fig. 11 the quantity ΔΩ∕Ω0, expressed inEq. (95), versus the age of solar-like stars for K to G spectral types. To compute this quantity, we used grids of the one-dimensional stellar evolution code STAREVOL (see Amard et al. 2019, for details of the code). In light of Fig. 11, the stars in the pre-main sequence (age ≲ 100 Myr) exhibit cylindrical rotation as they are fast rotating. During the main sequence, stars mostly feature solar-like rotation, while anti-solar-like rotation is observed at the end of the main sequence from 0.8 to 1.1 M. According to Fig. 11, a limit on the absolute value of the normalised shear contrast can be set to |ΔΩ|∕Ω0 < 0.5. However, as already stressed by Benomar et al. (2018), the latitudinal shear inferred by asteroseismology can be much larger than predicted by numerical simulations. This can actually be inferred by comparing the shear factors in their work (see Table S3 from the supplementary materials in Benomar et al. 2018) with ours given in Fig. 11. Moreover, according to their study, cylindrical and anti-solar differential rotation are hardly unambiguously detectable. Finally, one should recall that the Brun et al. (2017) scaling laws given in Eqs. (92) and (94) are derived for K and G spectral type stars only.

Since we now know values of the shear contrast, we can calculate the ‘shear’ Rossby number, Ro = U′∕(2Ω), given by thefollowing relationship: (96)

which has been derived from Eq. (5) by keeping only zero-order terms in y. Taking χ ≃ 0.3 as a representative value of the shear contrast for main sequence G and K stars, we find that the Rossby number is maximal when θ0 ≃ 55°, its maximum value being Ro ≃−0.17. In particular, Ro ≃−0.013 for θ0 = 10° and Ro ≃−0.076 for θ0 = 80°. These values are useful for interpreting wave flux action transmission at critical levels , considering Fig. 2 (we recall that at the co-rotation the wave action flux is fully transmitted). A downward (upward) propagating wave through () is: (1) totally absorbed provided that kz ≳ 0.1kx (αk ≳ 0.1) at θ0 = 10°, and that kzkx (αk ≳ 1) at θ0 = 80°, (2) strongly attenuated for kz ~ 0.1kx (αk ~ 0.1) at θ0 = 80°, and (3) fully transmitted given that kz ≪ 10−1kx (αk ≪ 10−1) for both inclinations.

These results also hold for anti-solar-like differential rotation since the transmission factor is a function of |Ro|. For larger values of |Ro| (i.e. for larger values of the shear contrast), waves are less damped at critical levels at a given αk = kzkx. The connection between this ratio of the vertical and azimuthal wavenumbers in the local model and an equivalent ratio of global wavenumbers in the spherical geometry is not straightforward. A first hint can be to state that kz ~ kr, where kr is the wavenumber in the global radialdirection, while kx ~ m∕(r0 sinθ0), where m is the azimuthal order of the considered mode of the tidal potential (when m≠ 0; Zahn et al. 1997). Then, we get αkkrr0 sinθ0m ≡ 2π r0λr × sinθ0m by introducing λr, the radial wavelength of the tidal wave. In the case where r0 > λrm, we should thus be in the regime where the tidal wave is attenuated.

thumbnail Fig. 11

Absolute values of the shear contrast normalised by the rotation at the pole against the age of K and G spectral type stars. Solid lines feature 0.3 < Rof < 0.9, dashed lines Rof ≲ 0.3, and transparent lines Rof > 1.

Open with DEXTER

5.2 Cylindrical differential rotation in Jupiter and Saturn

Mathis et al. (in prep.) have developed an equatorial model to examine inertial wave properties in the outer convective layers of giant gaseous planets, such as Saturn and Jupiter, which are subject to cylindrical differential rotation. In their model (built in cylindrical coordinates), they derived a Schrödinger-like differential equation for under the anelastic approximation, where ρ is the density, r the axial distance coordinate, and vr the axial velocity. For free inertial waves, their second-order differential equation is: (97)

where lΘ and m denote the equatorial and azimuthal wavenumbers, respectively, is the (linear) Doppler-shifted frequency, and κr is the ‘axial’ epicyclic frequency defined as . For cylindrical differential equation, the co-rotation resonance results in critical cylinders (see Baruteau & Rieutord 2013) characterised by a critical axial distance r = rc. The Taylor expansion of Eq. (97) at the first order around rc gives: (98)

by setting the local Rossby number in cylindrical coordinates. By writing , Eq. (98) becomes very similar to our ODE Eq. (22) when the box is located at the south pole and for a constant shear (i.e. with , f = −1, and ). We note that, when the local shear box model is located at the equator, the latitudinal coordinate y is directed along the (vertical) rotation axis, whereas when the box is at the poles, y is the axial distance. That is why a polar configuration of the box best reproduces the ‘equatorial’ model of Mathis et al. (in prep.). Moreover, the convention of a plus sign in the Doppler-shifted frequency explains why Eq. (98) is analogous to our wave propagation equation when the box is at the south pole rather than at the north pole.

In Fig. 12, we show Jupiter’s and Saturn’s local Rossby number from Mathis et al. (in prep.). Cylindrical differential rotation extends in the outer layer of the convective envelope of both planets, in agreement with the Juno and Cassini Grand Finale observations(Kaspi et al. 2017; Galanti et al. 2019, respectively). According to Fig. 3, where we notably plotted at the south pole (dark red and purple areas), we can assess the role of the critical level for wave transmission across the co-rotation, in terms of the wavenumber ratio αk. In Fig. 12, the Rossby number satisfies Ro ≲ 0.27 for Saturn and Ro ≲ 0.07 for Jupiter. Given this range of values, two regimes can be evidenced for waves and wave action fluxes through the co-rotation: (1) waves are strongly attenuated for lΘm (see also Fig. 4 for Ro = −|Ro| as the transmission factor in this figure is plotted for the north pole) and (2) waves can be over-reflected and over-transmitted for lΘm and can potentially lead to instabilities given specific boundary conditions.

To give an idea of the values these wavenumbers can take, we have listed in Table 4 three typical orbital states, where, in order, the asynchronous, eccentricity, and obliquity tides are supposed to be dominant (Ogilvie 2014). These states are described by the ‘spherical’ quadrupolar components of the dominant terms in the tidal potential, namely the degree l and the order m of the spherical harmonics. The analogy with the equatorial model is then made to get m, and lΘ is chosen toapproximate, as best as possible, the behaviour of the Legendre polynomial around the equator with a simple trigonometric function {exp[i(lΘΘ + ϕ)]}, where ϕ is the appropriate phase (Mathis et al., in prep.). To find the associated wave attenuation at co-rotation for the three main tides, one can use Figs. 3 and 4 for αk ≥ 1 and Ro ≲ 0.27 and look at the south pole (as Fig. 4 is plotted at the north pole, one has to take the opposite Rossby number). From Fig. 3, we can assess that waves excited by these tides are always in the so-called stable regime for these ranges of Ro and αk, which excludes an amplification of these waves. Moreover, in Fig. 4, we also observe that waves are completely absorbed at co-rotation for our given ranges of parameters. By consequence, waves excited by the asynchronous, inclined, or eccentric tides in Jupiter and Saturn are expected to transfer all their wave action flux to the mean flow at co-rotation.

thumbnail Fig. 12

Rossby numbers for Jupiter and Saturn as a function of the axial distance r.

Open with DEXTER
Table 4

Regime of wave transmission at co-rotation deduced from Figs. 3 and 4 for three orbital states of a satellite around a giant gaseous planet.

6 Conclusion and perspectives

The present study was motivated by the works of Baruteau & Rieutord (2013) and Guenel et al. (2016a,b), who showed that differential rotation can strongly affect the propagation and dissipation properties of (tidal) linear inertial waves. They considered different rotation profiles typical of stellar and planetary interiors and pointed out that tidal waves can deeply interact with zonal flows at co-rotation resonances, leading to intense wave energy dissipation, along with possible instabilities. In this paper, we have investigated the transmission of free inertial waves with latitudinal stratification and differential rotation at the co-rotation resonance (characterised by a zero Doppler-shifted wave frequency) and, more broadly, at critical levels (any singularities of the governing second-order wave propagation equation in the inviscid limit). For this purpose, we built a new local Cartesian box model with horizontal shear, modelling a small patch of the convective zone of a low-mass star or a giant planet. By considering the inclination of the local reference frame relative to the rotation axis, we have examined the effect on wave propagation through a critical level of a conical rotation profile at a general co-latitude when the box is tilted, or of a cylindrical rotation profile when the box is at the north or south poles. These rotation profiles are inspired by those observed or expected in the Sun, low-mass stars, and the giant gaseous planets in our Solar System. Three critical levels can be identified when the box is inclined relative to the rotation axis: the co-rotation resonance and two other critical levels that arise from the inclination between the gravity and the rotation vectors, which are defined by a Doppler-shifted frequency equal to plus or minus the latitudinal component of the rotation frequency. When the box is at the poles, critical levels are restricted to the co-rotation.

In order to diagnose the behaviour of a wave passing through a critical level for both aforementioned rotation profiles, we made use of an invariant called the wave action flux, which is independent of the latitudinal coordinate in a non-dissipative fluid flow. This invariant was used when the ‘directional’ flux of angular momentum (here latitudinal) cannot be constructed easily from the mean perturbed velocity, as is the case, for example, in Lindzen & Tung (1978) for Rossby waves in plane-parallel shear flows. The wave action flux has already been used in vertically stratified shear flows in the presence of rotation or magnetic fields to interpret the role of critical levels (Grimshaw 1975a, 1979; Andrews & McIntyre 1978; Mathis 2009; Mathis & de Brye 2012). Using the condition that this invariant is discontinuous at critical levels, we demonstrated in Sect. 3 that waves can be either fully transmitted, damped, or even amplified after passing through critical levels as a result of wave action flux exchanges. These different regimes of wave transmission are found with both conical and cylindrical rotation profiles; they depend on the critical level encountered, on the wave properties (e.g. the propagation direction and wavenumbers), and on the profile of the mean flow. Table 2 summarises the main analytical results.

We then compared our analytical results with a three-layer numerical model that comprises a shear zone, where the critical level is located, and two surrounding shear-free zones that allow incident, reflected, and transmitted waves. A difference with the analytical model is the introduction of a small dissipative force under the form of a Rayleigh friction (also called frictional force by Ogilvie 2009) to avoid strict singularities. This does not seem to affect the results since analytical and numerical results match quite well when using a power-law mean flow profile and varying the friction. This conclusion is also shared by the work of Alvan et al. (2013), who studied co-rotation resonances for gravity waves propagating in stratified and vertically shear flows.

Based on the analytical results, we discussed possible applications to stellar and planetary interiors in Sect. 5. We have estimated the rate of differential rotation in solar-like stars using the shear contrast (the rotation difference between the pole and theequator) and in giant gaseous planets through the local Rossby number (the ratio between the shear and the rotation frequency incylindrical coordinates). We find that throughout the lifetimes of K- and G-type stars, as well as for Jupiter and Saturn at the present time, a regime where inertial waves are strongly damped is largely preferred in the convective envelope of these objects. Similar conclusions were found by Alvan et al. (2013) for internal gravity waves through critical levels in the core of solar-like stars.

It is interesting to discuss different regimes of the wave transmission in terms of angular momentum transfer for cases of strong damping and wave amplification. First, we have to underline that the theoretical analysis presented in Sect. 3 (using the Frobenius method) does not adequately characterise wave (over-)reflection as similarly observed in the numerical section (Sect. 4). What we can access is the wave action flux on either side of the critical level. The analysis of the changes in the wave action flux across the critical level allows us to understand whether energy is deposited into or extracted from the mean flow, in line with the work carried out by Miles (1961), Booker & Bretherton (1967), Grimshaw (1975a), and Lindzen & Barker (1985). We have demonstrated that, in the presence of a locally conical differential rotation, a valve effect can be found for critical levels other than the co-rotation, analogous to the results of Acheson (1972) and Grimshaw (1975a) for hydromagnetic and gravito-inertial waves with a vertical shear. For these peculiar critical levels, waves can be attenuated when going in one direction, mainly featured by the sign of the rotation components in the box, or fully transmitted when going in the other direction. For cylindrical differential rotation, we have found a criterion analogous to the Miles-Howard theorem for stratified shear flows (Miles & Howard 1964), which, for inertial waves such as those studied in this work, can be formulated as: (99)

The above criterion depends on the shear Rossby number Ro = U′∕(2Ω) and the vertical (kz) and longitudinal (kx) wavenumbers. This last point is an important difference with the Miles-Howard criterion, which does not involve wavenumbers. For this reason, the analogy between Eq. (99) and the Miles-Howard stability criterion must be taken with care. We also stress that Eq. (99) is very different from the Rayleigh’s inflection point theorem for Rossby waves, which are a sub-class of inertial waves, when neglecting the vertical or the radial perturbed velocity (Bretherton 1966; Lindzen & Tung 1978). The Miles-Howard criterion allows us to distinguish between critical levels where strong wave attenuation is expected for Ri > 1∕4 (where Ri is the Richardson number) and those where over-reflection and over-transmission can lead to potential shear instabilities for Ri < 1∕4. Lindzen (1988) warns, however, that over-reflection and over-transmission are a necessary but not sufficient condition for shear instability. Such amplifications leading to the instability require peculiar conditions in a three-layer model, where the shear zone that features the critical level is surrounded by a region of incoming propagating waves and by a ‘sink’ zone to force waves to cross the evanescent shear zone. Special boundary conditions are necessary for the wave to successively return to the critical level and induce wave amplitude growth. Recent studies (see e.g. Carpenter et al. 2012, for a review) have revisited instabilities in stratified shear flows by studying multiple counter-propagating waves that can interact with one another to grow in amplitude with time (with conditions such as phase-locking). A parallel was drawn between over-reflection mechanisms and interacting counter-propagating waves by Harnik & Heifetz (2007) to describe baroclinic instabilities for Rossby waves.

Contrary to what our results predict, Baruteau & Rieutord (2013) did not observe any instabilities of inertial waves when using cylindrical differential rotation. Several reasons can be put forward to explain this discrepancy, such as boundary conditions (as discussed in the previous paragraph) or the values of the shear and horizontal wavenumbers, since the inertial waves may not be in the regime which allows instabilities, according to the criterion in Eq. (99), which needs to be adapted further to the global cylindrical geometry used in the work of Baruteau & Rieutord (2013). We stress that when exploring different power-laws for the mean flow profiles, over-reflection was not retrieved for a non-linear mean flow in cylindrical differential rotation. Guenel et al. (2016a,b) did observe instabilities with conical differential rotation, but only for sufficiently low viscosities, whereas our study shows little dependence on the friction and rather highlights possible over-transmission for non-linear flows. Lastly, we underline that a temporal analysis on the growth rate of perturbations should be undertaken to unravel instabilities, which has not been performed in this paper but has been in other separate papers (Park et al. 2020, 2021).

This ab initio analytical study is thus a first step towards understanding how inertial waves interact with a mean flow subject to latitudinal differential rotation at critical levels in the context of tidal dissipation in differentially rotating stars and planets. Possible feedbacks of the perturbed wave on equilibrium quantities and the mean flow are not taken into account in this study, nor are non-linearities in the perturbed hydrodynamical wave equations. Nonetheless, they should be considered in future studies since Barker & Ogilvie (2010) and Baruteau & Rieutord (2013) suggested important non-linear effects for inertial waves at co-rotation. Finally, magnetism may also play an important role in dissipating or redistributing angular momentum at critical levels through magnetic stresses (e.g. Wei 2016, 2018; Lin & Ogilvie 2018; Astoul et al. 2019).

Acknowledgements

We would like to thank the anonymous referee, as well as A. Barker, for the helpful comments and suggestions regarding our work. A.A., J.P., and S.M. acknowledge funding by the European Research Council through the ERC grant SPIRE 647383. The authors acknowledge the PLATO CNES funding at CEA/IRFU/DAp and IRAP. This research hasmade use of NASA’s Astrophysics Data System and of the software MATLAB version R2018a. Finally, we thank Quentin André and Jéremy Ahuir for fruitful discussions on the details of the analytical model.

Appendix A Derivation of the wave propagation equation

In this section, we detail the derivation of the ODE (cf. Eqs. (22) and (23)) in a dissipative medium and for forced inertial waves. Writing the perturbation variables as wave-like functions in x and z directions, the hydrodynamic equations, (Eqs. (16), (17), and (18)), are, respectively:

with s = σ + f a complex frequency that includes the Rayleigh friction frequency σf. We have removed the symbol ′ in each perturbed quantity, and it now refers to the derivative with respect to y. The projection of the curl of the perturbed momentum Eq. (A.1) on the (ex,   ey, ez) basis is:

We note that the term is of order ϵ and thus neglected. This enables us to define and use the quantity Π = pρ0, which we call reduced pressure in the present study. Using the continuity equation, the linear combination ikz (A.4) − ikx(A.6) reads: (A.7)

where . Thus, gives a second-order ordinary differential equation on the latitudinal velocity v, with a source term S: (A.8)

Coefficients are written as: (A.9)

Appendix B Derivation of the polarisation relations

The Navier-Stokes equation with Rayleigh friction, without forcing, and projected onto the Cartesian basis reads:

where we use Eq. (A.3) to replace the perturbed density. In order to obtain an equation for the reduced pressure and the latitudinal velocity perturbations only, we applied the following linear combination: (B.4)

which yields, using the continuity Eq. (A.2): (B.5)

To get the perturbed vertical velocity, one can use the linear combination , which yields, after some algebra: (B.6)

Furthermore, gives the perturbed longitudinal velocity: (B.7)

Appendix C The WKBJ approximation

Rapidly oscillating solutions are often studied within the WKBJ approximation (Press 1981). At the pole, the solution of the ODE (60) in this approximation takes the form: (C.1)

where A and B are the complex amplitudes of the wave function Ψ and κ(y) is the complex potential associated with the ODE. We can determine the validity domain of this approximation in the same way Alvan et al. (2013) did. In the WKBJ approximation, Ψ satisfies (C.2)

where for a constant shear. The WKBJ approximation is valid provided that: (i) Ro < 1 so that R > 0, and (ii) is negligible in front of . Then, we introduce the Liouville transformation: (C.3)

We deduce (C.4)

Equation (C.2) thus becomes: (C.5)

The WKBJ approximation states that |Φ|≪ 1. Given the definition of f, this leads to the validity condition |R|≫ 1∕4.

Appendix D Analytical properties in the no-shear regions (zones I and III)

D.1 Wave-like solutions

The ordinary differential equation without shear can be written as: (D.1)

For the following, we introduce s1 = ω + f in zone I, and s3 = ωkxΛ + f in zone III. In the no-flow region I, the analytic solution of Eq. (D.1) can be written in the form: (D.2)

which is the sum of an incident wave and a reflected wave of amplitudes AI and AR and wavenumbers: (D.3)

respectively. These relationships are given by the dispersion relation satisfied by the incident and reflected waves: (D.4)

Similarly, we write the latitudinal velocity in the uniform mean flow region III as: (D.5)

where AT is the amplitude of the transmitted wave and (D.6)

is its wavenumber, also given by the dispersion relation of the transmitted wave. It is important to stress that all varying parameters in the model, such as Λ, kx, kz, and ω, have been chosen such that kI, kR, and kT are not complexwithout friction. This prevents widely diverging waves in zones I or III.

Table D.1

Requirements on ω for the existence of critical points inside the shear region II depending on the sign of the shear Λ.

D.2 Wave action flux

As we did in Sect. 3.2.5, one can also determine the wave action flux (Eq. (41)) in region I, (D.7)

and in region III, (D.8)

with σ3 = ωkxΛ, and where we use the reduced pressure for all three waves: (D.9)

It is noteworthy that in the no-shear regions the latitudinal flux of energy is preserved (see also Eliassen & Palm 1961), and so is the wave action flux. Indeed, the characteristic frequencies (ω and σ3) in both zones (I andIII, respectively) are constants. Furthermore, as we stated in Sect. 3.2.5, the direction of the group velocity is given by . By consequence, the direction of energy propagation is constrained by the sign of the Doppler-shifted frequency σ3 in zone III, and by the sign of inertial frequency ω in zone I. Considering ω > 0, the incident wave is properly named since the group velocity is positive; likewise, the reflective wave has negative group velocity and is thus moving downwards. It is a little more complicated for the so-called transmitted wave. The sign of σ3 is directly related to the presence or the absence of a critical level inside zone II.

Table D.1 summarises which frequency has to be exited so that waves can meet one or several critical points in the shear region of range 0 < y < 1 for a linear mean flow profile.

Thus, if co-rotation is met in zone II with Λ > 0, we automatically have a negative Doppler-shifted frequency in zone III (i.e. σ3 < 0). Therefore, one must choose the − sign in the expression of kT (Eq. (D.6)) in order to construct a wave that moves away from the critical level. The same reasoning applies for the critical level . It is trickier for the singularity because , so σ3 can be either positive or negative in the interval [0, π∕2] depending on the values of and kxΛ.

References

  1. Abramowitz, M., & Stegun, I. A. 1972, Handbook of Mathematical Functions (New York: Dover Publications Inc.) [Google Scholar]
  2. Acheson, D. J. 1972, J. Fluid Mech., 53, 401 [Google Scholar]
  3. Alvan, L., Mathis, S., & Decressin, T. 2013, A&A, 553, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Amard, L., Palacios, A., Charbonnel, C., et al. 2019, A&A, 631, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. André, Q., Barker, A. J., & Mathis, S. 2017, A&A, 605, A117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Andrews, D. G., & McIntyre, M. E. 1978, J. Fluid Mech., 89, 647 [Google Scholar]
  7. Astoul, A., Mathis, S., Baruteau, C., et al. 2019, A&A, 631, A111 [EDP Sciences] [Google Scholar]
  8. Auclair-Desrotour, P., Le Poncin-Lafitte, C., & Mathis, S. 2014, SF2A-2014: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics, eds. J. Ballet, F. Martins, F. Bournaud, R. Monier, & C. Reylé, 199 [Google Scholar]
  9. Auclair Desrotour, P., Mathis, S., & Le Poncin-Lafitte, C. 2015, A&A, 581, A118 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Barker, A. J., & Ogilvie, G. I. 2009, MNRAS, 395, 2268 [Google Scholar]
  11. Barker, A. J., & Ogilvie, G. I. 2010, MNRAS, 404, 1849 [NASA ADS] [Google Scholar]
  12. Barnes, J. R., Collier Cameron, A., Donati, J. F., et al. 2005, MNRAS, 357, L1 [Google Scholar]
  13. Barnes, J. R., Jeffers, S. V., Haswell, C. A., et al. 2017, MNRAS, 471, 811 [Google Scholar]
  14. Baruteau, C., & Masset, F. 2008, ApJ, 672, 1054 [Google Scholar]
  15. Baruteau, C., & Rieutord, M. 2013, J. Fluid Mech., 719, 47 [Google Scholar]
  16. Bazot, M., Benomar, O., Christensen-Dalsgaard, J., et al. 2019, A&A, 623, A125 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Beaudoin, P., Strugarek, A., & Charbonneau, P. 2018, ApJ, 859, 61 [Google Scholar]
  18. Benbakoura, M., Réville, V., Brun, A. S., Le Poncin-Lafitte, C., & Mathis, S. 2019, A&A, 621, A124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Benomar, O., Bazot, M., Nielsen, M. B., et al. 2018, Science, 361, 1231 [Google Scholar]
  20. Bolmont, E., & Mathis, S. 2016, Celest. Mech. Dyn. Astron., 126, 275 [Google Scholar]
  21. Bolmont, E., Gallet, F., Mathis, S., et al. 2017, A&A, 604, A113 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Booker, J. R., & Bretherton, F. P. 1967, J. Fluid Mech., 27, 513 [Google Scholar]
  23. Bretherton, F. P. 1966, Quart. J. R. Meteorol. Soc., 92, 325 [Google Scholar]
  24. Bretherton, F. P., & Garrett, C. J. R. 1968, Proc. R. Soc. London Ser. A, 302, 529 [Google Scholar]
  25. Broad, A. S. 1995, Quart. J. R. Meteorol. Soc., 121, 1891 [Google Scholar]
  26. Brun, A. S., García, R. A., Houdek, G., Nandy, D., & Pinsonneault, M. 2015, Space Sci. Rev., 196, 303 [Google Scholar]
  27. Brun, A. S., Strugarek, A., Varela, J., et al. 2017, ApJ, 836, 192 [Google Scholar]
  28. Bryan, G. H. 1889, Phil. Trans. R. Soc. London Ser. A, 180, 187 [Google Scholar]
  29. Carpenter, J. R., Tedford, E. W., Heifetz, E., & Lawrence, G. A. 2012, Appl. Mech. Rev., 64, 061001 [Google Scholar]
  30. Cartan, E. 1922, Bull. Sci. Math., 46, 317 [Google Scholar]
  31. Damiani, C., & Mathis, S. 2018, A&A, 618, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Debras, F., & Chabrier, G. 2019, ApJ, 872, 100 [Google Scholar]
  33. Decressin, T., Mathis, S., Palacios, A., et al. 2009, A&A, 495, 271 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Duguid, C. D., Barker, A. J., & Jones, C. A. 2020, MNRAS, 491, 923 [Google Scholar]
  35. Eliassen, A., & Palm, E. 1961, On the Transfer of Energy in Stationary Mountain Waves (Det Norske Videnskaps-Akademi i Oslo. Geofysiske publikasjoner) (I kommisjon hos Aschehoug) [Google Scholar]
  36. Favier, B., Barker, A. J., Baruteau, C., & Ogilvie, G. I. 2014, MNRAS, 439, 845 [Google Scholar]
  37. Ford, E. B., & Rasio, F. A. 2006, ApJ, 638, L45 [Google Scholar]
  38. Fuller, J., Luan, J., & Quataert, E. 2016, MNRAS, 458, 3867 [Google Scholar]
  39. Galanti, E., Kaspi, Y., Miguel, Y., et al. 2019, Geophys. Res. Lett., 46, 616 [Google Scholar]
  40. Gallet, F., Bolmont, E., Mathis, S., Charbonnel, C., & Amard, L. 2017, A&A, 604, A112 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Gallet, F., Bolmont, E., Bouvier, J., Mathis, S., & Charbonnel, C. 2018, A&A, 619, A80 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. García, R. A., Turck-Chièze, S., Jiménez-Reyes, S. J., et al. 2007, Science, 316, 1591 [Google Scholar]
  43. Gastine, T., Wicht, J., & Aurnou, J. M. 2013, Icarus, 225, 156 [Google Scholar]
  44. Gerkema, T., Zimmerman, J. T. F., Maas, L. R. M., & van Haren, H. 2008, Rev. Geophys., 46, RG2004 [Google Scholar]
  45. Gliatto, M. T., & Held, I. M. 2020, J. Atm. Sci., 77, 859 [Google Scholar]
  46. Goldreich, P., & Nicholson, P. D. 1989, ApJ, 342, 1075 [Google Scholar]
  47. Goldreich, P., & Tremaine, S. 1979, ApJ, 233, 857 [Google Scholar]
  48. Greenspan, H. P. 1969, The theory of rotating fluids (Cambridge University Press) [Google Scholar]
  49. Grimshaw, R. H. J. 1975a, J. Fluid Mech., 70, 287 [Google Scholar]
  50. Grimshaw, R. H. J. 1975b, Tellus Ser. A, 27, 351 [Google Scholar]
  51. Grimshaw, R. 1979, Geophys. Astrophys. Fluid Dyn., 14, 303 [Google Scholar]
  52. Guenel, M., Baruteau, C., Mathis, S., & Rieutord, M. 2016a, A&A, 589, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Guenel, M., Mathis, S., Baruteau, C., & Rieutord, M. 2016b, ArXiv e-prints [arXiv:1612.05071] [Google Scholar]
  54. Harnik, N., & Heifetz, E. 2007, J. Atm. Sci., 64, 2238 [Google Scholar]
  55. Hut, P. 1980, A&A, 92, 167 [NASA ADS] [Google Scholar]
  56. Jackson, B., Greenberg, R., & Barnes, R. 2008, ApJ, 678, 1396 [Google Scholar]
  57. Jones, W. L. 1967, J. Fluid Mech., 30, 439 [Google Scholar]
  58. Jones, W. L. 1968, J. Fluid Mech., 34, 609 [Google Scholar]
  59. Jouve, L., & Ogilvie, G. I. 2014, J. Fluid Mech., 745, 223 [Google Scholar]
  60. Kaspi, Y., Guillot, T., Galanti, E., et al. 2017, Geophys. Res. Lett., 44, 5960 [Google Scholar]
  61. Kippenhahn, R., Weigert, A., & Weiss, A. 2012, Stellar Structure and Evolution (Berlin: Springer) [Google Scholar]
  62. Lai, D. 2012, MNRAS, 423, 486 [Google Scholar]
  63. Lainey, V., Jacobson, R. A., Tajeddine, R., et al. 2017, Icarus, 281, 286 [Google Scholar]
  64. Lainey, V., Casajus, L. G., Fuller, J., et al. 2020, Nat. Astron., 4, 1053 [Google Scholar]
  65. Latter, H. N., & Balbus, S. A. 2009, MNRAS, 399, 1058 [Google Scholar]
  66. Lin, Y., & Ogilvie, G. I. 2018, MNRAS, 474, 1644 [Google Scholar]
  67. Lindzen, R. S. 1988, Pure Appl. Geophys., 126, 103 [Google Scholar]
  68. Lindzen, R. S., & Barker, J. W. 1985, J. Fluid Mech., 151, 189 [Google Scholar]
  69. Lindzen, R. S., & Tung, K. K. 1978, J. Atm. Sci., 35, 1626 [Google Scholar]
  70. Luan, J., Fuller, J., & Quataert, E. 2018, MNRAS, 473, 5002 [Google Scholar]
  71. Maas, L. R. M., & Lam, F. P. A. 1995, J. Fluid Mech., 300, 1 [Google Scholar]
  72. Maslowe, S. A. 1986, Ann. Rev. Fluid Mech., 18, 405 [Google Scholar]
  73. Mathis, S. 2009, A&A, 506, 811 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Mathis, S. 2015, A&A, 580, L3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Mathis, S. 2019, EAS Publ. Ser., 82, 5 [Google Scholar]
  76. Mathis, S., & de Brye, N. 2012, A&A, 540, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  77. Mathis, S., Palacios, A., & Zahn, J. P. 2004, A&A, 425, 243 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  78. Mathis, S., Prat, V., Amard, L., et al. 2018, A&A, 620, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  79. Miles, J. W. 1961, J. Fluid Mech., 10, 496 [Google Scholar]
  80. Miles, J. W., & Howard, L. N. 1964, J. Fluid Mech., 20, 331 [Google Scholar]
  81. Militzer, B., Wahl, S., & Hubbard, W. B. 2019, ApJ, 879, 78 [Google Scholar]
  82. Morse, P. M., & Feshbach, H. 1953, Methods of Theoretical Physics (Feshbach Publishing) [Google Scholar]
  83. Ogilvie, G. I. 2009, MNRAS, 396, 794 [Google Scholar]
  84. Ogilvie, G. I. 2014, ARA&A, 52, 171 [Google Scholar]
  85. Ogilvie, G. I., & Lesur, G. 2012, MNRAS, 422, 1975 [Google Scholar]
  86. Ogilvie, G. I., & Lin, D. N. C. 2004, ApJ, 610, 477 [Google Scholar]
  87. Ogilvie, G. I., & Lin, D. N. C. 2007, ApJ, 661, 1180 [Google Scholar]
  88. Park, J., Prat, V., & Mathis, S. 2020, A&A, 635, A133 [EDP Sciences] [Google Scholar]
  89. Park, J., Prat, V., Mathis, S., & Bugnet, L. 2021, A&A, 646, A64 [EDP Sciences] [Google Scholar]
  90. Press, W. H. 1981, ApJ, 245, 286 [Google Scholar]
  91. Rieutord, M. 2015, Fluid Dynamics: An Introduction (Berlin: Springer) [Google Scholar]
  92. Rieutord, M., & Valdettaro, L. 1997, J. Fluid Mech., 341, 77 [Google Scholar]
  93. Rieutord, M., & Valdettaro, L. 2010, J. Fluid Mech., 643, 363 [Google Scholar]
  94. Rieutord, M., Georgeot, B., & Valdettaro, L. 2001, J. Fluid Mech., 435, 103 [Google Scholar]
  95. Ringot, O. 1998, PhD thesis (Paris 7), France [Google Scholar]
  96. Schmid, P., Henningson, D., & Jankowski, D. 2002, Appl. Mech. Rev., 55, B57 [Google Scholar]
  97. Schou, J., Antia, H. M., Basu, S., et al. 1998, ApJ, 505, 390 [Google Scholar]
  98. Shampine, L. F. & Reichelt, M. W. 1997, SIAM J. Sci. Comp., 18, 1 [Google Scholar]
  99. Thompson, M. J., Christensen-Dalsgaard, J., Miesch, M. S., & Toomre, J. 2003, ARA&A, 41, 599 [Google Scholar]
  100. Tsang, D., & Lai, D. 2009, MNRAS, 400, 470 [Google Scholar]
  101. Varela, J., Strugarek, A., & Brun, A. S. 2016, Adv. Space Res., 58, 1507 [Google Scholar]
  102. Watts, A. L., Andersson, N., & Williams, R. L. 2004, MNRAS, 350, 927 [Google Scholar]
  103. Wei, X. 2016, ApJ, 828, 30 [Google Scholar]
  104. Wei, X. 2018, ApJ, 854, 34 [Google Scholar]
  105. Yamanaka, M., & Tanaka, H. 1984, J. Meteorol. Soc. Jpn., 62, 1 [Google Scholar]
  106. Zahn, J. P. 1966, Ann. Astrophys., 29, 313 [NASA ADS] [Google Scholar]
  107. Zahn, J. P. 1977, A&A, 500, 121 [NASA ADS] [Google Scholar]
  108. Zahn, J. P., Talon, S., & Matias, J. 1997, A&A, 322, 320 [NASA ADS] [Google Scholar]

1

We note that a similar analysis can be undertaken at the south pole by using κ = 1 + U′.

2

The quantity Π is actually the enthalpy perturbation, but we will use the denomination ‘reduced pressure’ in the following.

3

We note that the velocity of energy density, Vg, in the latitudinal direction has been named the ‘group velocity’ for obvious physical reasons, but it differs from the group velocity defined in Sect. 3.2, which, unlike Vg, depends on latitudinal and vertical wavenumbers.

4

A singular point y0 of the second-order ODE v″(y) + p(y)v′(y) + q(y)v(y) = 0 is said to be regular when the function p(y)(yy0) and are analytical at y = y0.

All Tables

Table 1

Correspondence between local and global coordinate systems.

Table 2

Summary of analytical results at each critical level when the box is inclined or not with respect to the rotation axis, in the [0, π∕2] quadrant, and for positive wavenumbers.

Table 3

Reflection R and transmission T coefficients and the ratio of the wave action flux above and below a critical level for a linear, square, and cubic mean flow (n = 1,  2, and 3, respectively).

Table 4

Regime of wave transmission at co-rotation deduced from Figs. 3 and 4 for three orbital states of a satellite around a giant gaseous planet.

Table D.1

Requirements on ω for the existence of critical points inside the shear region II depending on the sign of the shear Λ.

All Figures

thumbnail Fig. 1

Sketch of the local Cartesian box in the convective region of a low-mass star or giant planet. Global spherical coordinates, including the depth r0, the inclination of the box θ0, and the co-latitude θ of a point of interest M inside the box, are shown to facilitate the analogy between the spherical and the Cartesian geometries.

Open with DEXTER
In the text
thumbnail Fig. 2

Transmission rate of a wave passing through any of the critical levels defined by as a functionof the absolute value of the shear Rossby number |Ro| = |U′| (where U′ is scaled by 2Ω0) and the ratioof the horizontal wavenumbers αk = kzkx for a co-latitude of the box θ0 = 10° (left panel) and θ0 = 80° (right panel).

Open with DEXTER
In the text
thumbnail Fig. 3

Diagram showing how compares to 1∕4 for two positions of the box (θ0 = 0 and π), as well as for a range of Rossby numbers Ro and a range of the vertical to longitudinal wavenumber ratio αk. The solid and dashed black lines mark where at the north and south poles, respectively. The purple domain shows where , and the white region where , regardless of whether the box is at the north or the south pole. In the dark red region, and and vice versa in the orange region.

Open with DEXTER
In the text
thumbnail Fig. 4

Transmission rate TF of the wave action flux across the co-rotation when the box is located at the pole for R > 1∕4. It is displayed against the Rossby number Ro and the ratioof the horizontal wavenumbers αk = kzkx. The forbidden region where R < 1∕4 is shown in white, and the dark red cone corresponds to values of TF that are lower than 10−9.

Open with DEXTER
In the text
thumbnail Fig. 5

Mean flow profiles used in the three-zone numerical model against y. In the no-shear regions I (y < 0) and III (y > 1), the mean flow is uniform and set respectively to U = 0 and U = Λ. The shearedregion II (grey-shaded) can have a linear, square, or cubic mean flow profile.

Open with DEXTER
In the text
thumbnail Fig. 6

Transmission coefficient (|T|, left panel) and reflection coefficient (|R|, right panel) when the box is at the pole. The coefficients are plotted in the regime R < 1∕4 (a possibly unstable case, see Sect. 3.4.4) as a function of the Rossby number Ro and the inertial frequency ω. The hatched areas do not feature critical points and correspond to regions where ω > kx Ro in our peculiar geometry (see Appendix D.2 for this particular matter). Vertical and longitudinal wave numbers are fixed: kx = 0.1 and kz = 0.1. Moreover, the contours that correspond to coefficients |R| and |T| that are equal to one are indicated by solid black lines. Crosses mark the set of parameters used in Fig. 9for the analysis of the behaviour of the velocity in the three-layer model.

Open with DEXTER
In the text
thumbnail Fig. 7

Same as Fig. 6, but for R > 1∕4. When Ro < 0 in these panels, hatched areas mark instances where ω > 0 (see Appendix D.2).

Open with DEXTER
In the text
thumbnail Fig. 8

Ratio of the wave action flux above and below the critical level against the Rossby number and the inertial frequency for R >1∕4 (left panel) and R < 1∕4 (right panel) when the box is at the pole. As in Fig. 7, kx = kz = 0.1 and hatchedzones represent areas without critical levels. Again, crosses indicate the set of parameters chosen to analyse the behaviour of the velocity in the three-layer model.

Open with DEXTER
In the text
thumbnail Fig. 9

Numerical outputs of the three-layer model. Top: real part of the latitudinal velocity v against y. The quantities vnum, vI, vR, vT, vF, and vW are the numerical, incident, reflected, transmitted, first-order Frobenius, and Whittaker velocities, respectively. Bottom: wave action flux against y. The quantities , , , and are the numerical, transmitted, incident and reflected, and Frobenius wave action fluxes, respectively. For all panels, θ0 = 0, the mean flow is linear in the grey-shaded shear regions, and the critical level is marked by dashed lines. The horizontal wave numbers are set to kx = kz = 0.1. From left to right: (i) ω = 0.02 and Ro = 0.3 (R > 1∕4), (ii) ω = 0.002 and Ro = 0.8 (R > 1∕4), and (iii) ω = 0.09 and Ro = 1.8 (R < 1∕4).

Open with DEXTER
In the text
thumbnail Fig. 10

Same quantities as in Fig. 9 (when the box is at the pole), but for a box tilted by 10 degrees relative to the pole and for different values of the shear, wavenumbers, and inertial frequencies. We note that, unlike Fig. 9, we do not have analytical solutions in the shear region. In all panels, the horizontal wavenumbers are set to kx = kz = 1 and the shear is Ro = 0.3. From left to right: inertial frequencies set to ω = 0.31, ω = 0.16, and ω = 0.1. Third panels: we indicate y0 and y with dasheddotted and vertical dashed lines, respectively.

Open with DEXTER
In the text
thumbnail Fig. 11

Absolute values of the shear contrast normalised by the rotation at the pole against the age of K and G spectral type stars. Solid lines feature 0.3 < Rof < 0.9, dashed lines Rof ≲ 0.3, and transparent lines Rof > 1.

Open with DEXTER
In the text
thumbnail Fig. 12

Rossby numbers for Jupiter and Saturn as a function of the axial distance r.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.