Issue 
A&A
Volume 637, May 2020



Article Number  A65  
Number of page(s)  9  
Section  The Sun and the Heliosphere  
DOI  https://doi.org/10.1051/00046361/201936251  
Published online  15 May 2020 
Rossby modes in slowly rotating stars: depth dependence in distorted polytropes with uniform rotation
^{1}
MaxPlanckInstitut für Sonnensystemforschung, JustusvonLiebigWeg 3, 37077 Göttingen, Germany
email: damiani@mps.mpg.de
^{2}
Institut für Astrophysik, GeorgAugustUniversität Göttingen, FriedrichHundPlatz 1, 37077 Göttingen, Germany
^{3}
Center for Space Science, NYUAD Institute, New York University Abu Dhabi, PO Box 129188, Abu Dhabi, UAE
Received:
5
July
2019
Accepted:
5
March
2020
Context. Largescale Rossby waves have recently been discovered based on measurements of horizontal surface and nearsurface solar flows.
Aims. We are interested in understanding why it is only equatorial modes that are observed and in modelling the radial structure of the observed modes. To this aim, we have characterised the radial eigenfunctions of r modes for slowly rotating polytropes in uniform rotation.
Methods. We followed Provost et al. (1981, A&A, 94, 126) and considered a linear perturbation theory to describe quasitoroidal stellar adiabatic oscillations in the inviscid case. We used perturbation theory to write the solutions to the fourth order in the rotational frequency of the star. We numerically solved the eigenvalue problem, concentrating on the type of behaviour exhibited where the stratification is nearly adiabatic.
Results. We find that for freesurface boundary conditions on a spheroid of nonvanishing surface density, r modes can only exist for ℓ = m spherical harmonics in the inviscid case and we compute their depth dependence and frequencies to leading order. For quasiadiabatic stratification, the sectoral modes with no radial nodes are the only modes which are almost toroidal and the depth dependence of the corresponding horizontal motion scales as r^{m}. For all r modes, except the zero radial order sectoral ones, nonadiabatic stratification plays a crucial role in the radial force balance.
Conclusions. The lack of quasitoroidal solutions when stratification is close to neutral, except for the sectoral modes without nodes in radius, follows from the need for both horizontal and radial force balance. In the absence of super or subadiabatic stratification and viscosity, both the horizontal and radial parts of the force balance independently determine the pressure perturbation. The only quasitoroidal cases in which these constraints on the pressure perturbation are consistent are the special cases where ℓ = m and the horizontal displacement scales with r^{m}.
Key words: Sun: oscillations / methods: analytical / stars: oscillations / stars: rotation / stars: interiors
© C. Damiani et al. 2020
Open 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.
Open Access funding provided by Max Planck Society.
1. Introduction
Rossby waves, which are largescale waves of radial vorticity with retrograde phase speed, have recently been discovered based on measurements of horizontal surface and nearsurface solar flows (Löptien et al. 2018, and confirmed by Liang et al. 2019). These observed waves have frequencies near those of sectoral traditional Rossby waves in a uniformly rotating fluid system (e.g. LonguetHiggins 1964), corresponding to sectoral spherical harmonics of azimuthal order 3 ≤ m ≤ 15. Löptien et al. found that the amplitudes of these Rossby waves do not depend strongly on depth down to 21 Mm below the photosphere but they could not further characterise the radial dependence of the eigenfunctions. Assuming the motion is incompressible, they argued that viscous damping is the reason why they observe only sectoral Rossby modes in the Sun. Here, we are interested in understanding, in the case of a more realistic stellar stratification, why it is only these modes that are observed and we aim to model the radial structure of the observed modes. We limit our attention to the retrogradepropagating Rossby waves discussed in Löptien et al., which have a dispersion relation close to that of traditional (nonmagnetic) Rossby waves, so we do not discuss magnetic Rossby waves (e.g. Zaqarashvili et al. 2010; McIntosh et al. 2017; Dikpati et al. 2018).
The restoring force for traditional Rossby waves is the Coriolis force. Rossby waves have been studied extensively in the geophysical context (see e.g. the textbook by Vallis 2006) with a special focus on their horizontal motion and applications to the Earth’s atmosphere and oceans (Rossby 1939; Dickinson 1978), as well as to the atmospheres of Jupiter and Venus (e.g. Allison 1990; Covey & Schubert 1982; Nara et al. 2019).
In the stellar context, waves analogous to planetary Rossby waves are known as r modes (or quasitoroidal modes, see e.g. Papaloizou & Pringle 1978; Unno et al. 1989). They have been considered in the photosphere of the Sun (starting with the speculative work of Plaskett 1966), as well as near the base of the convection zone (see the series of papers starting with Gilman 1969) where the stratification is assumed to be subadiabatic. Wolff & Blizard (1986) studied the properties of the r modes in the convective zone of the Sun but did not predict any restriction on the existence of nonsectoral modes for uniform rotation.
In general, it can be shown that for lowfrequency nonradial oscillations of a rotating star, the spheroidal components associated to the spherical harmonics of degree ℓ couple with the toroidal components of adjacent degrees ℓ ± 1. Furthermore, these toroidal components ℓ ± 1 themselves couple with the spheroidal components associated with ℓ and ℓ ± 2. Thus, without significant simplifications, a nonradial oscillation mode in a rotating star is given by an infinite sum of terms proportional to spherical harmonics with different degrees ℓ for a given azimuthal index m (Zahn 1966; Berthomieu et al. 1978). In a numerical analysis, a truncation of the series is inevitable. Previous works have opted for drastic truncation, retaining only the first two terms (see e.g. Lee & Saio 1987), which may affect the results significantly.
Two approaches have been considered in order to remove this difficulty. One of these is the socalled traditional approximation, in which the horizontal component of the rotation vector is neglected (Lee & Saio 1989, 1997). Then the Coriolis force associated with radial motion and the radial component of the Coriolis force associated with horizontal motion are both neglected. Alternatively, the solution can be sought using asymptotic expansion relative to a small parameter proportional to the rotation frequency (Provost et al. 1981; Smeyers et al. 1981). The former approximation is valid locally in regions of the star where both the rotation frequency and the pulsation angular frequency in the corotating frame are significantly smaller than the BruntVäisälä frequency (Lee & Saio 1989; Townsend 2003), whereas the latter only requires slow stellar rotation to be valid.
In this paper, we are interested in the lowdegree modes observed in a slowly rotating star, so we will consider the framework of Provost et al. (1981, hereafter P81), concentrating on the cases where the stratification is close to, but not exactly, adiabatic. Our results are not merely an addition to those of P81, but they provide better insight into the nature of r modes, and arguably amend those of P81. Assuming that the stellar interior is inviscid and the motions are adiabatic, we show that the ℓ = m mode with no radial nodes is the only almost toroidal Rossby mode which can be present for uniform rotation. The corresponding eigenfunctions scale as r^{m}.
2. Digest of the analysis of Provost et al. (1981)
This paper follows the formalism developed in P81, with minor changes in the notation. In this section, we summarise the main points of their method that we find to be significant for understanding our results. P81 considers a rotating star, in a corotating reference frame of basis vectors (), with the origin at the star’s centre of mass and spherical coordinates (r, ϑ, φ), where r is the radial distance to the origin, ϑ is the polar angle, and φ is the azimuthal angle. Rotation is assumed to be uniform, with angular frequency parallel to . The densities are in units of , the pressures are in units of , the times are in units of , and the lengths are in units of R_{⋆}, where R_{⋆} and M_{⋆} are the radius and mass of the star and G the universal gravitational constant. The dimensionless angular frequency of the star is then denoted as:
For the Sun^{1}, ε = Ω_{⊙}τ_{⊙} = 4.5 × 10^{−3}. In a uniformly rotating star, the isobaric and isopsynic surfaces coincide with the level surfaces of constant total potential (gravitational and centrifugal). This is known as the PoincaréWavre theorem and holds whatever the equation of state of the gas. For slow rotators (ε ≪ 1), those surfaces can be expressed through a distortion term of the order of ε^{2} and a function α determined by the internal structure. These surface levels can be used to implicitly define a set of curvilinear coordinates (x, θ, ϕ),
where the new coordinate x is constant on surfaces of constant density and pressure. The surface of the star is an isobaric surface x = constant and we select a given normalisation such that x = 1 at the surface. With respect to the new variable x, the equilibrium pressure p(x) and density ρ(x) are then independent of θ and ϕ. The oscillations are treated as a small perturbation around this static equilibrium state. To study linear modes of oscillations, the temporal and longitudinal structure of all perturbed quantities and the displacement ξ are assumed to be proportional to^{2},
The equations governing the small amplitude, periodic, adiabatic oscillations of a uniformly rotating star are obtained by writing the linearised equations for the conservation of angular momentum, mass, and energy:
where p′ and ρ′ are the Eulerian perturbations of pressure and density and is the first adiabatic exponent. Here, we have neglected the perturbation of the gravitational potential (Cowling’s approximation) and have assumed the flows are adiabatic and the viscosity is negligible. Cowling’s approximation was shown to be justified in most cases by P81.
Saio (1982) argues that the adiabatic assumption is justified in stellar radiative zones except near the boundaries. The assumption that the flow is inviscid is also likely to be justifiable in radiative zones; however, in stellar convection zones, nonadiabatic mixing by turbulent convective motions and turbulent viscosity are likely to be relevant on timescales shorter than or comparable to the rotation period of the star. Notwithstanding, we chose to follow P81 in considering the inviscid and adiabatic case because it reveals interesting elements of basic physics. In addition, we expect that modes which rely for their existence on viscosity or nonadiabatic processes will decay faster than modes which exist in the adiabatic, inviscid case. Equations (6) through (8) form a system of partial differential equations that requires appropriate boundary conditions to constitute a wellposed boundary value problem. We consider the boundary conditions at the centre and the surface of the rotating star. At the centre, the displacement must remain finite. The conservation of momentum across the nonspherical surface of the star requires that the Lagrangian pressure perturbation δp = 0 at the surface.
2.1. Series expansion in terms of ε
P81 considers solutions for ξ = (ξ_{r}, ξ_{θ}, ξ_{ϕ}),p′,ρ′,σ, that are solutions for the system represented in Eqs. (6) through (8) for cases where the rotation rate is small. This motivates an expansion of the form:
Our methodology is identical to that of P81 but our notation is slightly different. P81 use the indices j = 0, 1, … only for nonnull terms in the final expansion: for example, in their notation (their Eq. (2)). In our notation, we index all terms. Thus, we write σ = σ_{0} + εσ_{1} + ε^{2}σ_{2} + ε^{3}σ_{3} + …. Our notation is closer to that of Smeyers et al. (1981). As is implicit in the notation of P81, symmetry arguments lead to the conclusion that the coefficients of the even powers of ε in the expansion of σ must be 0, so this difference with respect to P81 is purely one of notation. Similarly, an inspection of Eqs. (6) through (8) reveals that the even and odd coefficients in the expansion of ξ_{θ}, ξ_{ϕ}, ξ_{r}, p′, and ρ′ decouple. The solutions for the odd terms in the expansion of these quantities are trivial if we know the solution for the even powers (they correspond to the solutions for the expansion keeping only the even terms where all quantities except σ are multiplied by ε).
P81 deals with with the generalisation of modes that are purely toroidal for nonrotating stars. The displacement vector of these toroidal modes satisfies ∇ ⋅ ξ = 0 and ξ_{r} = 0. They are also characterised by σ^{2} = 0. In the case where rotation is present, these toroidal modes become nontrivial and develop characteristics similar to Rossby waves in the Earth’s atmosphere and oceans. They are often referred to as r modes, after the seminal work of Papaloizou & Pringle (1978), or quasitoroidal mode according to the nomenclature of P81. They are called quasitoroidal because, to zerothorder, they have the same properties as the toroidal modes (). To summarise, the quasitoroidal modes are nonradial modes of lowfrequency, whose radial displacement is small compared to their horizontal motion.
P81 performs the equivalent of a fourthorder expansion in terms of the small parameter ε in Eqs. (6) through (8). At the zeroth order, the system (6)–(8) only retains zerothorder quantities and is reduced to:
, and
which are the nonrotating toroidal modes. At higherorder, we have to consider the value of the Ledoux discriminant,
which is a measure of convective instability as will be discussed further in Sect. 3. The case A = 0 will be discussed later. If A ≠ 0, the secondorder approximation of the conservation of momentum ( components) yields:
where
is the unperturbed gravity. The radial part of the conservation of momentum Eq. (12) involves the component of the displacement that is normal to isopotential surfaces, . Equation (12) is P81’s (7a). Combining our Eqs. (13) and (14) gives their Eq. (6a), which only involves zerothorder quantities.
Equations (13), (14), and (10) express the conservation of total vertical angular momentum and of mass to zero order. By elimination of p^{′(2)}, they form a classical Legendre equation. This allows us to find exact solutions for σ_{1} and the angular dependence of and in terms of ℓ, m, along with the associated Legendre polynomials. The radial dependence of the eigenfunctions is not determined at this order. The solutions are Rossby waves on surfaces of constant x, as in Eq. (9a) in P81:
Here, we are interested in retrograde Rossby waves, so we are only going to consider m > 0. At the fourthorder level of the approximation, the θ and ϕ components of the conservation of momentum can be combined to eliminate the pressure perturbation terms p^{′(4)} and p^{′(2)}, yielding Eq. (7c) from P81. Finally, the system of equations is closed by taking the secondorder expansion of the continuity equation, yielding Provost et al.’s Eq. (7d).
After some manipulation, the aforementioned closed system of equations can be reduced to a single differential equation for the amplitude of the horizontal displacement, Eq. (11) of P81,
where λ_{1}, λ_{2}, λ_{3}, and λ_{4} are functions of m and ℓ and are defined in P81. When A ≠ 0, for x ∈ (0, 1) and with appropriate boundary conditions, Eq. (19) is a Sturm–Liouville eigenvalue problem.
2.2. Boundary condition at the centre
Near the centre x ∼ 0, the quantities that appear in Eq. (19) are in order of magnitudes as follows:
where the subscript c denotes values at the centre. Therefore, through an asymptotic analysis of Eq. (19) in the vicinity of x = 0, we have:
Using the Frobenius method (Bender & Orszag 1999), we look for a solution of the form in the vicinity of x = 0. According to (25), the coefficient of the lowest power of x must satisfy:
which has solutions β = ℓ and β = −ℓ−1. The only solution that ensures the regularity of is in the vicinity of x = 0. Hence we use the following boundary condition,
and
2.3. Boundary condition at the surface
At the free surface the Lagrangian pressure perturbation,
must vanish. This is a very stringent condition for the existence of nonsectoral Rossby waves in the inviscid case, and its implications depend on whether ρ/A vanishes (or not) at the surface.
2.3.1. Boundary conditions for a star where the surface has ρ/A = 0
P81 considers a complete polytrope^{3} as a model of a star. The stellar surface for such models is defined by ρ = p = 0.
In general, for a complete polytrope or otherwise, when ρ/A vanishes at the surface, δp = 0 is met at the surface as long as and its derivative are regular. This is the boundary condition that was used in P81.
For the particular case of a polytropic stellar model, characterised by a polytropic index n, the behaviour of the quantities that appear in (19) can be approximated by Taylor expansion of , the solution to the LaneEmden equation, near x = 1 as follows:
Through an asymptotic analysis of Eq. (19) in the vicinity of x = 1, we have:
where
We note that λ_{1}, λ_{2} and λ_{3} are simple, but tediously long, rational functions of ℓ and m. They are given explicitly in P81.
Thus, the condition of regularity for at the surface is:
Hereafter, the notation in bracket with subscript S indicates that the value is taken at the surface.
2.3.2. Boundary conditions when ρ/A ≠ 0 at the surface
Truncated polytropes have been used as models of stars that feature an atmosphere (e.g. Hendry 1993; Bogdan & Cally 1995). In these models, the polytrope is truncated at some location, with nonzero pressure and density, that represents the stellar surface. In the limit that the density of the atmosphere is small compared to the surface density, the boundary condition corresponds to a free surface (with vanishing Lagrangian pressure perturbation at the surface). Such truncated polytropes with a free surface have previously been used in the study of helioseismic acoustic waves (e.g. Bogdan & Cally 1995). The analysis for Rossby waves proceeds differently according to whether ℓ = m or not.
Cases where ℓ = m. When ℓ = m, the freesurface boundary yields one condition:
In the cases where {ρ/A}_{S} ≠ 0, the boundary condition is reduced to:
Cases where ℓ≠ m. When ℓ ≠ m, the freesurface boundary condition involve two associated Legendre polynomials and , whose coefficients must both vanish and, thus, this yields two conditions^{4} :
In general, ρ/A does not vanish at the surface of stars, Eqs. (40) and (41) can only be met simultaneously by requiring both and at the surface. These two boundary conditions, together with the boundary condition at the centre in Eq. (27), fully specify the problem and the only solution is then the trivial solution, .
Consequently, for nonvanishing density at the surface, there is no nontrivial solution when ℓ ≠ m. The sectoral mode is the only quasitoroidal mode that can satisfy the freesurface boundary condition in the inviscid case. Let us note that the nonsectoral modes found by Wolff & Blizard (1986), using P81’s derivations applied to the Sun, are obtained by imposing at x = 0.999, which is inconsistent with a freesurface boundary condition. We remark as well that the nopenetration boundary condition is with:
so that Wolff & Blizard (1986) also appears to be inconsistent with the nopenetration boundary condition.
2.4. Comparison of boundary conditions of polytropes and truncated polytropes
It is interesting to compare the boundary conditions for the complete and truncated polytropes. For ℓ = m these are Eqs. (37) and (39), respectively. The quantity Q(1) is given by Eq. (36) evaluated at the surface and depends on m. For the special cases where ℓ = m, γ = 5/3 and n = 3/2, we find that Q(1)/(n + 1) = m. The values of Q/(n + 1) evaluated at the surface of the full polytrope (x = 1) are given in Table 1 for some of the other cases studied in this paper with ℓ = m = 3. In Table 1, we see that Q/(n + 1) is approximately equal to m near n = 1.5 (it is exactly m at n = 1.5). A consequence of this is that in the neighbourhood of n = 1.5, the boundary condition for the truncated polytrope is equivalent to that for the complete polytrope. This indicates that the results have some robustness to the details of the model.
Value of Q/(n + 1) evaluated at x = 1 for ℓ = m = 3 for different polytropic indices.
For ℓ ≠ m the situation is more complicated and no solutions exist for the truncated polytrope. This is also the case for the complete polytrope when n = 1.5, (see Appendix I in P81). For other values of n, the solution for the complete polytrope of course exists, and is as described in P81.
3. Results
To solve the boundaryvalue problem for Eq. (19), we used an implementation of a fourthorder collocation algorithm based on control of residuals, provided by the function solve_bvp from the scipy.integrate python module. We follow P81 in considering polytropes, characterised by a polytropic index n. In the stellar context, a polytrope is a gas spheroid in gravitational equilibrium, where the pressure is related to the density by the relation,
where K is a constant of proportionality. By definition of the free surface, a complete polytrope has p = 0 and ρ = 0 at the surface. For all polytropes, the Ledoux discriminant A is a monotonic function of x and it is everywhere positive when n < 1.5 (superadiabatic stratification). The polytrope with index n = 1.5 has the interesting property that A = 0 everywhere in the star (adiabatic stratification). Finally, for n > 1.5, A is everywhere negative (subadiabatic stratification. The Ledoux discriminant is the argument of the criterion for convection which develops or not whether A > 0 or A < 0 (Ledoux & Walraven 1958). So we also refer to the convectively unstable case when n < 1.5, neutrally stratified case when n = 1.5 and the radiative case when n > 1.5.
The shape function α(x) in Eq. (2) can be derived for the distorted polytrope (Chandrasekhar 1933) and is obtained by solving the LaneEmden equation modified for rotation. We obtain the full polytrope by setting the surface x = 1 at the first zero of the LaneEmden function and the truncated polytrope by setting the surface at x = 0.999.
To test our solver, we computed the eigenvalues for n = 1 and n = 3 given in Tables 1 and 2 of P81 to the same decimal place accuracy, for a complete polytrope using the boundary conditions in Eqs. (27), (28), and (37). We also reproduce the radial dependence of for the same (ℓ,m, k) as in their Figs. 1–3.
We also computed the eigensolutions for the truncated polytropes when ℓ = m, using the boundary conditions in Eqs. (27), (28), and (39). As is to be expected from Table 1, for sectoral modes, the eigensolutions for the full polytrope and the same polytrope truncated very close to the surface are nearly identical close to n = 1.5, and the differences even for n = 1 and n = 3 are small and not distinguishable in the figures that follow. Hence, we only show solutions for a polytrope truncated at x = 0.999.
For cases when ℓ ≠ m, we refer to P81 for the solutions in a complete polytrope. Again we note that there are no solutions satisfying the freesurface boundary conditions for truncated polytropes.
3.1. Nonadiabatic stratification is essential to the normal force balance for all sectoral rmodes except those with no radial nodes
To get a better understanding of the effect of stratification on the r modes, we also solve the problem for 1 ≤ n ≤ 3. As an example, Fig. 1 shows the eigenfrequencies found for ℓ = 3 and m = 3 as a function of the polytropic index n for several values of the radial order k. We see that for all modes, with the notable exception of the sectoral mode with k = 0, the eigenfrequencies become increasingly large as n gets close to 1.5 (where A = 0). Let us stress here that the derivation of Eq. (19) is obtained by a singular perturbation method. Thus the asymptotic expansion (9) does not necessarily converge, and the solution is valid when σ_{3}≪σ_{1}/ε^{2}, which means, here, ε^{3}σ_{3}/(2πτ_{⊙})≪10^{5} pHz for ℓ = m = 3 (see Fig. 1).
Fig. 1. Thirdorder term in the frequency expansion σ_{3}ε^{3} as a function of the polytropic index n, for the sectoral r modes ℓ = m = 3 and several radial orders k. The corresponding frequency term in the incompressible case, which has only ℓ = m, k = 0 solutions, is given by a blue dotted line (see text). The frequency is given in pHz, using the relevant solar quantities, and is displayed in symmetric log scaling, with linear scaling between ±6 pHz. The quantity τ_{⊙} is our unit of time and corresponds to the Sun’s dynamical time scale. 
Solving Eq. (19) allows the derivation of all the perturbed quantities of our problem. Figure 2 shows the radial force balance associated to r modes given by Eq. (12) for the polytrope n = 1.49 and for ℓ = m = 3. In the momentum equation, the Coriolis term is balanced by a nontrivial combination of and . For the modes with k ≠ 0 (Fig. 2, left) the nonadiabatic stratification plays an essential role in the radial force balance. This turns out to be true for all the modes (i.e. also those with ℓ ≠ m) except the ℓ = m, k = 0 modes (discussed in Sect. 4). For the case with ℓ = m and k = 0, the radial force balance essentially lies between the Coriolis term and (Fig. 2, right) – the term involving the nonadiabatic stratification essentially plays no role. This is a special property of the ℓ = m, k = 0, r modes.
Fig. 2. Terms involved in the normal force balance at θ = 0.2, according to Eq. (12) for the slightly subadiabatic polytropic index n = 1.49, in the case of the (ℓ,m, k) = (3, 3, 1) mode (left) and the (3,3,0) sectoral mode of zero radial order (right). All the quantities are dimensionless but they have been normalised consistently with Fig. 4. 
3.2. The depth dependence of the sectoral modes of zero radial order is x^{m} for quasiadiabatic stratification
Figure 3 shows the radial structure for (ℓ,m, k) = (3, 3, 0) for different values of the polytropic index n. As shown in P81, when A = 0, the only nontrivial modes are the sectoral modes with ℓ = m. The solution must then have the form , which has no radial nodes. However, in this case there is no finite value of σ_{3} that can give a finite radial displacement at the surface of the complete polytrope, since in this case the solution is divergent at the surface (P81). This problem does not exist for truncated polytropes.
Fig. 3. Solutions to Eq. (19) for the sectoral modes ℓ = m = 3 of zero radial order for different polytropic indices. The functions have been normalised in order to have in all cases. 
In the limit of n → 1.5, as we can see in Fig. 3, we find that the depth dependence of the sectoral mode with k = 0 is very close to x^{m}. Since this mode is the only one allowed to exist in the case A = 0, we find that there is no discontinuity of solutions for this mode and the solutions slowly depart from an x^{m} dependence as the stratification departs from neutral.
The solution near n = 1.5 (where A = 0) has a radial dependence proportional to x^{m}. This is also the form of the incompressible Rossby wave (Bryan 1889; Provost et al. 1981), which has no dependence on the stratification, as will be discussed in Section 4.
3.3. Symmetries of the eigenfunctions about n = 1.5
Figures 4 and 5 show eigenfunctions for the (3, 3, 1) and (3, 3, 0) modes, both for n = 1.49 and n = 1.51. These values of n were chosen because the background is convectively unstable in one case and stable in the other. Here, instead of the radial component of the displacement, we plot the more relevant component of the displacement that is normal to isopotential surfaces.
Fig. 4. Meridional cuts (x, θ) showing the leading order terms, for solar rotation, of the latitudinal σ_{1}ξ_{θ} and azimuthal σ_{1}ξ_{ϕ} flow velocity, the relative normal displacement ξ_{x}/H_{p}, the relative pressure p′/p, and density ρ′/ρ perturbations (from left to right respectively, as labelled in the colour bar) for the (ℓ,m, k) = (3, 3, 1) mode and for n = 1.49 (top) and n = 1.51 (bottom). H_{p} is the pressure height scale. The radial dependence of all the quantities has been normalised consistently so that at phase 0, a dimensionless velocity of 1 on this scale corresponds to a velocity of 1 m s^{−1} on the Sun. 
Fig. 5. Same as Fig. 4 but for the (ℓ,m, k) = (3, 3, 0) mode. From left to right: latitudinal and azimuthal flow velocity, relative normal displacement, relative pressure, and density perturbations for n = 1.49 (top) and n = 1.51 (bottom). 
The solution for the case with k ≠ 0 (Fig. 4) demonstrates that most of the physical quantities vary only slightly across n = 1.5 except for ξ_{x}, which flips the sign. The reason for this flip in the sign can be inferred from Eq. (12), which is further discussed in Sect. 4. Conversely, in the case ℓ = m and k = 0 (Fig. 5), all the quantities vary only slightly across n = 1.5 and, in particular, the sign for ξ_{x} does not change.
This difference results from the fact that the radial force balance in the cases with ℓ ≠ m or k ≠ 0 the super– or sub–adiabaticity plays an essential role and ξ_{x} changes sign with A. By way of contrast, the sectoral mode of zero radial order does not depend essentially on the adiabaticity and ξ_{x} varies smoothly with n.
4. Discussion
The solutions with ℓ = m and k = 0 are qualitatively different from all other solutions. We show the mode with (ℓ,m, k) = (3, 3, 0) as an example. In the neighbourhood of n = 1.5 (corresponding to A ∼ 0), all values of ℓ, m, and k admit solutions in a complete distorted polytrope, as discussed above. Only the sectoral modes (ℓ = m) would be admissible if the density at the surface did not vanish. For the cases with ℓ = m and k = 0, the term remains small and the solution is not sensitive to A; the solution does not depend on nonadiabatic stratification. All other solutions are baroclinic modes and the term is approximately constant in the neighbourhood of n = 1.5, with ξ_{r} being arbitrarily large in the limit n → 1.5 (A → 0) and flipping sign at n = 1.5. Similar reasoning explains the behaviour shown in Fig. 1 for the eigenfrequencies.
The radial displacements associated with solutions where the nonadiabatic stratification is essential (i.e. where the contribution to the radial force balance from Agξ_{r} is substantial) tend to infinity as n approaches 1.5 (where A → 0). These solutions cease to be quasitoroidal. Hence the only solutions which are valid^{5} near n = 1.5 are those with ℓ = m and k = 0.
The lack of quasitoroidal solutions as n approaches 1.5, except when ℓ = m and k = 0, is not a consequence of the chosen expansion. It follows from the statement that the system needs to be in both horizontal and radial force balance. In the absence of super or subadiabatic stratification and viscosity, both the horizontal and radial force balances independently determine the pressure perturbation. The only case in which the two determinations of the pressure perturbation are consistent and quasitoroidal are the special cases where ℓ = m, k = 0 and the horizontal displacement scales with x^{m}. It is here that we make contact with the example Rossby waves in an incompressible, unstratified spherical shell discussed by Löptien et al. (2018).
The existence of nontrivial solutions only in the ℓ = m case with an r^{m} radial dependence is a finding also described in P81 in reference to the incompressible case with arbitrary stratification. This is understood by considering that the horizontal force balance sets the horizontal structure of the pressure perturbation independently of A or incompressibilitity. In the incompressible case (or in the case with A = 0) and for quasitoroidal motion, the pressure perturbation must alone balance the radial component of the Coriolis force. This is only possible in the case, ℓ = m, and it results in an r^{m} dependence, with no dependence on stratification.
The expression of the thirdorder term σ_{3} of the frequency expansion for the ℓ = m mode in the incompressible case for nonspherical shapes have been obtained by Bryan (1889) and P81,
Figure 1 shows these values of σ_{3} for m = 3 as a function of n as a blue dotted line^{6}; they are the same order of magnitude as the eigenfrequency associated to the sectoral mode of zero radial order. This is to be expected since here as we are seeking solutions that are quasitoroidal, which means that the flow is divergence free to zerothorder.
All these solutions are approximately the same, although the problem is formulated differently in each case (compressible vs incompressible, polytrope versus neutrally stable stratification). In the incompressible case or the neutrally stratified case, the radial force balance is independent of the normal displacement. In this case, the problem admits only solutions that have ℓ = m and k = 0 and an x^{m} dependence. For the quasitoroidal modes, the solutions are again approximately the same because in this case also the toroidal components of the motion are the same (spherical harmonics characterised by ℓ, m). This, in turn, determines the pressure perturbation, which then must balance the radial component of the Coriolis force if the motions are to remain quasitoroidal. This only happens for ℓ = m and k = 0. For ℓ ≠ m or k ≠ 0, the modes develop substantial radial velocity whenever the stratification is close to, but not strictly, neutral.
We also found the eigenvalues and eigenfunctions for the case where we artificially set α = 0, which then treats the star as a perfect sphere. Consistent with P81, we found the changes to σ_{3} were substantial relative to σ_{3} (i.e. can have the wrong sign), but since σ_{3}ε^{3} is small, this amounts to a change of the order of 10 pHz for the solar rotation rate. More important are the changes in the radial structure of the eigenfunctions, which can be seen by comparing Fig. 6 with the equivalent in Fig. 3. This would suggest that the distortion of the background star due to rotation should be included when determining the radial structure of the eigenfunctions, especially for nonadiabatic stratifications.
Fig. 6. Same as Fig. 3 for the (ℓ,m, k) = (3, 3, 0) mode except with spherical geometry assumed by setting α = 0. The radial structure of the eigenfunctions is affected by omission of the shape distortion when n ≠ 1.5. 
5. Conclusions
In this paper, we used polytropes to understand the fundamental properties of quasitoroidal modes for slowly and uniformly rotating stars.
The sectoral r modes of zero radial order are qualitatively different from the other r modes in that they do not rely on nonadiabatic stratification to balance the radial component of the Coriolis force. This is critical when the stratification is close to neutral (for polytropes, this is in the neighbourhood of n = 1.5). In this neighbourhood, the modes with ℓ ≠ m or k ≠ 0 all involve large radial displacements and are no longer nearly toroidal. The ℓ = m, k = 0 modes retain a small radial displacement through the star (except possibly at the surface if the surface has p = ρ = T = 0 as pointed out by P81). The depth dependence of the horizontal displacement is close to x^{m}, as it is the only solution allowing both horizontal and vertical force balance, in the absence of viscosity and the lack of a buoyant contribution from nonadiabatic stratification.
Consequently, in the case of the Sun, we speculate that only the ℓ = m, k = 0 quasitoroidal modes can exist in the convection zone, which is very close to adiabatically stratified (corresponding approximately to the polytropic index n = 1.49). They are presumably the modes that are observed at the solar surface. In Fig. 7, we show the kinetic energy density associated with those modes on a meridional cut. In Fig. 8, we plot the corresponding radial dependence at the equator. These figures suggest that solar Rossby waves have diagnostic potential because different modes have different radial and latitudinal distribution of the kinetic energy density. It can be noticed that the kinetic energy density of the modes peaks in the interior (not at the surface as p modes would), and that the m = 3 mode has a kinetic energy density that peaks near x = 0.75.
Fig. 7. Meridional cuts (x, θ) showing the normalised kinetic energy density of the sectoral modes of zero radial order, for n = 1.49, and solar rotation for different values of m as labelled. The dashed line indicates the position of the base of the convection zone in the Sun, at x = 0.7. The contour levels are marked by black lines in the colour bar. 
Fig. 8. Radial dependence at the equator of the normalised kinetic energy density of the sectoral modes of zero radial order shown in Fig. 7. The dashed line indicates the position of the base of the convection zone in the Sun, at x = 0.7. 
However, the results of this study cannot be directly applied to the Solar case since the Sun is not a uniformly rotating inviscid polytrope. In particular, latitudinal differential rotation will modify the latitudinal eigenfunctions (Gizon et al., in prep.). Also, the inclusion of radial differential rotation has not been considered here. Furthermore, the Sun has a convectively stable radiative interior beneath its convective envelope that requires a careful consideration of the matching conditions.
We follow the same sign conventions as in Löptien et al. (2018), this means that our frequencies have the opposite sign as those of P81.
Polytropes will be described in more details in Sect. 3.
The same boundary conditions were already given in Smeyers et al. (1981) Eqs. (65) and (66). There is however a sign error in their Eq. (65), which is corrected in our Eq. (41).
In Eq. (44) only the α coefficient is a function of n, because the shape distortion depends on density distribution.
Acknowledgments
We thank the referee for useful comments which helped improve the manuscript. We acknowledge partial support from the German Aerospace Center (Deutsches Zentrum fúr Luft und Raumfahrt) under PLATO Data Center grant 50OO1501. We also acknowledge partial support from ERC Synergy Grant WHOLESUN 810218. This research made use of NumPy (Van Der Walt et al. 2011), SciPy (Jones et al. 2001) and matplotlib, a Python library for publication quality graphics (Hunter 2007).
References
 Allison, M. 1990, Icarus, 83, 282 [NASA ADS] [CrossRef] [Google Scholar]
 Bender, C. M., & Orszag, S. A. 1999, Advanced Mathematical Methods for Scientists and Engineers I: asymptotic Methods and Perturbation Theory (Berlin, Heidelberg: SpringerVerlag) [Google Scholar]
 Berthomieu, G., Gonczi, G., Graff, P., Provost, J., & Rocca, A. 1978, A&A, 70, A597 [NASA ADS] [Google Scholar]
 Bogdan, T. J., & Cally, P. S. 1995, ApJ, 453, 919 [NASA ADS] [CrossRef] [Google Scholar]
 Bryan, G. H. 1889, Philos. Trans. T. Soc. London Series A, 180, 187 [NASA ADS] [CrossRef] [Google Scholar]
 Chandrasekhar, S. 1933, MNRAS, 93, 390 [NASA ADS] [CrossRef] [Google Scholar]
 Covey, C., & Schubert, G. 1982, J. Atmos. Sci., 39, 2397 [NASA ADS] [CrossRef] [Google Scholar]
 Dickinson, R. E. 1978, Ann. Rev. Fluid Mech., 10, 159 [NASA ADS] [CrossRef] [Google Scholar]
 Dikpati, M., Belucz, B., Gilman, P. A., & McIntosh, S. W. 2018, ApJ, 862, 159 [NASA ADS] [CrossRef] [Google Scholar]
 Gilman, P. A. 1969, Sol. Phys., 8, 316 [NASA ADS] [CrossRef] [Google Scholar]
 Hendry, A. W. 1993, Am. J. Phys., 61, 906 [NASA ADS] [CrossRef] [Google Scholar]
 Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
 Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open source scientific tools for Python [Google Scholar]
 Ledoux, P., & Walraven, T. 1958, Handbuch der Physik, 51, 353 [NASA ADS] [CrossRef] [Google Scholar]
 Lee, U., & Saio, H. 1987, MNRAS, 224, 513 [NASA ADS] [CrossRef] [Google Scholar]
 Lee, U., & Saio, H. 1989, MNRAS, 237, 875 [NASA ADS] [Google Scholar]
 Lee, U., & Saio, H. 1997, ApJ, 491, 839 [NASA ADS] [CrossRef] [Google Scholar]
 Liang, Z.C., Gizon, L., Birch, A. C., & Duvall, T. L. 2019, A&A, 626, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 LonguetHiggins, M. S. 1964, Proc. R. Soc. London Ser. A, 279, 446 [NASA ADS] [CrossRef] [Google Scholar]
 Löptien, B., Gizon, L., Birch, A. C., et al. 2018, Nat. Astron., 2, 568 [NASA ADS] [CrossRef] [Google Scholar]
 McIntosh, S. W., Cramer, W. J., Pichardo Marcano, M., & Leamon, R. J. 2017, Nat. Astron., 1, 0086 [CrossRef] [Google Scholar]
 Nara, Y., Imamura, T., Murakami, S., et al. 2019, J. Geophys. Res. (Planets), 124, 1143 [NASA ADS] [Google Scholar]
 Papaloizou, J., & Pringle, J. E. 1978, MNRAS, 182, 423 [NASA ADS] [CrossRef] [Google Scholar]
 Plaskett, H. H. 1966, MNRAS, 131, 407 [NASA ADS] [CrossRef] [Google Scholar]
 Provost, J., Berthomieu, G., & Rocca, A. 1981, A&A, 94, 126 [NASA ADS] [Google Scholar]
 Rossby, C.G. 1939, J. Marine Res., 2, 38 [CrossRef] [Google Scholar]
 Saio, H. 1982, ApJ, 256, 717 [NASA ADS] [CrossRef] [Google Scholar]
 Smeyers, P., Craeynest, D., & Martens, L. 1981, Ap&SS, 78, 483 [NASA ADS] [CrossRef] [Google Scholar]
 Townsend, R. H. D. 2003, MNRAS, 340, 1020 [NASA ADS] [CrossRef] [Google Scholar]
 Unno, W., Osaki, Y., Ando, H., Saio, H., & Shibahashi, H. 1989, Nonradial Oscillations of Stars (Tokyo: University of Tokyo Press) [Google Scholar]
 Vallis, G. K. 2006, Atmospheric and Oceanic Fluid Dynamics: Fundamentals and Largescale Circulation (Cambridge University Press) [CrossRef] [Google Scholar]
 Van Der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Comput. Sci. Eng., 13, 22 [CrossRef] [Google Scholar]
 Wolff, C. L., & Blizard, J. B. 1986, Sol. Phys., 105, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Zahn, J. P. 1966, Ann. Astrophys., 29, 313 [NASA ADS] [Google Scholar]
 Zaqarashvili, T. V., Carbonell, M., Oliver, R., & Ballester, J. L. 2010, ApJ, 709, 749 [NASA ADS] [CrossRef] [Google Scholar]
All Tables
Value of Q/(n + 1) evaluated at x = 1 for ℓ = m = 3 for different polytropic indices.
All Figures
Fig. 1. Thirdorder term in the frequency expansion σ_{3}ε^{3} as a function of the polytropic index n, for the sectoral r modes ℓ = m = 3 and several radial orders k. The corresponding frequency term in the incompressible case, which has only ℓ = m, k = 0 solutions, is given by a blue dotted line (see text). The frequency is given in pHz, using the relevant solar quantities, and is displayed in symmetric log scaling, with linear scaling between ±6 pHz. The quantity τ_{⊙} is our unit of time and corresponds to the Sun’s dynamical time scale. 

In the text 
Fig. 2. Terms involved in the normal force balance at θ = 0.2, according to Eq. (12) for the slightly subadiabatic polytropic index n = 1.49, in the case of the (ℓ,m, k) = (3, 3, 1) mode (left) and the (3,3,0) sectoral mode of zero radial order (right). All the quantities are dimensionless but they have been normalised consistently with Fig. 4. 

In the text 
Fig. 3. Solutions to Eq. (19) for the sectoral modes ℓ = m = 3 of zero radial order for different polytropic indices. The functions have been normalised in order to have in all cases. 

In the text 
Fig. 4. Meridional cuts (x, θ) showing the leading order terms, for solar rotation, of the latitudinal σ_{1}ξ_{θ} and azimuthal σ_{1}ξ_{ϕ} flow velocity, the relative normal displacement ξ_{x}/H_{p}, the relative pressure p′/p, and density ρ′/ρ perturbations (from left to right respectively, as labelled in the colour bar) for the (ℓ,m, k) = (3, 3, 1) mode and for n = 1.49 (top) and n = 1.51 (bottom). H_{p} is the pressure height scale. The radial dependence of all the quantities has been normalised consistently so that at phase 0, a dimensionless velocity of 1 on this scale corresponds to a velocity of 1 m s^{−1} on the Sun. 

In the text 
Fig. 5. Same as Fig. 4 but for the (ℓ,m, k) = (3, 3, 0) mode. From left to right: latitudinal and azimuthal flow velocity, relative normal displacement, relative pressure, and density perturbations for n = 1.49 (top) and n = 1.51 (bottom). 

In the text 
Fig. 6. Same as Fig. 3 for the (ℓ,m, k) = (3, 3, 0) mode except with spherical geometry assumed by setting α = 0. The radial structure of the eigenfunctions is affected by omission of the shape distortion when n ≠ 1.5. 

In the text 
Fig. 7. Meridional cuts (x, θ) showing the normalised kinetic energy density of the sectoral modes of zero radial order, for n = 1.49, and solar rotation for different values of m as labelled. The dashed line indicates the position of the base of the convection zone in the Sun, at x = 0.7. The contour levels are marked by black lines in the colour bar. 

In the text 
Fig. 8. Radial dependence at the equator of the normalised kinetic energy density of the sectoral modes of zero radial order shown in Fig. 7. The dashed line indicates the position of the base of the convection zone in the Sun, at x = 0.7. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext 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 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.