© H. Dhouib et al. 2021

1. Introduction

The traditional approximation of rotation (TAR) was introduced for the first time in geophysics by Eckart (1960) to study the dynamics of the shallow Earth atmosphere and oceans. Then, Berthomieu et al. (1978) used it for the first time in an astrophysical context. Afterwards, Bildsten et al. (1996) and Lee & Saio (1997) applied it to stellar pulsations to study low-frequency nonradial oscillations in rotating stars, more precisely the gravito-inertial waves (GIWs). These waves (often called g modes) propagate under the combined action of the buoyancy force and the Coriolis acceleration. A subclass of these waves are the Rossby waves (often called r modes) whose dynamics are driven by the conservation of the vorticity of the fluid.

The TAR is applicable for low-frequency waves propagating in strongly stratified zones. In this case, a hierarchy of frequencies is verified. Namely, the buoyancy force is stronger than the Coriolis force (i.e., 2Ω ≪ N, where Ω is the angular velocity and N is the Brunt-Vaïsälä frequency) in the direction of stable entropy or chemical stratification and the Brunt-Vaïsälä frequency is larger than the frequency of the waves in the co-rotating frame ω (i.e., ω ≪ N). Considering these two assumptions, we can neglect the horizontal projection of the rotation vector Ω_H = Ω sin θ in the linearised hydrodynamic equations, where θ is the colatitude. Then, the vertical component of the Coriolis acceleration can be neglected and the wave velocity is almost horizontal within the description of the dynamics of GIWs. In this framework, the linearised hydrodynamic system, which is coupled in the general case (Dintrans et al. 1999; Dintrans & Rieutord 2000), becomes separable as in the non-rotating case and can be rewritten in the form of the Laplace tidal equation (Laplace 1799).

When the GIWs are observed with high-precision photometry, they allow us to study the internal structure and dynamics of stars with asteroseismology (Aerts et al. 2010; Aerts 2021). GIWs also trigger transport of angular momentum and mixing inside the stars, modifying their evolution (e.g., Mathis et al. 2008; Mathis 2009; Rogers et al. 2013; Rogers 2015; Aerts et al. 2019; Neiner et al. 2020). In this context, GIWs which propagate in stably stratified stellar radiation regions are of major interest. First, they are detected and analysed with asteroseismology in intermediate- and high-mass stars (e.g., Neiner et al. 2012; Van Reeth et al. 2015) allowing the characterisation of the internal stellar rotation (Van Reeth et al. 2016, 2018). In addition, GIWs are one of the candidates to explain the small rotation contrast observed between the stellar surface and core (Mathis 2009; Rogers 2015; Aerts et al. 2019). So, the TAR and its flexibility allow us to derive powerful seismic diagnostics, based on the period spacings between consecutive radial order g modes in uniformly and differentially rotating spherical stars (Bouabid et al. 2013; Ouazzani et al. 2017; Van Reeth et al. 2018). This enables us to determine the properties of the chemical stratification and the rotation rate near the convective core of early-type stars, which constrain our modelling of stellar structure, evolution and internal angular momentum transport (e.g., Pedersen et al. 2018; Aerts et al. 2018, 2019; Ouazzani et al. 2019).

In addition to the two conditions described above, three other assumptions are made when applying the TAR. First, the rotation is assumed to be uniform. Second, the star is assumed to be spherical, in other words, the centrifugal deformation of the star is neglected (i.e., , where Ω_K is the Keplerian critical (breakup) angular velocity, and G, M, and R are the universal constant of gravity, the mass of the star, and the stellar radius, respectively). Third, the Cowling approximation (Cowling 1941) is assumed where the fluctuations of the gravitational potential caused by the redistribution of mass by waves motion is neglected.

In addition, the anelastic approximation (Spiegel & Veronis 1960) is assumed where the high-frequency acoustic modes can be filtered out. Finally, for asymptotic low frequency modes which are rapidly oscillating in the vertical direction, we can assume the Jeffreys, Wentzel, Kramers, Brillouin (JWKB) approximation (e.g., Fröman & Fröman 1965).

The assumption of uniform rotation was abandoned by Ogilvie & Lin (2004) and Mathis (2009), who took into account the effects of general differential rotation (both in radius and latitude in a spherical star) on low-frequency GIWs within the framework of the TAR. The results of this new formalism have been successfully applied in Van Reeth et al. (2018) to derive the variation of the asymptotic period spacing in the case of a weak radial differential rotation and evaluate the sensitivity of GIWs to the effect of differential rotation.

The assumption of spherical symmetry was partially abandoned by Mathis & Prat (2019), who generalised the TAR for slightly deformed moderately rotating stars by considering the effects of the centrifugal acceleration using a perturbative approach. First, they derived the centrifugal perturbation of the spherical hydrostatic balance by determining the deformation of an isobar in the case of moderate uniform rotation. Afterwards, using this deformation they defined a mapping between the spherical coordinates system and a spheroidal coordinates system which takes into account this weak centrifugal perturbation of the star. Then, by keeping only first-order terms in the deformation, they re-derived the linearised hydrodynamic equations (Saio 1981; Lee & Baraffe 1995) and deduced the perturbed Laplace tidal equation within this new coordinates system. Then, they carried out a numerical exploration of the eigenvalues and horizontal eigenfunctions (the so-called Hough functions; Hough 1898) of the perturbed Laplace tidal equation. Finally, they derived the asymptotic expression for the frequencies of GIWs, including the effect of the centrifugal acceleration. This perturbative framework was then optimised for practical intense seismic forward modelling (Aerts 2021) using one-dimensional stellar structure models by Henneco et al. (2021).

In this work, we consider the general case of rapidly and uniformly rotating stars and planets which will generalise for the first time the previous work of Lee & Saio (1997) and Mathis & Prat (2019). The shape of such objects becomes a strongly deformed spheroid because of the action of the centrifugal acceleration which cannot be modelled as a perturbation in this case. We study if the TAR could be generalised to include the effects of this deformation for any stellar (planetary) rotation including those close to their break-up velocity. The derived formalism could have several key applications. For instance, it could be used to build new seismic diagnostic and to study angular momentum transport (e.g., Lee & Saio 1993; Mathis et al. 2008; Mathis 2009; Lee et al. 2014) and tidal dissipation (e.g., Ogilvie & Lin 2004, 2007; Braviner & Ogilvie 2014) induced by low-frequency GIWs in deformed stars and planets.

To do so, we first introduce in Sect. 2 a new spheroidal coordinate system which follows the shape of a deformed star for any rotation. It is based on the work by Bonazzola et al. (1998) and was first used by Rieutord et al. (2005) to develop rotating polytropic models. Then Lignières et al. (2006) employed it to study acoustic oscillations in strongly deformed rotating stars. Afterwards, Reese et al. (2006, 2009) and Ouazzani et al. (2012) make use of it to develop the two-dimensional oscillation program (TOP) and the adiabatic code of oscillation including rotation (ACOR), respectively. They compute non-radial pulsations of rotating stars. Finally, Espinosa Lara & Rieutord (2013) utilised it to develop the Evolution STEllaire en Rotation (ESTER) code which computes the structure of an early-type star including its differential rotation and the associated meridional circulation. Here, we first consider the case of uniform rotation to distinguish the effects of the deformation from those of differential rotation associated to the centrifugal and Coriolis accelerations, respectively. From there, we derive the system of linearised hydrodynamic equations in spheroidal coordinates. In Sect. 3, we apply the TAR in this new coordinate system by adopting the adequate assumptions. In Sect. 4, we can then rewrite the oscillation equations in the form of a new generalised Laplace tidal equation (GLTE) and deduce the asymptotic frequencies of low-frequency GIWs. In Sect. 5, we identify the validity domain of the TAR when using ESTER models and carry out a numerical exploration of the eigenvalues and the eigenfunctions of the GLTE. In Sect. 6, we quantify the centrifugal deformation effect of the star on the asymptotic period spacing pattern and we discuss its detectability. In Sect. 7, we verify the validity of the hierarchy of the terms in the linearised hydrodynamics equations imposed by the TAR in deformed stars. Finally, we discuss and summarise our work and results in Sect. 8.

2. Hydrodynamic equations in uniformly rotating deformed stars and planets

2.1. Spheroidal geometry

Because of the flattened shape of a rapidly rotating star, we use the spheroidal coordinate system (ζ, θ, φ) proposed by Bonazzola et al. (1998), where ζ is the radial coordinate (the pseudo-radius), θ the colatitude and φ the azimuth. It was designed not only to fit the shape of the star at the surface but also to get closer to the spherical shape at the centre as illustrated in Fig. 1. Following Rieutord et al. (2013), this new coordinate system can be linked to the usual spherical one (r, θ, φ) via the following mapping

(1)

Fig. 1. Spheroidal coordinate system with n = 13 subdomains used to compute the equilibrium model of a star rotating at 60% of its break-up velocity corresponding to the case of the 3 M⊙ ESTER model (with the central fraction in Hydrogen Xc = 0.7) studied in Sect. 5. ζ and θ are the pseudo-radius and the colatitude, respectively. ηi and ΔRi(θ) are defined in our mapping (Eq. (1)).

where we split the spheroidal domain 𝒟 in n subdomains 𝒟_{i ∈ ⟦0, n − 1⟧} ∈ [R_i(θ),R_i + 1(θ)] where R_{i ∈ ⟦0, n⟧}(θ) are series of functions, such that R_n(θ) = R_s(θ) is the outer boundary and R₀(θ) = 0 is the centre. Additionally, η_i = R_i(θ = 0) are the polar radii of the interfaces between the subdomains, Δη_i = η_i + 1 − η_i, ΔR_i(θ) = R_i + 1(θ)−R_i(θ) and ξ = (ζ − η_i)/Δη_i. The functions A_{i ∈ ⟦1, n − 1⟧}(ξ) = − 2ξ³ + 3ξ² and A₀(ξ) = − 1.5ξ⁵ + 2.5ξ³ and the constants a_i (in practice a_i = 1) are chosen to satisfy the boundary conditions between the different subdomains so that the mapping is continuous and derivable. The other coordinates (θ and φ) remain unchanged. Since we choose the polar radius as the unit length R_s(θ = 0) = 1, ζ = 1 coincides with the surface (i.e., r(ζ = 1, θ) = R_s(θ)), while r = 0 is obtained at the centre when ζ = 0 (i.e., r(ζ = 0, θ) = 0). Along the rotation axis (θ = 0), the deformed radial coordinate ζ is the distance from the centre normalised by the polar radius and when the centrifugal deformation is zero, ζ is equivalent to the radius, hence the name pseudo-radius.

Once the spheroidal coordinate system is established, it is mandatory to specify the corresponding basis. To do so, we define the spheroidal basis (a_ζ, a_θ, a_φ) from the natural covariant basis (E_ζ, E_θ, E_φ) defined as E_i = ∂_i(re_r)

(2)

where r_ζ ≡ ∂_ζr, r_θ ≡ ∂_θr and (e_r, e_θ, e_φ) is the usual spherical unit-vector basis. This basis has the advantage of being reduced to the usual spherical unit-vector basis when we approach the spherical geometry (R_s(θ) tends to the constant function 1).

2.2. Linearised hydrodynamic equations in spheroidal coordinates

Once the geometry is fixed, we can derive the complete adiabatic inviscid system of linearised equations to treat the wave dynamics in a uniformly rotating, strongly deformed star or planet. First, the linearised momentum equation for an inviscid fluid is written as (we refer the reader to Appendix A.1 for a detailed derivation)

(3)

(4)

(5)

where v^ζ, v^θ and v^φ are the contravariant components of the velocity field v = vⁱa_i and Ω = Ω(cos θe_r − sin θe_θ) is the uniform angular velocity. ρ, Φ and P are the fluid density, gravitational potential and pressure, respectively. Each of these scalar quantities has been expanded as

where X₀ is the hydrostatic component of X and X′ the wave’s associated linear fluctuation. When assuming the Cowling (1941) approximation, the fluctuation of the gravitational potential is neglected (Φ′ = 0).

Subsequently, the linearised continuity equation is obtained (we refer the reader to Appendix A.2 for a detailed derivation)

(6)

Then, the linearised energy equation in the adiabatic limit is derived (we focus in this work on adiabatic oscillations as in studies of the TAR in spherical stars)

(7)

where Γ₁ = (∂lnP₀/∂lnρ₀)_S (S being the macroscopic entropy) is the adiabatic exponent, is the background effective gravity which includes the centrifugal acceleration and N² is the squared Brunt–Väisälä frequency given by

(8)

Finally, since in this work we are interested in studying standing eigenmodes of oscillation we consider a discrete spectrum of eigenfrequencies rather than an integral over a continuous spectrum. So, we expand the wave’s velocity and scalar quantities fluctuations (X′≡{ρ′,P′}) on discrete Fourier series both in time and in azimuth

(9)

(10)

where m is the azimutal order (m is an integer because we are looking for a linear solution of a wave problem whose axisymmetric background depends only on ζ and θ) and ωⁱⁿ is the wave eigenfrequency in an inertial reference frame (we don’t put an index for the eigenfrequencies in order to make the notations simpler). In a uniformly rotating region, the waves are Doppler shifted due to the rotation. We can define the wave eigenfrequency ω in the co-rotating reference frame as

(11)

3. Generalised TAR

3.1. Approximation on the stratification profile: N²(ζ, θ) ≈ N²(ζ)

To obtain a separable system of equations by applying the TAR, a first necessary condition is that the Brunt-Väisälä frequency N² depends mainly on the pseudo-radius ζ, that is to say it has a small variation with respect to the colatitude θ. This also implies that the hydrostatic pressure P₀ and the hydrostatic density ρ₀ depend mainly on ζ so the linearised energy equation simplifies into

(12)

In order to compute the relative error made by adopting this approximation and to define its validity domain, we compare in Sect. 5.2 an approximate value of N² where it depends only on ζ to the exact one depending on ζ and θ by fixing the value of θ and varying the angular velocity of the star using two-dimensional ESTER stellar models.

3.2. The TAR with centrifugal acceleration

Using this first approximation, and assuming the same hierarchies of frequencies (2Ω ≪ N and ω ≪ N) and velocity scales (|v^ζ|≪{|v^θ|,|v^φ|}) that in the spherical and in the weakly deformed cases (this will be discussed in details in Sect. 7), we now built the generalised framework for the TAR in the case of a uniformly and rapidly rotating strongly deformed star (planet).

By adopting the approximation N²(ζ, θ) ≈ N²(ζ) and the Cowling approximation, we can rewrite the radial momentum Eq. (3) as

(13)

where is the normalised pressure. The propagation of low-frequency GIWs can be studied within the anelastic approximation in which acoustic waves are filtered out (Dintrans & Rieutord 2000; Mathis 2009). Therefore, Eqs. (13) and (12) can be simplified accordingly by neglecting the terms and . Within the TAR, we focus on low-frequency waves (i.e., ω ≪ N) propagating in strongly stratified regions (i.e., 2Ω ≪ N). In this case, we can neglect the vertical component of the Coriolis acceleration because it is dominated by the buoyancy force and the vertical wave velocity since |v^ζ|≪{|v^θ|,|v^φ|} because of the strong stable stratification. Thus, the radial momentum equation can be simplified into

(14)

where

(15)

Subsequently, we examine the latitudinal component of the momentum equation (Eq. (4)), which reduces, for the same reasons, to

(16)

where

(17)

(18)

Finally, the azimuthal component of the momentum equation (Eq. (5)) simplifies onto

(19)

where

(20)

We thus obtain a system of equations for strongly deformed rotating stars which has the same mathematical form that in the case of spherically symmetric stars (Lee & Saio 1997) and weakly deformed stars (Mathis & Prat 2019), and where we thus still manage to partially decouple the vertical and horizontal components of the velocity. By solving the system formed by Eqs. (16) and (19) we can express u^θ and u^φ as a function of as follows

(21)

(22)

where we introduce

(23)

and the spin parameter

(24)

The structure of the new equations that we obtain in the spheroidal case is similar to the one in the usual spherical case (e.g., Lee & Saio 1997; Townsend 2003; Pantillon et al. 2007). The most important difference is that the coefficients 𝒜, ℬ, 𝒞, 𝒟 and ℰ are function of ζ and θ through r(ζ, θ) and their derivatives. As this is done in Lee & Saio (1997) and Mathis & Prat (2019), we introduce the reduced latitudinal coordinate x = cos θ, so that Eqs. (21) and (22) transform to

(25)

(26)

4. Dynamics of low-frequency gravito-inertial waves

Our goal in this section is to derive the generalised Laplace tidal equation (GLTE) for the normalised pressure , which allows us to compute the frequencies and periods of low-frequency GIWs and to build the corresponding seismic diagnostics following the method by Bouabid et al. (2013), Van Reeth et al. (2018), and Mathis & Prat (2019). Applying the anelastic approximation and the approximation N²(ζ, θ) ≈ N²(ζ) in the continuity equation (Eq. (6)) allows us to neglect the terms ∂_tρ′ and ∂_θρ₀ respectively, so it simplifies into

(27)

4.1. JWKB approximation

Assuming that ω ≪ N, each scalar field and each component of u can be expanded using the two-dimensional Jeffreys-Wentzel-Kramers-Brillouin (JWKB) approximation which was introduced for strongly stratified differentially rotating fluids by Mathis (2009). In this case, the vertical wave number k_V (≡k_ζ) is assumed to be very large, and the associated wavelength is thus very small. Therefore, the spatial variation of the wave in the pseudo-radial direction is very fast compared to that of the hydrostatic background given by ρ₀, P₀ and g_eff. So the spatial structure of the waves can be described by the product of a rapidly oscillating plane-like wave function in the pseudo-radial direction multiplied by a bi-dimensional slowly varying envelope:

(28)

(29)

with j ≡ {ζ, θ, φ}, k is the index of a latitudinal eigenmode (cf. Sect. 4.3) and A_νkm is the amplitude of the wave. Substituting the expansion given in Eqs. (28) and (29) into Eqs. (14), (25) and (26), the final pseudo-radial, latitudinal and azimuthal components of the velocity are obtained:

(30)

(31)

(32)

4.2. Approximation: 𝒜(ζ, θ) ≈ 𝒜(ζ)

To be able to derive a generalised Laplace tidal equation as in the case of slightly deformed spheres (Mathis & Prat 2019) and to introduce its horizontal eigenvalues Λ_νkm, which should depend only on ζ, we have to do a partial separation between the pseudo-radial and latitudinal variables in the pseudo-radial component of the velocity . In other words, we have to write Eq. (30) as the product of a function which depends only on ζ and the normalised pressure w_νkm. To do so, we have to assume that the coefficient 𝒜(ζ, θ) depends mainly on ζ. To evaluate the relative error made by adopting this approximation and to define its validity domain, we compare in Sect. 5.2 an approximate value of 𝒜 where it depends only on ζ to the exact one depending on ζ and θ by fixing the value of θ and varying the angular velocity of the star using two-dimensional ESTER stellar models.

4.3. Generalised Laplace tidal equation (GLTE)

Substituting the radial component of the velocity (Eq. (30)) into the continuity equation (Eq. (27)) we get:

(33)

where the JWKB approximation allows us to neglect ζ∂_ζρ₀/ρ₀u^ζ and 2u^ζ in front of the dominant term ζ∂_ζu^ζ. Using Eqs. (31) and (32), we obtain the GLTE for the normalised pressure w_νkm

(34)

When the centrifugal acceleration is taken into account, the eigenvalues Λ_νkm and the generalised Hough functions (eigenfunctions) w_νkm of the GLTE vary with the pseudo-radius ζ. We choose to define our latitudinal ordering number k to enumerate, for each (ν, m), the infinite set of solutions as in Lee & Saio (1997) and Mathis & Prat (2019) by considering the eigenvalues and eigenfunctions at the centre where they are not affected by the centrifugal acceleration since our mapping (Eq. (1)) is such that r(ζ → 0, θ)→0.

We can check the validity of our model analytically by calculating the coefficients of the GLTE (Eq. (34)) in the spherical case (cf. Table 1) and by substituting these values into the GLTE. We then recover the standard Laplace tidal equation (SLTE) of Lee & Saio (1997):

(35)

In this case, the eigenvalues are independent of the pseudo-radius ζ. We note that at the limit ν = 0, which corresponds to the non-rotating case, the SLTE reduces to the differential equation for the associated Legendre polynomials P_ℓ, m of degree ℓ and order m with and w_0km(x) = a_ℓ, mP_ℓ, m(x) where ℓ = |m|+k and

is the normalisation coefficient.

In our work, we solve for each spheroid only one linear ordinary differential equation (the generalised Laplace tidal equation) which depends on the colatitude and in a parametric way on the pseudo-radius. Indeed, since we focus on low-frequency gravito-inertial waves, which are rapidly oscillating along the pseudo-radial direction, we adopt the JWKB approximation. It allows us to describe the spatial structure of the waves by the product of a rapidly oscillating plane-like wave function in the pseudo-radial direction multiplied by a bi-dimensional slowly varying envelope (Eqs. (28) and (29)) which is obtained by solving the GLTE (Eq. (34)) for each spheroidal shell for the normalised pressure and polarisation relationships for the different components of the velocity (Eqs. (30)–(32)). Therefore, we avoid having an ordinary differential equation to solve in the pseudo-radial direction. This method has been introduced and used by Mathis (2009) and Van Reeth et al. (2018) in the case of the generalisation of the TAR to the case of spherical differentially rotating stars and by Mathis & Prat (2019) and Henneco et al. (2021) in the case of the generalisation of the TAR for slightly deformed rotating stars. Finally, solving the GLTE allows us to compute the horizontal eigenvalues and the corresponding eigenfrequencies.

4.4. Asymptotic frequency and period spacing of low-frequency GIWs

From Eqs. (33) and (34), we can identify the dispersion relation for low-frequency GIWs within the TAR for strongly deformed stars

(36)

To be able to derive the eigenfrequencies of low-frequency GIWs, we applied the quantisation relation in the vertical (pseudo-radial) direction used by Mathis (2009) and Mathis & Prat (2019) and introduced by Unno et al. (1989), Gough et al. (1993), and Christensen-Dalsgaard (1997)

(37)

where ζ₁ and ζ₂ are the turning points of the Brunt–Väisälä frequency N and n is the radial order.

So, by substituting the dispersion relation (Eq. (36)) into the vertical quantisation one (Eq. (37)), we get the asymptotic expression for the frequencies of low-frequency GIWs

(38)

and the corresponding period

(39)

with

(40)

We can thus compute the period spacing ΔP_km = P_n + 1km − P_nkm. We note that by applying a first-order Taylor development we can generalise the result obtained by Bouabid et al. (2013) in the spherical case

(41)

which will potentially allow us to probe the internal rotation of a large number of rapidly rotating stars taking into account their possibly strong flattening by the centrifugal acceleration and perform their detailed seismic modelling.

This result generalises and is very similar to the one obtained in the spherical case (Lee & Saio 1997; Bouabid et al. 2013) and in the weakly deformed case using a perturbative approach (Mathis & Prat 2019). The major interest of the eigenfrequencies, periods and periods spacing (Eqs. (38)–(41)) is that they can be easily evaluated for a large number of rapidly rotating stars and compared to observations which is a great asset to perform detailed seismic modelling. This flexibility is what makes the TAR a very interesting treatment.

5. Application to rapidly rotating early-type stars

In order to implement our equations and solve the GLTE, we must first compute an equilibrium model which gives the structure of the rotating star. The most advanced model of stellar evolution that we have and which is adapted to the formalism that we developed here is the ESTER model (Espinosa Lara & Rieutord 2013). It is unique of its kind because it calculates the effect of the centrifugal acceleration in a non-perturbative and a self-consistent way at any rotation rate.

5.1. ESTER models

The ESTER code is a 2D stellar structure code which takes into account the effects of rotation in a non-perturbative and coherent manner. The ESTER code makes it possible to calculate the equilibrium state of a star with a mass greater than 2 M_⊙ by solving the stationary equations of hydrodynamics for pressure, density, temperature and angular velocity for the entire volume. These realistic two-dimensional models of rapidly rotating stars are computed with the actual baroclinic meridional flows and differential rotation for intermediate-mass and massive stars. In this work, as a first step we work only with a uniform rotation (the mean value of the differential rotation over the pseudo-radius and the colatitude). The inclusion of differential rotation in the TAR formalism will be done in Paper II. The code uses spectral methods in the radial and latitudinal directions to transform the original system of non-linear partial differential equations in a system of non-linear algebraic equations (Chebyshev polynomials in the pseudo-radial direction and with spherical harmonics in the horizontal direction) (Rieutord et al. 2013).

In this study, we use 3 M_⊙ stellar models with a hydrogen mass fraction in the core X_c = 0.7 at [0%,90%] of the Keplerian break-up rotation rate Ω_K. Such a model is representative of observed stars near the zero-age main-sequence with asteroseismic probes of their angular velocity (Van Reeth et al. 2015, 2016, 2018; Pápics et al. 2015, 2017; Li et al. 2019a,b).

5.2. Domain of validity of the TAR

In this section, we study the validity domain of the two approximations (3.1 and 4.2) that are necessary to build the TAR in the case of rapidly rotating deformed bodies:

(42)

(43)

To visualise the problem, we first illustrate in Figs. 2 and 3 the functions N²(ζ, θ) and 𝒜(ζ, θ) computed with an ESTER model rotating at Ω/Ω_K = 60% (we use r instead of ζ in these figures so we can show the spheroidal shape of the star caused by the centrifugal deformation). We can see from Fig. 2 that the convective region of the star is located between ζ = 0 and ζ = 0.153 and that the pseudo-radial variation of the Brunt–Väisälä frequency is the highest at the surface. From Fig. 3 we can see that the coefficient 𝒜 varies the most, according to the pseudo-radius, at the centre of the star where we approach its singularity. In order to be able to apply our approximations (Eqs.(42) and (43)), we are interested in quantifying the latitudinal variation of this two quantities.

Fig. 2. Squared Brunt–Väisälä frequency N2(ζ, θ) profile for Ω/ΩK = 60%.

Fig. 3. Coefficient 𝒜(ζ, θ) profile for Ω/ΩK = 60%.

So to evaluate the committed error by adopting Eqs. (42) and (43) we represent in Figs. 4 and 5 the Brunt–Väisälä frequency and the coefficient 𝒜 pseudo-radial profiles for different values of the colatitude θ_f, at a fixed rotation rate Ω = 0.6Ω_K. In the bottom panels, we show the corresponding relative error

(44)

Fig. 4. Profile of the squared Brunt–Väisälä frequency N2 and the relative error of the adopted approximation (Eq. (42)) as a function of ζ at different colatitudes θf using an ESTER model rotating at 60% of the Keplerian breakup rotation rate (The light orange area indicates the margin of error which we allow and the red hatched area represents the convective region of the star).

Fig. 5. Profile of the coefficient 𝒜 and the relative error δ𝒜θf of the adopted approximation (Eq. (43)) as a function of ζ at different colatitudes θf using an ESTER model rotating at 60% of the Keplerian break-up rotation rate (The light orange area indicates the margin of error which we allow).

between the exact value X(ζ, θ_f) given by the model and the approximated model X_approx(ζ) given by the weighted average of the exact value over the colatitude θ

(45)

with

where θ_f is a fixed value of the colatitude θ.

Then, in order to determine the validity domain of these approximations as a function of the rotation rate and the pseudo-radius, we vary the rotation rate of the star Ω/Ω_K (we use here 10 sample points for the rotation rate) and we calculate for each case the associate maximum relative error max_iδX_θi(ζ). Afterwards, we fix a maximum error rate equal to 10% and we deduce the pseudo-radius limit ζ_limit where the maximum relative error exceeds this threshold and we decide that the approximations become invalid. In Fig. 6, we display the pseudo-radius limit ζ_limit as a function of the rotation rate Ω/Ω_K. This curve defines the validity domains of the two approximations (Eqs. (42) and (43)) which consequently define the validity domain of the TAR.

Fig. 6. Validity domain of the TAR within the framework of ESTER models (3 M⊙, Xc = 0.7) as a function of the rotation rate Ω/ΩK and the pseudo-radius ζ (with 90% degree of confidence).

We notice that the validity domains of the two approximations are more and more restricted (the pseudo-radius limit ζ_limit decrease) by increasing the rate of rotation and approaching the critical angular velocity. The tendency of the reduction of the pseudo-radius limit with the increasing of the rotation rate is expected since the action of the centrifugal acceleration is stronger further away from the axis of rotation. Thus, to be able to apply the TAR in the whole space domain, we shouldn’t exceed a rotation rate Ω = 0.2Ω_K. It is therefore possible to apply the generalised TAR to early-type stars (3 M_⊙, X_c = 0.7) with an external radiative zone rotating up to 20% of the Keplerian critical angular velocity. This limit is lower than the one proposed by Mathis & Prat (2019) which was Ω = 0.4Ω_K. The limit proposed here should be more realistic since it is derived from a two-dimensional model which takes into account the centrifugal acceleration in a non-perturbative way whereas Mathis & Prat (2019) treated it in a perturbative manner. In the case where ζ_limit < 1, generalised TAR should also be applicable to the radiative core of rapidly rotating (young) solar-type stars for which ζ_rad < ζ_limit, where ζ_rad is the pseudo radius of the radiation-convection interface. So, this formalism can be applied to all rapidly rotating bodies (stars of all types and planets) as long as we have a 2D model with access to the quantities that we used to derive our equations.

5.3. Eigenvalues and Hough functions

We solve the GLTE for different pseudo-radii, spin parameters, and rotation rates within the defined validity domain using the implementation based on Chebyshev polynomials presented and discussed by Wang et al. (2016) and then used by Mathis & Prat (2019) and Henneco et al. (2021).

Figure 7 shows the eigenvalues Λ_νkm of the GLTE at ζ = ζ_limit as a function of the spin parameter ν for m = 1 at two different rotation rates: Ω = 0.2Ω_K the maximum rotation rate for which the generalised TAR is valid in all the radiative region (ζ_limit = 1) and Ω = 0.9Ω_K the highest rotation rate that we consider in this study for which the generalised TAR is valid in a restricted space domain (ζ_limit = 0.57). Since Ω points in the direction of θ = 0 and the oscillations are proportional to e^{i(ωt − mφ)}, a prograde (retrograde) oscillation has a positive (negative) value of the product mν. Therefore, since m = 1 the prograde and retrograde modes correspond respectively to the positive and negative values of ν. We see here that the centrifugal deformation of the star causes a modest gradual shift in the eigenvalues as already identified in previous studies where the centrifugal acceleration has been treated in a perturbative way (Mathis & Prat 2019; Henneco et al. 2021).

Fig. 7. Spectrum of the GLTE as a function of the spin parameter ν at ζ = ζlimit and m = 1 for Ω = 0.2ΩK (above) and Ω = 0.9ΩK (below). Blue (respectively, orange) dots correspond to even (respectively, odd) eigenfunctions.

We can distinguish two families of eigenvalues. Namely, the first one is gravity-like solutions (Λ_νkm ≥ 0) which exist for any value of ν and we attach positive k’s to them. They correspond to gravity waves (g modes) modified by rotation. The second one is Rossby-like solutions which exist only when |ν|> 1 and we attach negative k’s to them. They only propagate in rotating stars. They correspond to Rossby modes (r modes) if they are retrograde and have positive eigenvalues (Λ_νkm > 0) (Saio et al. 2018) and to overstable convective modes if they are prograde and have negative eigenvalues (Λ_νkm < 0). These modes are able to propagate in convective regions stabilised by the rotation under the action of the Coriolis acceleration (Lee & Saio 1997, 2020; Lee 2019). Even though the TAR is not a valid approximation in convective zones (Ogilvie & Lin 2004), we do find such propagative modes as solutions from our approach. Contrary to the spherical case we obtain here eigenvalues which depend slightly on the pseudo-radius ζ as shown in Fig. 8 where we represent the spectrum of the GLTE as a function of the pseudo-radius at two different rotation rates Ω/Ω_K = {20%,90%} for m = 1 and ν = 10.

Fig. 8. Spectrum of the GLTE as a function of the pseudo-radius for m = 1 and ν = 10 at Ω/ΩK = 0.2 (above) and Ω/ΩK = 0.9 (below). Blue (respectively, orange) dots correspond to even (respectively, odd) Hough functions.

We focus now on the generalised Hough functions w_νkm (normalised to unity) which vary with the pseudo-radius ζ and the horizontal coordinate x. This dependence is illustrated at ν = 10 in Fig. 9 at two different rotation rates Ω/Ω_K = {20%,90%} for gravity like-solutions (prograde dipole {k = 0, m = 1} modes) and Rossby like-solutions (retrograde Rossby modes with {k = −2, m = −1}). We can see that increasing the rotation rate from 20% to 90% of the critical angular velocity does hardly affect the Hough functions. This is caused by the fact that increasing the rotation rate implies decreasing the pseudo-radius limit (ζ_limit) so we can no longer trust the computed Hough functions near the surface where the centrifugal effect is maximal. We can see clearly in this figure that the centrifugal acceleration slightly influences and modifies the eigenfunctions of the GLTE by introducing a new dependency on the pseudo-radial coordinate ζ, whereas in the spherically symmetric case they were only dependant on the latitudinal coordinate θ. In fact, the gravity-like solutions and the Rossby-like solutions shift as the distance from the border of the radiative zone of the star (ζ = 0.153) to its surface (ζ = 1) increases. More specifically, we observe a non-monotonous behaviour in Rossby-like solutions: these solutions migrate at first onwards, away from the equator (x = 0), causing a broadening of its shape. Next, by increasing the reduced latitudinal coordinate x this behaviour becomes less and less visible until it changes to the opposite one where these solutions migrate inwards, towards the equator, causing a narrowing of its shape. However, in gravity-like solutions we found a monotonous behaviour for which these solutions migrate inwards, towards the equator, causing a narrowing of its shape. Contrary to the non-monotonous behaviour found in Rossby-like solutions, Henneco et al. (2021) found using the perturbative TAR a monotonous one with a migration towards the equator in the two types of solutions.

Fig. 9. Generalised Hough functions (normalised to unity) as a function of the horizontal coordinate x at different pseudo-radii from ζ = 0.153 (blue) to ζ = ζlimit (yellow) for ν = 10, Ω/ΩK = 0.2 (above), and Ω/ΩK = 0.9 (below). The solid lines correspond to gravity-like solutions with {k = 1, m = 1}, while the dotted lines correspond to Rossby-like solutions with {k = −2, m = −1}.

6. Asymptotic seismic diagnosis

6.1. Asymptotic period spacing pattern

To compute the period spacing patterns, we first calculated the asymptotic frequencies by following the method developed by Henneco et al. (2021) (the detailed steps of this computation can be found in their AppendixB). First, we compute for each radial order n of a given mode (k, m), the corresponding ν_nkm = 2Ω/ω_nkm using the implicit equation of the asymptotic frequencies of low-frequency GIWs (Eq. (38)). Then, by taking their inverse and multiplying them with 2Ω, we retrieve the asymptotic frequencies ω_nkm. The corresponding asymptotic frequencies in the inertial (observer’s) frame are then found using the angular Doppler shift (Eq. (11))

(46)

Now, we can easily calculate the asymptotic periods in the inertial frame () and the corresponding period spacing which is shown in Fig. 10 for the prograde dipole sectoral {k = 0, m = 1} and the retrograde Rossby {k = −2, m = −1} modes using the standard TAR (in spherical geometry) and the generalised TAR. The asymptotic periods within the standard TAR are computed using the following asymptotic frequencies deduced from the SLTE (Eq. (35))

(47)

Fig. 10. Period spacing pattern in the inertial frame for {k = 0, m = 1} (above) and {k = −2, m = −1} (below) modes at Ω/ΩK = 20% within the standard (red line) and the generalised (blue line) TAR frameworks. The bottom panel for each mode shows the differences with respect to the standard TAR period spacing pattern (the sudden fluctuations in the period spacing pattern are caused by the numerical noise which is introduced by the numerical derivatives of the mapping with respect to ζ and θ used in the resolution of the GLTE).

By comparing the period spacing pattern computed using the standard TAR and the generalised TAR in each of bottom panels of Fig. 10, we find that the spacing values of the {k = 0, m = 1} and {k = −2, m = −1} modes increase with the decrease in radial orders (low (high) radial orders equivalent to short (long) pulsation periods for g modes and long (short) pulsation periods for r modes) under the influence of the centrifugal deformation. The effect of the centrifugal acceleration on the period spacing that we found here is different from the one discussed in Henneco et al. (2021) who found an increase in the period spacing of the two modes. Another notable difference between the generalised TAR and the standard (Lee & Saio 1997) and the perturbative (Henneco et al. 2021) TAR is in the domain of validity and applicability, where the general treatment adopted in this work should be more reliable than the standard and the perturbative ones. Figure 10 reveals the difference in the period spacing values computed with the standard and the generalised formalisms and also offers a comparison of this difference with typical observational uncertainties of measured period spacing values from modern space photometric observations. In fact, Van Reeth et al. (2015) shows that for a sample of 40 γ Dor stars the error margins of the period spacings have values up to 25 seconds for the prograde dipole modes (larger than our period spacing differences) and between 20 and 50 seconds for the retrograde Rossby modes (smaller than our period spacing differences).

6.2. Detectability of the centrifugal deformation effect

To evaluate the detectability of the centrifugal effect with space-based photometric observations, we first compute the frequency differences Δf_centrifugal between asymptotic frequencies calculated in the standard TAR (TAR_s) and those computed in the generalised TAR (TAR_g) through Eq. (38):

(48)

Then, by comparing the obtained frequency differences with the frequency resolutions (f_res = 1/T_obs) of nominal 4-year Kepler (Borucki et al. 2009) and TESS CVZ (Transiting Exoplanet Survey Satellite Continuous Viewing Zone; Ricker et al. 2014) light curves covering quasi-continuously observation times of T_obs = 4 yrs and T_obs = 351 d, respectively, we can deduce the radial orders n for which the frequency differences are, at least in principle, expected to be detectable. In Fig. 11, we display these results for the {k = 0, m = 1} and {k = −2, m = −1} modes for a stellar model rotating at 0.2 Ω_K. We can see that the centrifugal effect for the {k = 0, m = 1} mode is detectable for all radial orders using Kepler observations, but detectable only for n < 41 using TESS observations (yellow dashed line in Fig. 11). Indeed even if for n ∈ [20, 40] using nominal TESS CVZ light curves and n ∈ [75, 100] using nominal Kepler light curves we have Δf_centrifugal/f_res ≳ 1, so we would expect that the detectability of the centrifugal effect for these radial orders is in principle guaranteed. However, since the frequency differences is nearly equal to the resolution one (Δf_centrifugal ≈ f_res), this effect become difficult to detect in this case. Regarding the {k = −2, m = −1} mode, we can see that the centrifugal effect is less detectable in comparison with the {k = 0, m = 1} mode but remains in principle observable (yet difficult to observe because Δf_centrifugal ≈ f_res) using Kepler observations for small radial orders (n ∈ [5, 8]) and for n ∈ [16, 45].

Fig. 11. Detectability of {k = 0, m = 1} (above) and {k = −2, m = −1} (below) modes as a function of the radial order n at a rotation rate Ω/ΩK = 20% based on the frequency resolution of TESS (brown line) and Kepler (indigo line).

Our results in Figs. 10 and 11 are formal comparisons in the sense that uncalibrated astrophysical effects occurring in stellar models (e.g., the presence of rotational mixing and atomic diffusion effects) usually imply larger uncertainties in the prediction of the period spacing values than the differences shown in these figures, as outlined in Aerts et al. (2018). For this reason, our theoretical findings in the current work offer a key result in terms of reliability of future asteroseismic analyses of observed pulsators. The centrifugal effect is not significant when assessed within the validity domain of the generalised TAR which is consistent with the results of Ouazzani et al. (2017) and Henneco et al. (2021).

7. Evaluation of the terms hierarchy imposed by the TAR within uniformly rotating deformed stars

The hierarchy of terms imposed by the TAR can be summarised by the following frequency hierarchy:

(49)

(50)

which ensures the other hierarchies. In fact, the vertical wave vector is larger than the horizontal one (|k_H|≪|k_V|) which comes from the dispersion relation of GIWs (Lee & Saio 1997)

(51)

In addition, the latitudinal component of the rotation vector in the momentum equation Ω_H = Ω sin θ is neglected for all colatitudes because (2Ω ⋅ k) ≈ 2Ω_Vk_V. And finally, the wave velocity is almost horizontal (|v^ζ|≪{|v^θ|,|v^φ|}). Indeed, by adopting the anelastic approximation we obtain

(52)

The anelastic approximation filters out acoustic waves which have higher frequencies, and is justified by the fact that we study only low frequency waves (Eq. (50)) within the framework of the TAR treatment.

So in order to verify whether the TAR is still valid in uniformly rotating deformed stars, we discuss the frequency hierarchy (Eqs. (49) and (50)) using the Brunt–Väisälä frequency and rotation profiles from ESTER models and using the asymptotic frequencies calculated in Sect. 6.

7.1. The strong stratification assumption: 2Ω ≪ N

We evaluate the term N/2Ω using the Brunt–Väisälä frequency and the angular velocity computed with ESTER models (3 M_⊙, X_c = 0.7) for rotation rates Ω/Ω_K ∈ [0.1,0.9]. As shown in Fig. 12, the strong stratification assumption is valid only in the radiative zone far from the border between convection and radiation (ζ ≥ 0.2 at Ω/Ω_K ≤ 0.2). Beyond this critical value, the Brunt–Väisälä frequency and the frequency of rotation have a very close order of magnitude so we can no longer adopt the strong stratification approximation for stars close to the zero-age main-sequence. As stars evolve, however, they build up a strong gradient of molecular weight, increasing the value of N(ζ) in the region near their convective core. Therefore, the strong stratification assumption becomes better fulfilled as stars evolve along the main sequence because the spin-up of their core remains modest (Li et al. 2020; Aerts 2021), while N(ζ) increases.

Fig. 12. N/2Ω term as a function of the pseudo-radius ζ at different rotation rates Ω/ΩK (the hatched area represents the convective region of the star).

7.2. The low frequency assumption: ω ≪ N

We evaluate the term N/ω using the Brunt–Väisälä frequency computed by ESTER models and the asymptotic frequencies for {k = 0, m = 1} and {k = −2, m = −1} modes calculated in Sect. 6. As shown in Fig. 13, the low frequency assumption is valid for {k = −2, m = −1} mode for all radial orders and pseudo-radii while for {k = 0, m = 1} mode this assumption is valid for all radial orders only if we are far from the transition layer between the convective core and the radiative envelope. If this is not the case, we have to work with high radial orders (n > 20). Below this critical value, the Brunt–Väisälä frequency and the wave frequency have a very close order of magnitude so we can no longer adopt the low frequency approximation.

Fig. 13. N/ω term as a function of the pseudo-radius ζ and the radial order n for {k = 0, m = 1} (top panel) and {k = −2, m = −1} (bottom panel) modes at Ω/ΩK = 20%.

8. Discussion and conclusions

The aim of this work was to study the possibility of carrying out a new generalisation of the TAR that takes into account the centrifugal acceleration in the case of strongly deformed, rapidly rotating stars and planets. We approached this issue by deriving the generalised Laplace Tidal equation in a general spheroidal coordinate system. We relied on two necessary assumptions (described by Eqs. (42) and (43)) that are applicable in a specific domain, which will then define the validity domain of the generalised TAR. The equation that we derive has a similar form to the one obtained when the TAR is applied to weakly rotating spherical stars (Lee & Saio 1997), differentially rotating spherical stars (Mathis 2009; Van Reeth et al. 2018) and moderately rotating, weakly deformed stars (Mathis & Prat 2019). This new formalism can be applied in the radiative region of all types of stars and planets.

Using this theoretical framework with 2D ESTER models we can define the validity domain of the TAR with a chosen confidence level as a function of the rotation rate and the pseudo-radius. We apply this general formalism to rapidly rotating early-type stars (M = 3 M_⊙,  X_c = 0.7), where we found that the generalised TAR is only applicable for all rotation rates for values of the pseudo-radius smaller than 0.57 (ζ_limit < 0.57). To calculate the asymptotic pulsation frequencies we then have to integrate the dispersion relation on the pulsation mode cavity (the radiative zone) (Eq. (37)) where ζ ∈ [ζ₁, ζ₂] (in our considered case ζ₁ = 0.153 and ζ₂ = 1). This condition will limit the rotation rates where the TAR is valid. We can therefore deduce a rotation rate limit (Ω_limit = 0.2Ω_K in our case) where the TAR ceases to be applicable with a fixed level of confidence (90% in our study).

We found a difference between the lsimit of the TAR applicability defined in Mathis & Prat (2019) which was Ω_limit = 0.4Ω_K and the one defined here. The current new limit should be more realistic and accurate since it is first derived from a 2D model which takes the impact of the centrifugal acceleration into account in a non-perturbative way. Next, the validity domain of the adopted approximations is studied carefully in more details and with a level of confidence of 90%.

The centrifugal acceleration affects two aspects within our new generalised theoretical formalism compared to the standard application of the TAR (Lee & Saio 1997). On the one hand, the stellar structure as a whole becomes deformed which causes the usual radial coordinate r, the hydrostatic background and thus the Brunt-Vaïsälä frequency to have a dependency not only on the pseudo-radius ζ but also on the colatitude θ. This introduces a new challenge in the decoupling of the vertical and the horizontal dynamics and the applicability of the quantisation relation unless we adopt adequate approximations which will restrict the validity domain of our formalism. On the other hand, the dynamics of the GIWs represented by the GLTE, is altered and get a radial dependence compared to the one in the spherical case represented by the SLTE. Then, we apply our new formalism to a rapidly rotating early-type star (M = 3 M_⊙,  X_c = 0.7,  Ω/Ω_K = 20%). In this case, we observe a non-monotonous geometrical behaviour in the Rossby-like solutions and a monotonous one (only an inward migration) in the gravity-like solutions. In fact, near the equator the Rossby-like solutions migrate onwards, causing a broadening of its shape. But far from the equator this behaviour becomes less and less visible until it changes to the opposite one where these solutions migrate inwards, causing a narrowing of its shape. We found also that the centrifugal acceleration increases the period spacing values for low radial orders (equivalent to short pulsation periods for g modes and long pulsation periods for r modes) and decrease slightly for high radial orders (equivalent to long pulsation periods for g modes and short pulsation periods for r modes). But yet the effects of the centrifugal acceleration (within the validity domain of the generalised TAR) are of the order of or smaller than the effects of other missing stellar physical processes in our models (e.g., accurate rotational mixing, atomic diffusion, etc.). This can mask the effect of the centrifugal acceleration in forward asteroseismic modelling analyses, even when it is theoretically detectable.

This new generalisation of the TAR can also be used to study the tidal dissipation induced by low-frequency GIWs in rapidly rotating deformed stars and planets (Braviner & Ogilvie 2014) and the angular momentum transport with a formalism directly implantable in ESTER models. In addition, we can apply our formalism to different 2D ESTER model grids and then compare the obtained results with direct computations using 2D oscillation codes TOP (Reese et al. 2021) and ACOR (Ouazzani et al. 2017).

The next steps will be first to take into account the general differential rotation in the deformed case. So far, Ogilvie & Lin (2004) and Mathis (2009) included the effects of general differential rotation on low-frequency GIWs within the framework of the TAR for spherical stars. Then, their formalism has been successfully applied in Van Reeth et al. (2018) to derive the variation of the asymptotic period spacing in the case of a weak radial differential rotation and evaluate the sensitivity of GIWs to the effect of differential rotation in spherical stars. Next, we will take into account the magnetic field in a non-perturbative way (Mathis & de Brye 2011, 2012). So far, Prat et al. (2019, 2020) and Van Beeck et al. (2020) have focused on the case where magnetic fields are weak enough to be treated within a perturbative treatment to study the effects of a magnetic field on the seismic parameters of g modes which become magneto-gravito inertial modes.

So in the near future, magneto-gravito inertial waves have to be studied in the general case of rapidly and differentially rotating stars and planets. As a first step, in a forthcoming paper (Dhouib et al. 2021, Paper II) we will examine if the TAR can be generalised to the case of strongly deformed differentially rotating stars and planets.

Acknowledgments

We thank the referee, Pr. H.Shibahashi, for his constructive report that allowed us to improve the quality of our work and of the article. We thank the ESTER code developers for their efforts and for making their codes publicly available. We are also grateful to Pr. C.Aerts who kindly commented on the manuscript and suggested improvements. H.D., V.P. and S.M. acknowledge support from the European Research Council through ERC grant SPIRE 647383. H.D. and S.M. acknowledge support from the CNES PLATO grant at CEA/DAp. T.V.R. gratefully acknowledges support from the Research Foundation Flanders (FWO) under grant agreement No. 12ZB620N.

References

Appendix A: Derivation of the linearised momentum and continuity equations in the spheroidal basis

A.1. The momentum equation

The linearised momentum equation for an inviscid uniformly rotating fluid is written as

(A.1)

Using the spheroidal basis (a_ζ, a_θ, a_φ) defined in Eq. (2) we can express the velocity field as follows

(A.2)

with , , , (E^ζ, E^θ, E^φ) is the dual basis (the contravariant basis) and g_ij is the covariant metric tensor defined as follows

(A.3)

We write also the uniform angular velocity in the natural covariant basis as

(A.4)

(A.5)

We then express the Coriolis acceleration in the spheroidal basis as follows

(A.6)

where ε^ilk is the Levi-Civita tensor

(A.7)

with and

(A.8)

Finally, we express the gradient of the scalar function f = {P′,P₀, Φ′} in the spheroidal basis as

(A.9)

With these definitions, it is now possible to give an explicit expression of the momentum equation (Eq. (A.1)) in the spheroidal basis

(A.10)

for i = ζ, i = θ and i = φ we get the pseudo-radial (Eq. (3)), the latitudinal (Eq. (4)), and the azimuthal (Eq. (5)) components of the momentum equation, respectively.

A.2. The continuity equation

The linearised continuity equation is given by

(A.11)

then by expressing the divergence in the spheroidal basis as

(A.12)

we get the Eq. (6).

Appendix B: Visual presentation of the terms involved in the derivation of the GLTE

We represent in Fig. B.1 the profiles of the 8 terms (𝒜 is already represented in the top panel of Fig. 5) involved in the derivation of the GLTE (Eq. (34)) in the general case of spheroidal geometry presented in the Table 1 as a function of the pseudo-radius ζ at different reduced colatitudes x. In the main text (Sect. 5.2), we mainly focused on the coefficient 𝒜 and the Brunt-Vaïsälä frequency N because they are the important quantities to be able to build the generalised TAR formalism. The oscillating behaviour of some terms is caused by the fact that the mapping (Eq. (1)) is different in each subdomain of the model which is continuous and derivable, but not twice differentiable with respect to ζ between subdomains.

Fig. B.1. Profiles of the coefficients ℬ, 𝒞, 𝒟, ℰ, ∂xℬ, ∂x𝒞, ∂x𝒟 and involved in the derivation of the GLTE (Eq. (34)) as a function of ζ at different reduced colatitudes x using an ESTER model (M = 3 M⊙,  Xc = 0.7) rotating at 20% of the Keplerian break-up rotation rate for a spin parameter ν = 10.

All Tables

All Figures

	Fig. 1. Spheroidal coordinate system with n = 13 subdomains used to compute the equilibrium model of a star rotating at 60% of its break-up velocity corresponding to the case of the 3 M⊙ ESTER model (with the central fraction in Hydrogen Xc = 0.7) studied in Sect. 5. ζ and θ are the pseudo-radius and the colatitude, respectively. ηi and ΔRi(θ) are defined in our mapping (Eq. (1)).
In the text

	Fig. 2. Squared Brunt–Väisälä frequency N2(ζ, θ) profile for Ω/ΩK = 60%.
In the text

	Fig. 3. Coefficient 𝒜(ζ, θ) profile for Ω/ΩK = 60%.
In the text

	Fig. 4. Profile of the squared Brunt–Väisälä frequency N2 and the relative error of the adopted approximation (Eq. (42)) as a function of ζ at different colatitudes θf using an ESTER model rotating at 60% of the Keplerian breakup rotation rate (The light orange area indicates the margin of error which we allow and the red hatched area represents the convective region of the star).
In the text

	Fig. 5. Profile of the coefficient 𝒜 and the relative error δ𝒜θf of the adopted approximation (Eq. (43)) as a function of ζ at different colatitudes θf using an ESTER model rotating at 60% of the Keplerian break-up rotation rate (The light orange area indicates the margin of error which we allow).
In the text

	Fig. 6. Validity domain of the TAR within the framework of ESTER models (3 M⊙, Xc = 0.7) as a function of the rotation rate Ω/ΩK and the pseudo-radius ζ (with 90% degree of confidence).
In the text

	Fig. 7. Spectrum of the GLTE as a function of the spin parameter ν at ζ = ζlimit and m = 1 for Ω = 0.2ΩK (above) and Ω = 0.9ΩK (below). Blue (respectively, orange) dots correspond to even (respectively, odd) eigenfunctions.
In the text

	Fig. 8. Spectrum of the GLTE as a function of the pseudo-radius for m = 1 and ν = 10 at Ω/ΩK = 0.2 (above) and Ω/ΩK = 0.9 (below). Blue (respectively, orange) dots correspond to even (respectively, odd) Hough functions.
In the text

	Fig. 9. Generalised Hough functions (normalised to unity) as a function of the horizontal coordinate x at different pseudo-radii from ζ = 0.153 (blue) to ζ = ζlimit (yellow) for ν = 10, Ω/ΩK = 0.2 (above), and Ω/ΩK = 0.9 (below). The solid lines correspond to gravity-like solutions with {k = 1, m = 1}, while the dotted lines correspond to Rossby-like solutions with {k = −2, m = −1}.
In the text

Fig. 10. Period spacing pattern in the inertial frame for {k = 0, m = 1} (above) and {k = −2, m = −1} (below) modes at Ω/ΩK = 20% within the standard (red line) and the generalised (blue line) TAR frameworks. The bottom panel for each mode shows the differences with respect to the standard TAR period spacing pattern (the sudden fluctuations in the period spacing pattern are caused by the numerical noise which is introduced by the numerical derivatives of the mapping with respect to ζ and θ used in the resolution of the GLTE).

In the text

	Fig. 11. Detectability of {k = 0, m = 1} (above) and {k = −2, m = −1} (below) modes as a function of the radial order n at a rotation rate Ω/ΩK = 20% based on the frequency resolution of TESS (brown line) and Kepler (indigo line).
In the text

	Fig. 12. N/2Ω term as a function of the pseudo-radius ζ at different rotation rates Ω/ΩK (the hatched area represents the convective region of the star).
In the text

	Fig. 13. N/ω term as a function of the pseudo-radius ζ and the radial order n for {k = 0, m = 1} (top panel) and {k = −2, m = −1} (bottom panel) modes at Ω/ΩK = 20%.
In the text

	Fig. B.1. Profiles of the coefficients ℬ, 𝒞, 𝒟, ℰ, ∂xℬ, ∂x𝒞, ∂x𝒟 and involved in the derivation of the GLTE (Eq. (34)) as a function of ζ at different reduced colatitudes x using an ESTER model (M = 3 M⊙,  Xc = 0.7) rotating at 20% of the Keplerian break-up rotation rate for a spin parameter ν = 10.
In the text