Open Access
Issue
A&A
Volume 646, February 2021
Article Number A64
Number of page(s) 19
Section Stellar structure and evolution
DOI https://doi.org/10.1051/0004-6361/202038654
Published online 08 February 2021

© J. Park et al. 2021

Licence Creative Commons
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.

1. Introduction

Stellar rotation is one of the key physical processes to build a modern picture of stellar evolution (e.g., Maeder 2009, and references therein). Indeed, it triggers the transport of angular momentum and of chemicals, which drives the rotational and chemical evolution of stars, respectively (e.g., Zahn 1992; Maeder & Zahn 1998; Mathis & Zahn 2004). This has a major impact on the late stages of their evolution (e.g., Hirschi et al. 2004), their magnetism (e.g., Brun & Browning 2017), their winds and mass losses (e.g., Ud-Doula et al. 2009; Matt et al. 2015), and the interactions with their planetary and galactic environment (e.g., Gallet et al. 2017; Strugarek et al. 2017).

A robust ab-initio evaluation of the strength of each (magneto-)hydrodynamical mechanism that transports momentum and chemicals is thus mandatory to understand astrophysical observations. In particular, our knowledge of the internal rotation of stars and its evolution has been revolutionized thanks to space-based asteroseismology with the Kepler space mission (NASA; Aerts et al. 2019, and references therein). It has been proven that stars mostly host weak differential rotation in the whole Hertzsprung-Russell diagram, such as our Sun. Stellar interiors are thus the seat of efficient mechanisms that transport angular momentum all along their evolution. These mechanisms have not been identified yet even if several candidates have been proposed, such as stable and unstable magnetic fields (e.g., Moss 1992; Charbonneau & MacGregor 1993; Spruit 1999, 2002; Fuller et al. 2019), stochastically-excited internal gravity waves (e.g., Talon & Charbonnel 2005; Rogers 2015; Pinçon et al. 2017), and mixed gravito-acoustic modes (Belkacem et al. 2015a,b). In this framework, improving our knowledge of the hydrodynamical turbulent transport induced by the instabilities of the stellar differential rotation is mandatory since it has been recently proposed as another potential efficient mechanism to transport angular momentum (Barker et al. 2020; Garaud 2020) while it constitutes one of the cornerstones of the theory of the rotational transport and mixing along the evolution of stars (Zahn 1992).

Since stellar radiation zones are stably stratified, rotating regions, it is expected that the turbulent transport triggered there by the instabilities of the vertical and horizontal shear of the differential rotation should be anisotropic (Zahn 1992). As pointed out in Mathis et al. (2018) and Park et al. (2020), it is the vertical shear instabilities that have received important attention in the literature, in particular by taking into account the impact of the high thermal diffusion in stellar radiation zones (Zahn et al. 1983; Lignières 1999) and the interactions with the horizontal turbulence induced by the horizontal shear (Talon & Zahn 1997). The obtained prescriptions, mainly derived using phenomenological modelings, have been broadly implemented in state-of-the-art evolution models of rotating stars (e.g., Ekström et al. 2012; Marques et al. 2013; Amard et al. 2019). They are now tested using direct numerical simulations devoted to the stellar regime (Prat & Lignières 2013, 2014; Garaud et al. 2017; Prat et al. 2016; Gagnier & Garaud 2018; Kulenthirarajah & Garaud 2018). The horizontal turbulence induced by horizontal gradients of the differential rotation has only been examined in the stellar regime in a few works, which again use phenomenological arguments (Zahn 1992; Maeder 2003), results from lab experiments studying differentially rotating flows (Richard & Zahn 1999; Mathis et al. 2004), and first devoted numerical simulations (Cope et al. 2020). The systematic study of the combined effect of stable stratification, rotation, and thermal diffusion in the stellar context has only been recently undertaken.

In Park et al. (2020), we thus examined the behavior of shear instabilities sustained by a horizontal shear as a function of stratification, rotation, and thermal diffusion. However, in this first work, we neglected the horizontal projection of the rotation vector and the corresponding terms of the Coriolis acceleration, following the traditional approximation of rotation (TAR) which is often used to describe geophysical and astrophysical flows in stably stratified, rotating regions (Eckart 1960). However, recent works have disputed this approximation as it can fail in some cases. For instance, the dynamics of near-inertial waves is strongly influenced by the nontraditional effects in such a way that properties of the wave reflection and associated mixing are largely modified (Gerkema & Shrira 2005; Gerkema et al. 2008). Moreover, couplings between gravito-inertial waves, inertial waves, and wave-induced turbulence in stars can only be properly treated when the full Coriolis acceleration is considered (Mathis et al. 2014). Finally, Zeitlin (2018) demonstrated analytically that the instability of a linear shear flow can be significantly modified if the full Coriolis acceleration is considered. Nontraditional effects can be particularly important for stellar structure configurations where the Coriolis acceleration can compete with the Archimedean force in the direction of both entropy and chemical stratification, for instance during the formation of the radiative core of pre-main-sequence, low-mass stars or in the radiative envelope of rapidly rotating upper main-sequence stars. These regimes should be treated properly to build robust one- or two-dimensional (1D or 2D) secular models of the evolution of rotating stars (e.g., Ekström et al. 2012; Amard et al. 2019; Gagnier et al. 2019).

In this paper, we, therefore, continue our previous work that examined the effect of thermal diffusion on horizontal shear instabilities in stably stratified rotating stellar radiative zones (Park et al. 2020), but we now consider the full Coriolis acceleration. In particular, we investigate how the full Coriolis acceleration modifies the dynamics of two types of shear instabilities: the inflectional and inertial instabilities. In Sect. 2, we formulate linear stability equations derived from the Navier-Stokes equations when the fluid is stably stratified and thermally diffusive in the rotating nontraditional f-plane where both vertical and horizontal components of the rotation vector are taken into account. We consider a base shear flow in a hyperbolic tangent form as a canonical shear profile used in previous studies (Schmid & Henningson 2001; Deloncle et al. 2007; Griffiths 2008; Arobone & Sarkar 2012). In Sect. 3, we provide general numerical results on the inflectional and inertial instabilities. In Sect. 4, we use the Wentzel-Kramers-Brillouing-Jeffreys (WKBJ) approximation in the asymptotic limit of large vertical wavenumbers to investigate how the inertial instability is modified by nontraditional effects in the asymptotic nondiffusive and high-diffusivity cases. In Sect. 5, we analyze in detail numerical results in wide ranges of parameters for both the inflectional and inertial instabilities. In Sects. 6 and 7, we propose expressions for the horizontal turbulent viscosity and the characteristic time of the turbulent transport, which can be possibly applied to the 1D- and 2D-secular evolution models of rotating stars. Finally, in Sect. 8, we provide our conclusions and discussions about the nontraditional effects on the horizontal shear instabilities, and we propose the perspectives of this work for astrophysical flows.

2. Problem formulation

2.1. Navier-Stokes equations and base steady state

We consider the Navier-Stokes and heat transport equations under the Boussinesq approximation in a rotating frame. We define the local Cartesian coordinates (x, y, z) with x the longitudinal coordinate, y the latitudinal coordinate, and z the vertical coordinate:

(1)

(2)

(3)

where u = (u,v,w) is the velocity, p is the pressure, Θ is the temperature deviation from the reference temperature T0, f = (0,fh, 0,fv, 0) is the Coriolis vector with the horizontal (latitudinal) Coriolis component fh, 0 = 2Ω0 sin θ and the vertical Coriolis component fv, 0 = 2Ω0 cos θ where is the stellar rotation rate, θ is the colatitude, ρ0 is the reference density, g = (0, 0, −g) is the gravity, ν0 is the reference viscosity, κ0 is the reference thermal diffusivity, αT is the thermal expansion coefficient that assumes a linear relation between the density deviation and the temperature deviation Θ, and ∇2 = ∂2/∂x2 + ∂2/∂y2 + ∂2/∂z2 denotes the Laplacian operator. Figure 1 illustrates the local coordinate system of the horizontal shear flow in the radiative zone when the nontraditional f-plane is considered at a given colatitude θ. For the base state, we consider a canonical example of the steady horizontal shear flow U = (U(y),0,0) in a hyperbolic tangent form:

(4)

thumbnail Fig. 1.

(a) Schematic of the radiative and convective zones, colored as yellow and orange, respectively, for the configuration of low-mass stars rotating with angular speed Ω0. (b) Horizontal shear flow U(y) in a local nontraditional f-plane at a colatitude θ in the radiative zone; zc denotes the transition altitude between the radiative and convective zones.

Open with DEXTER

where U0 and L0 are the reference velocity and length scales, respectively. This hyperbolic tangent profile has an inflection point at y = 0. Therefore, it has been considered in many previous stability studies (see e.g., Schmid & Henningson 2001; Deloncle et al. 2007; Arobone & Sarkar 2012) to investigate the inflectional instability. We also adopt this profile to compare with other studies and extend our understanding of the inflectional instability for a horizontal shear when the full Coriolis acceleration is taken into account. In the presence of the horizontal Coriolis parameter fh, 0, such a base flow is steady if a thermal-wind balance between U(y) and the base temperature is satisfied as follows:

(5)

Such a thermal-wind balance has been adopted for cases where the rotation vector and the base shear are not aligned in stratified fluids; for instance, ageostrophic instability in vertical shear flows in stratified fluids in the traditional f-plane (Wang et al. 2014). In addition to the thermal-wind balance, the base temperature is considered to be stably stratified with a linearly increasing profile in z; therefore, it has the form

(6)

where is the difference in base temperature along the vertical distance Δz at a given y.

2.2. Linearized equations

We consider the velocity perturbation , the pressure perturbation where P is the base pressure, and the temperature perturbation . We nondimensionalize variables in Eqs. (1)–(3) by considering the length scale L0, the velocity scale U0, the time scale as t0 = L0/U0, the pressure scale as , and the temperature scale as . For infinitesimally small perturbations, we have the following nondimensional linearized equations:

(7)

(8)

(9)

(10)

(11)

where the prime symbol (′) denotes the derivative with respect to y, f = 2Ω cos θ and are the dimensionless vertical and horizontal Coriolis parameters, respectively, Ω is the dimensionless stellar rotation rate, Re is the Reynolds number

(12)

N is the dimensionless Brunt-Väisälä frequency1

(13)

and Pe is the Péclet number

(14)

In Eq. (11), the term on the left-hand side appears due to the thermal-wind balance (5).

To derive linear stability equations, we apply a normal mode expansion to perturbations as

(15)

where , and are the mode shapes of velocity, pressure, and temperature perturbations, respectively, i2 = −1, kx is the streamwise (longitudinal) wavenumber, kz is the vertical wavenumber, σ = σr + iσi is the complex growth rate where the real part σr is the growth rate and the imaginary part σi is the temporal frequency, and c.c. denotes the complex conjugate. Using the normal mode expansion, we obtain the following linear stability equations

(16)

(17)

(18)

(19)

(20)

where with . By eliminating using the continuity Eq. (16), Eqs. (16)–(20) can be gathered into a matrix form as an eigenvalue problem:

(21)

where 𝒜 and ℬ are the operator matrices expressed as

(22)

(23)

where

(24)

(25)

(26)

(27)

We discretize numerically the operators 𝒜 and ℬ in the y-direction using the rational Chebyshev function that maps the Chebyshev domain ycheb ∈ ( − 1, 1) onto the physical space y ∈ ( − ∞, ∞) via the mapping where Lmap is the mapping factor (Deloncle et al. 2007; Park et al. 2020). To distinguish physical and numerically-converged modes from spurious modes, we use a convergence criterion based on the residual of coefficients on the Chebyshev functions proposed by Fabre & Jacquin (2004). This technique allows us to separate highly-oscillatory spurious modes from physical modes that decays smoothly at boundaries as y → ±∞. Furthermore, we chose the number of collocation points in the y-direction from 100 to 200, which was sufficient in our study to confirm the convergence of physical modes. We impose vanishing boundary conditions as y → ±∞ by suppressing terms in the first and last rows of the operator matrices (Antkowiak 2005; Park 2012). Numerical results of the horizontal shear instability are compared and validated with results of Deloncle et al. (2007), Arobone & Sarkar (2012), Park et al. (2020) in stratified and rotating fluids for the traditional case when .

While the stability of horizontal shear flows in stratified-rotating fluids in the traditional f-plane has one symmetry condition on σ with ±kx and ±kz (Deloncle et al. 2007; Park et al. 2020), an analogous symmetry condition does not exist due to the horizontal Coriolis parameter . Instead, there exist two separate symmetry conditions on σ in terms of kx and kz as follows:

(28)

and

(29)

where the asterisk (*) denotes the complex conjugate. Therefore, in this paper, we investigate both negative and positive kz for kx ≥ 0.

2.3. Simplified equations for

Equations (16)–(20) can be simplified into ordinary differential equations (ODEs) for in particular cases. For instance, in the inviscid and nondiffusive limits (Re → ∞ and Pe → ∞), we can derive the single second-order ODE for as follows:

(30)

where

(31)

s = −iσ + kxU is the Doppler-shifted frequency, ,

(32)

and

(33)

For finite Pe in the inviscid limit, we can find for kx = 0 the single fourth-order ODE as follows:

(34)

where

(35)

(36)

(37)

(38)

In Sect. 4, the above single ODEs (30) and (34) will be used for asymptotic analyses with the WKBJ approximation for large kz to derive asymptotic dispersion relations for the inertial instability.

3. General stability results

In this section, we investigate examples of numerical stability results to see how the full Coriolis acceleration modifies the instabilities of the horizontal shear flow. Figure 2 shows spectra of the eigenvalue σ for various sets of (kx, kz) at Re = ∞,  N = 1, and Pe = 1 for Coriolis components f = 0.5 and (i.e., 2Ω ≃ 2.062 and θ ≃ 76°). In the eigenvalue spectra, the real part of the complex growth rate σr determines the stability while the imaginary part σi denotes the frequency of the mode. In Fig. 2a for (kx, kz) = (0.25, 0), we see different types of eigenvalues: widely distributed stable eigenvalues, clusters of neutral eigenvalues, and one unstable eigenvalue. The neutral or stable modes are continuous modes, which are not exponentially decreasing with |y| → ∞, or gravito-inertial waves that can possess critical points (Astoul et al. 2020). However, this paper only investigates unstable modes that are the most likely to cause turbulence with horizontal shear.

thumbnail Fig. 2.

Eigenvalue spectra in the space (σi, σr) at Re = ∞,  N = 1, Pe = 1, f = 0.5, and (i.e., Ω = 1.031 and θ = 76°) for (a) (kx, kz) = (0.25, 0), (b) (kx, kz) = (0, 6), and (c) (kx, kz) = (0.25, 6). Triangles denote the maximum growth rates.

Open with DEXTER

The eigenfunction of the most unstable mode in Fig. 2a is plotted in Fig. 3 in panels a and b. The mode shape has maxima around y ≃ ±2.96, and the mode displayed in the physical space (x, y) shows that the velocity v(x, y) at z = 0 has an inclined structure in the direction opposite to the shear. This feature is a characteristic of the inflectional instability mode as previously studied for stratified fluids in the traditional f-plane (Arobone & Sarkar 2012; Park et al. 2020).

thumbnail Fig. 3.

Eigenfunctions of the most unstable modes in Fig. 2 at Re = ∞,  N = 1, Pe = 1, f = 0.5, and : (a) the mode shape and (b) velocity v in the physical space (x, y) at z = 0 for (kx, kz) = (0.25, 0), and σ = 0.0712, (c) the mode shape and (d) velocity v in the physical space (y, z) at x = 0 for (kx, kz) = (0, 6), and σ = 0.250. In (a) and (c), black, blue, and red lines denote the absolute value, real part, and imaginary part of the mode shape , respectively.

Open with DEXTER

Figure 2b shows an eigenvalue spectrum that has nearly-neutral clustered eigenvalues, four neutral modes with frequency |σi| > 2, a stable eigenvalue, and one unstable eigenvalue. We note that f = 0.5 belongs to the inertial instability regime for kz >  0 in the traditional approximation, and the instability is not from the inflectional point of the base flow since kx = 0. For the unstable eigenvalue, we plot the corresponding mode in panels c and d of Fig. 3. We see that the mode shape is oscillatory in the region |y| < 10, and this wavelike behavior centered around y = 0 reminds us of the inertial instability mode (Park et al. 2020). We also observe clearly in Fig. 3d a characteristic alternating pattern of v when plotted in the physical space (y, z).

In Fig. 2c, we plot the eigenvalue spectra for (kx, kz) = (0.25, 6). They have weakly unstable or stable clustered eigenvalues due to the critical point where the Doppler-shifted frequency becomes zero (i.e., s = 0), separate neutral and stable eigenvalues, and one unstable eigenvalue. The mode arisen by the critical point is beyond the scope of this study, and it will be covered by the follow-up work (Astoul et al. 2020). At the wavenumbers kz = 6 and kx = 0.25, the inflectional instability is not present and we only have the inertial instability with the mode shape (not shown) similar to that in Fig. 3c.

Figure 4 shows contours of the growth rate σr for the most unstable mode in the parameter space (kx, kz) for different sets of parameters at N = 1 and Re = ∞. We found that the most unstable modes have zero temporal frequency σi; therefore, we plot only the real part of the maximum growth rate σr. In Fig. 4a, we consider the nondiffusive limit Pe = ∞ in an inertially-stable regime at f = 0 according to the criterion in the traditional approximation (Arobone & Sarkar 2012). If we consider (i.e., 2Ω = 0), only the inflectional instability exists and its maximum growth rate σmax = 0.1897 is attained at kx = 0.445 and kz = 0 (see also, Deloncle et al. 2007). As increases (i.e., on the equator for rotating stars with 2Ω >  0), we see that the regime of the inflectional instability is widened in the parameter space (kx, kz), and the higher maximum growth rate σmax ≃ 0.277 is attained at kz = 0 and at a higher value of kx ≃ 0.54.

thumbnail Fig. 4.

Contours of the maximum growth rate max(σr) in the parameter space of (kx, kz) (solid lines) for Re = ∞ and N = 1 for (a) , (b) , and (c) . For (a) and (b), this corresponds to (Ω, θ) = (1, 90° ), and (Ω, θ) = (1.031, 76°) for (c). Black and blue solid lines are contours of the maximum growth rate of the inflectional and inertial instabilities, respectively, and gray dashed lines in the background are contours of the maximum growth rate for the same parameters but in the traditional f-plane approximation where .

Open with DEXTER

In Fig. 4b, we now consider a finite Péclet number Pe = 1. At (i.e., 2Ω = 0), the instability is stabilized as the thermal diffusivity increases (see also, Park et al. 2020) and the unstable regime in the parameter space (kx, kz) is slightly shrunk compared to that at Pe = ∞ in Fig. 4a. The maximum growth rate of the inflectional instability at Pe = 1 and still remains as σmax = 0.1897 at kx = 0.445. As increases, we see that the inflectionally-unstable regime at is even more shrunk than the unstable regime at . The maximum growth rate is decreased to σmax ≃ 0.0717 and is found at (kx, kz) ≃ (0.26, 0.02). This implies that a slightly three-dimensional inflectional instability is now more unstable than the two-dimensional inflectional instability at finite Pe and . What is also very interesting in Fig. 4b is that we observe other unstable regimes. The growth rate contours for large |kz| are reminiscent of those for the inertial instability, which has a maximum growth rate as |kz| → ∞ at kx = 0. It is remarkable to observe this inertial instability at f = 0 (i.e., in the inertially-stable regime in the traditional approximation) when the horizontal Coriolis component becomes positive and the fluid is thermally diffusive.

We now consider in Fig. 4c the case f = 0.5, which belongs to the inertially-unstable regime in the traditional f-plane. We see a clear difference between growth-rate contours at (i.e., on the pole with 2Ω = 0.5) and those at (i.e., 2Ω ≃ 2.062 and θ ≃ 76°). While the growth rate smoothly increases as |kz| increases in the traditional case at , the regime of the inflectional instability is well separated from the regime of the inertial instability for . The inflectional instability at has a maximum growth rate σmax ≃ 0.0725 around (kx, kz) ≃ (0.26, 0.07) while the inertial instability has a maximum growth rate as |kz| → ∞ at kx = 0.

In the traditional f-plane approximation, the inflectional instability reaches its maximum growth rate σmax ≃ 0.1897 in the inviscid limit Re = ∞ at kx ≃ 0.445 and kz = 0 (also shown by the gray dashed line in Fig. 5). In this case, the maximum growth rate σmax is independent from f and N (Park et al. 2020). But without the traditional approximation, the growth rate of the inflectional instability now becomes dependent of N and as shown in Fig. 5 at Pe = ∞ and kz = 0. We verified numerically that the most unstable mode is reached for a finite kx at kz = 0 in the parameter space of (kx, kz) in the inviscid and nondiffusive limits (i.e., Re = Pe = ∞). In these limits, we have the following second-order ODE for

(39)

thumbnail Fig. 5.

Inviscid growth rate of the inflectional instability for various parameter sets of and N at f = 0 (i.e., θ = 90°), Pe = ∞, and kz = 0.

Open with DEXTER

We clearly see that Eq. (39) is still independent from f but it depends on N if . For a fixed N, we see in Fig. 5 that the maximum growth rate increases with , while it decreases with N at a fixed . This implies that the stratification stabilizes the inflectional instability, as similarly observed in the traditional approximation (Park et al. 2020). While the wavenumber range of the inflectional instability is 0 <  kx <  1 at , it is noticeable that the instability can sustain for kx >  1 as increases.

Figure 6 shows growth rate curves for various values of Pe at N = 1, f = 0, kx = 0, and (i.e., on the equator with 2Ω = 2). Although f = 0 implies that it is inertially stable in the traditional approximation, it is still unstable for finite Pe at and the growth rate increase as Pe increases. At a fixed Pe, the growth rate increases with kz and the maximum growth rate is reached as kz → ∞. Furthermore, the growth rate curves converge to a single curve as Pe → 0.

thumbnail Fig. 6.

Inviscid growth rate of the inertial instability for various Péclet numbers at N = 1 and kx = 0 for f = 0 and , i.e., (Ω, θ) = (1, 90° ).

Open with DEXTER

Preliminary results presented in this section tell us that nontraditional effects can significantly modify the properties of the inflectional and inertial instabilities. For instance, the unstable regime of the inertial instability is changed as increases. The maximum growth rate of the inflectional instability increases with and depends on N. Similarly to the traditional case at , thermal diffusion destabilizes the inertial instability when . In the following section, we discuss more thoroughly how the inertial instability depends on physical parameters such as f, , N, or Pe using the WKBJ approximation by taking the asymptotic limit kz → ∞ for the inviscid case at Re = ∞. We derive analytical expressions of the dispersion relations for the inertial instability in the nondiffusive (Pe = ∞) and highly-diffusive (Pe → 0) cases to understand how the inertial instability is modified as the horizontal Coriolis parameter increases.

4. Asymptotic description of the inertial instability

Our previous work in Park et al. (2020) with the traditional approximation () was successful to perform detailed analyses on the inertial instability employing the WKBJ approximation in the inviscid limit Re = ∞. For instance, asymptotic dispersion relations were explicitly proposed for large kz in two limits: Pe → 0 and Pe → ∞. From the dispersion relations, we found that the unstable regime of the inertial instability is 0 <  f <  1 and the maximum growth rate is reached as kz → ∞. In the traditional approximation, σmax is related to the epicyclic frequency ωep, which satisfies and leads to the Solberg-Høiland criterion for stability when (Solberg 1936; Høiland 1941). Besides, σmax is analogous to the growth rate suggested in the limit Pr → 0 for the Goldreich-Schubert-Fricke (GSF) instability (Goldreich & Schubert 1967; Fricke 1968) induced by the horizontal and vertical shear (for more details, we refer the reader to Knobloch & Spruit 1982; Maeder et al. 2013; Barker et al. 2019). Besides, the maximum growth rate σmax is attained independently of the stratification and thermal diffusivity (i.e., N and Pe) while the first-order term of the growth rate strongly depends on them. However, this argument is not applicable without the traditional approximation as demonstrated by numerical results in the previous section. Thus, it is imperative to investigate the inertial instability in the nontraditional case to see how the properties of the instability are changed as the horizontal Coriolis parameter increases.

In the following subsections, we perform the WKBJ analysis for large kz at kx = 0 in the two limits of the Péclet number Pe: the nondiffusive case as Pe → ∞ and the higly-diffusive case as Pe → 0. This analysis is an extension of the work by Park et al. (2020) in the traditional f-plane approximation. For the two cases, we describe below the properties of the inertial instability such as dispersion relations, the maximum growth rate, or regimes of the instability in the parameter space.

4.1. WKBJ formulation in the nondiffusive limit Pe→∞

In this subsection, the thermally nondiffusive case with Pe → ∞ is considered. We verified numerically that the inviscid maximum growth rate of the inertial instability is found as |kz| → ∞ at kx = 0, and it has a zero temporal frequency (i.e., σi = 0). For simplicity, we focus on this most unstable case at kx = 0 and consider only kz >  0 due to the symmetry σ(0, kz) = σ(0, −kz) at kx = 0. The second-order ODE for in Eq. (30) at kx = 0 becomes

(40)

We note that the term multiplied by the first derivative has the order of kz and is comparable with all the other terms when the WKBJ approximation is applied for large kz. For convenience, we introduce a new variable prior to the WKBJ analysis

(41)

Then Eq. (40) becomes

(42)

which is now in the form of the Poincaré equation. Applying to Eq. (42) the WKBJ approximation of for large kz:

(43)

we get

(44)

where

(45)

Since , the exponential behavior of determined by the sign of S0 depends on the sign of . On the one hand, the solution is evanescent if :

(46)

where A1 and A2 are constants, and . On the other hand, the solution is wavelike if :

(47)

where B1 and B2 are constants.

The wavelike or evanescent solution behavior changes at turning points yt where . To find yt, it is useful to define the function σt(y) such that

(48)

(see also, Park & Billant 2012, 2013; Park et al. 2017). By finding σt(y) = σ at a given σ, we can find turning points where . An example of σt(y) is shown in Fig. 7 for f = 0.5, , and N = 1. The gray area denotes where the WKBJ solution is wavelike (i.e., ) while the white area denotes where the solution is evanescent (i.e., ). At a given instability growth rate such that |σ| < σt, max, there exist two turning points, one for y >  0 and the other for y <  0. For instance, if the growth rate is σ = 0.250, there exist two turning points: one at yt+ = 1.365 and the other at yt = −1.365. In this case, we can construct an eigenfunction such that the solution decays exponentially as |y| → ∞ while it is wavelike between the two turning points yt±. If the growth rate is greater than the maximum of σt (i.e., |σ| > |σt, max|), the solution is evanescent everywhere with no turning point and we cannot construct an eigenfunction satisfying the decaying boundary conditions as |y| → ∞. Therefore, to construct the eigenfunction, the growth rate has to lie in the range 0 <  |σ| < σt, max where

(49)

thumbnail Fig. 7.

Function σt(y) (solid lines) at N = 1 for (i.e., Ω = 1.031 and θ = 76°). Gray and white areas denote the regimes where is positive and negative, respectively. The dashed line represents an example of the growth rate σ = 0.250.

Open with DEXTER

By performing a turning point analysis in the growth rate range 0 <  σr <  σt, max, we can further derive an asymptotic dispersion relation as done in Park et al. (2020). We consider first the evanescent WKBJ solution outside the turning point y >  yt+:

(50)

where A is a constant. As y approaches the turning point yt+, the WKBJ solution (50) is not valid and a local solution is needed to find the WKBJ solution below the turning point y <  yt+. We use a new scaled coordinate where ϵ is a small parameter defined as , is the derivative of at yt+, and we assume at leading order that around the turning point yt+. The following local equation is obtained

(51)

Solutions of the local equation (51) are the Airy functions: , where a1 and b1 are constants (Abramowitz & Stegun 1972). From the asymptotic behavior of the Airy functions as and of the WKBJ solution as y → yt+, we obtain the matched WKBJ solution in the region yt <  y <  yt+:

(52)

where constants C± satisfy

(53)

Similarly, we consider the WKBJ solution below the lower turning point y <  yt decaying exponentially as y → −∞:

(54)

After the local analysis around yt, we find the matched wavelike solution in the range yt <  y <  yt+:

(55)

where constants B± satisfy

(56)

Finally, matching the wavelike solutions (52) and (55) between two turning points leads to the following dispersion relation in the quantized form

(57)

where m is the positive integer denoting the mode number. This quantized dispersion relation becomes equivalent to that in the traditional approximation as becomes zero (Park et al. 2020).

The right-hand side of Eq. (57) is always fixed at a finite m so that the integral on the left-hand side of Eq. (57) should go to zero as . This implies that the two turning points converge to zero as . Using this property, we can find a more explicit dispersion relation in the Taylor expansion form for the growth rate σ as

(58)

where

(59)

and

(60)

Since σ1 is positive, the first mode with m = 1 is the most unstable mode. There are also higher-order modes with m >  1 as investigated in Park et al. (2020), but their growth rate is smaller than that of the first mode. Thus, we hereafter consider only the m = 1 mode for asymptotic results. We note that σ = σ0 is the maximum growth rate σmax as kz → ∞, and σ0 equals to σt, max of Eq. (49). While the maximum growth rate in the traditional f-plane approximation is , which depends only on f, we see that the maximum growth rate σ0 now depends on f, , and N altogether. For instance, at a fixed N, the maximum growth rate σ0 increases with as shown in Fig. 8a. Compared to the inertially-unstable regime 0 <  f <  1 in the traditional approximation, we take σmax real and positive and find the following range of the inertial instability for fixed N and :

(61)

thumbnail Fig. 8.

(a) Contours of the maximum growth rate σ0 from Eq. (59) in the parameter space for N = 2. Black dashed lines represent the upper limit of the instability from Eq. (61) for different N, and the white dashed line represents fmax from (64) where the maximum growth rate is attained. (b) Contours of σ0 from Eq. (59) in the parameter space at f = 1. (c) Contours of the first-order term σ1 in the parameter space for N = 2 and the first branch at m = 1.

Open with DEXTER

If we use the relations f = 2Ω cos θ and , we have the corresponding range

(62)

which is equivalent to

(63)

if 2Ω >  0.

Moreover, at the vertical Coriolis parameter f = fmax

(64)

we can find the maximum growth rate at fixed and N as

(65)

In Fig. 8b, we see that the maximum growth rate σ0 at fixed f and decreases as N increases, which highlights the stabilizing role of the stratification. We also display in Fig. 8c the dependence of the first-order term σ1 on f and at N = 2. In this case, σ1 increases with f at a fixed and goes to infinity as f reaches its upper limit since σ0 is in the denominator of σ1 and σ0 → 0 as .

Figure 9 shows the growth rate σ versus the vertical wavenumber kz for different values of and N at f = 1 and kx = 0 in the inviscid and nondiffusive limits. The vertical Coriolis parameter f = 1 belongs to the inertially stable regime under the traditional approximation, but we see that it becomes unstable if we consider the nontraditional effects with the horizontal component , and the growth rate increases with , as predicted from the WKBJ analysis. We notice that numerical stability results are in good agreement with asymptotic predictions from Eq. (58), especially for large kz. We also verify that the stratification stabilizes the inertial instability, and the onset of the inertial instability appears at a higher kz as N increases.

thumbnail Fig. 9.

Numerical results (lines) and asymptotic results (circles) from (58) for various values of and N at f = 1, kx = 0, Re = ∞, and Pe = ∞.

Open with DEXTER

4.2. The WKBJ formulation in the limit Pe→0

In this subsection, we now turn our attention to the highly-diffusive case in the limit of Pe → 0. This limit is relevant for stellar interiors. To analyze the inertial instability, we consider the case with kx = 0 and use the fourth-order ODE (34) for the WKBJ analysis. As performed in the previous subsection, we introduce for convenience the variable

(66)

to understand more clearly the exponential or wavelike behavior around turning points. The fourth-order ODE (34) can be rewritten in terms of as follows:

(67)

where the right-hand side of order will be neglected in the WKBJ analysis for large kz. We apply the WKBJ approximation

(68)

to the fourth-order ODE (67), and we find and four solutions for S0 where

(69)

The WKBJ solution with is

(70)

where A3 and A4 are constants. This solution implies that thus simply increases or decreases exponentially as |y| → ∞ without turning points. Therefore, an eigenfunction satisfying the decaying boundary conditions as |y| → ∞ cannot be constructed. Thus, we impose A3 = A4 = 0. The WKBJ solution with can be either evanescent or wavelike depending on the sign of where

(71)

On the one hand, the WKBJ solution is evanescent if :

(72)

where C1 and C2 are constants, and . On the other hand, the WKBJ solution is wavelike if :

(73)

where D1 and D2 are constants. In these expressions, we consider only the leading order term S0 and neglect higher-order terms.

As conducted in the previous subsection, to derive asymptotic dispersion relations, we need to perform analyses around turning points where . To construct an eigenfunction, we first note that is negative as y → ∞ since . We thus have the exponentially decaying WKBJ solution for as follows:

(74)

Around the turning point where , we expand , where is the derivative of at the turning point, , and is a small parameter satisfying . Using this expansion in Eq. (67), we obtain the local equation around :

(75)

whose solutions are the Airy functions , where c1 and d1 are constants. Matching the asymptotic behavior of the Airy functions as and of the WKBJ solution (74) as , we have the following WKBJ solution in the range between the two turning points:

(76)

For y → −∞, we consider the decaying WKBJ solution for as

(77)

From the local solution behavior around the turning point , we can match the two WKBJ solutions (76) and (77). It leads to the following dispersion relation in a quantized form:

(78)

where m0 is the positive integer for the mode number.

Furthermore, we apply the Taylor expansion for σ

(79)

to the quantized dispersion relation (78), and we get

(80)

and

(81)

Since σ1, 0 is positive, we consider hereafter the first mode with m0 = 1 for the asymptotic results of the most unstable mode. In Fig. 10, we display the growth rate σ as a function of the vertical wavenumber kz for various values of f and at small Pe = 0.1 for kx = 0, N = 1, and Re = ∞, and we compare the numerical results with the asymptotic predictions from Eq. (79). We see that they are in very good agreement, especially as kz → ∞.

thumbnail Fig. 10.

Numerical results (lines) of the growth rate σ as a function of kz for various parameter sets of : (1,1) (black, Ω = 0.707 and θ = 45°), (1,2) (blue, Ω = 1.118 and θ = 63.4°), (1,3) (red, Ω = 1.581 and θ = 71.6°), (−0.5, 1) (green, Ω = 0.559 and θ = 116.6°), (−1, 1) (yellow, Ω = 0.707 and θ = 135°), at Pe = 0.1, kx = 0, N = 1, and Re = ∞, and asymptotic results (circles) from Eq. (79) in the limit Pe → 0.

Open with DEXTER

Figure 11 shows how the maximum growth rate σ0, 0 depends on the Coriolis parameters f and . It is very remarkable that for , σ0, 0 is positive for any value of f including negative (f <  0) and large (f >  1) values, which are the two ranges outside the inertially unstable regime in the traditional approximation. The maximum value of σ0, 0 at a fixed is always attained at f = fmax = 0.5 with the corresponding value σmax, 0 = 0.5.

thumbnail Fig. 11.

Contours of σ0, 0 in the parameter space of .

Open with DEXTER

Using the expressions of the Coriolis parameters f = 2Ω cos θ and , we can see how the maximum growth rate σ0, 0 depends on the rotation rate Ω and colatitude θ. Figure 12 shows the maximum growth rate σ0, 0 plotted against the absolute Coriolis parameter 2Ω for various colatitudes. At the northern pole (θ = 0°), the case considered in the traditional approximation, the inertial instability only exists in the range 0 <  2Ω <  1. On the other hand, in the northern hemisphere 0° < θ ≤ 90°, we see that the inertial instability exists in any values of 2Ω >  0, and the growth rate reaches its maximum σmax = 0.5 at 2Ωmax = 0.5/cos θ. After the peak, σmax decreases with 2Ω and asymptotes to a constant value as 2Ω → ∞. By considering the limit 2Ω → ∞, we can find the asymptote of σmax as follows:

(82)

thumbnail Fig. 12.

Maximum growth rate σmax = σ0, 0 of Eq. (80) as a function of the absolute Coriolis parameter 2Ω at various colatitudes θ. Solid lines and dashed lines denote the results in the northern and southern hemispheres, respectively.

Open with DEXTER

In the southern hemisphere 90° < θ <  180°, the growth rate always increases with 2Ω and it reaches the maximum sin θ/2 as 2Ω → ∞. At the southern pole (θ = 180°), we find no inertial instability.

As a summary, it is remarkable that the inertial instability for the high-diffusivity case (Pe → 0) exists in the ranges

(83)

or

(84)

This colatitude range 0° ≤θ <  180° is much wider than the range of Eq. (63), which belongs to the northern hemisphere in the nondiffusive case (Pe = ∞). Since the high-thermal-diffusivity regime is relevant for stellar interiors, the inertial instability due to the horizontal shear can be an important source of turbulence in stars at any colatitude except at the southern pole. This unstable regime in the limit Pe → 0 is equivalent to the regime of the GSF instability that occurs when the rotation profile varies along the rotation axis for small Pe2.

5. Detailed parametric investigation

The detailed WKBJ analysis performed in the previous section provided explicit expressions for the dispersion relation of the inertial instability in two cases: one with no thermal diffusion and the other with high thermal diffusivity, without the traditional approximation. However, we do not fully understand how the inertial and inflectional instabilities are modified in other regimes such as finite Pe, finite Re, etc. In this section, we present more broad and detailed numerical results to understand the parametric behavior of the inflectional and inertial instabilities on Pe, Re, N, f, and .

5.1. Maximum growth rate of the inflectional instability: effects of N, and Pe

In the traditional approximation (i.e., ), the inflectional instability always has the maximum growth rate σmax ≃ 0.1897 in the inviscid limit Re = ∞ at kx ≃ 0.445 and kz = 0, independently from the Coriolis parameter f, the Brunt-Väisälä frequency N, and the Péclet number Pe (Deloncle et al. 2007; Arobone & Sarkar 2012; Park et al. 2020). But as shown in Figs. 4 and 5, the maximum growth rate of the inflectional instability depends on the values of N and Pe if . For instance, we plot in Fig. 13 the maximum growth rate σmax of the inflectional instability over the parameter space (kx, kz) at f = 0 (i.e., on the equator) in the inviscid and nondiffusive limits (i.e., Re = ∞ and Pe = ∞). We see how the maximum growth rate σmax is modified by the stratification for various values of . We found numerically that the maximum growth rate σmax of the inflectional instability is attained for a finite kx at kz = 0 for the nondiffusive case Pe = ∞, and that it is independent of f as explained by Eq. (39). While the maximum growth rate σmax = 0.1897 is invariant for , we found that σmax at a fixed is reduced for weakly stratified fluids with N ≪ 0.7 while σmax surpasses 0.1897 for stratified fluids with N ≳ 0.7. The maximum of σmax over N, namely max(σmax), is thus greater than 0.1897 and it increases with .

thumbnail Fig. 13.

Maximum growth rate σmax of the inflectional instability as a function of the Brunt-Väisälä frequency N for different values of at Re = Pe = ∞.

Open with DEXTER

Figure 14 shows the maximum growth rate of the inflectional instability over N as a function of for various values of Pe at kz = 0 and f = 0. We see that at a finite Pe, max(σmax) remains constant around 0.1897 for small and then increases with . At Pe = 1, we see that the maximum remains constant in the range . The onset of the growth-rate increase is delayed as Pe decreases; therefore, we can say that the inflectional instability is stabilized as the thermal diffusivity increases (i.e., as Pe decreases). This stabilization of the inflectional instability was also reported in the traditional approximation (see e.g., Park et al. 2020), thus we can say that the inflectional instability is always stabilized by the thermal diffusivity for . And we can expect that the inflectional instability will not sustain in stars due to their high thermal diffusivity. As a summary, we display in Table 1 the parametric dependence for the inflectional instability in the nontraditional case (for other parameters f and N in the traditional case, refer to Deloncle et al. 2007; Arobone & Sarkar 2012; Park et al. 2020).

thumbnail Fig. 14.

Maximum of σmax of the inflectional instability over N as a function of for different values of Pe (solid lines) at kz = 0 and f = 0.

Open with DEXTER

Table 1.

Summary table for the inflectional-instability growth rate by its variation with parameters , Pe, and Re

5.2. Maximum growth rate of the inertial instability at finite Pe

Figure 15 shows plots of the growth rate as a function of kz for different values of Pe at N = 1 and kx = 0. Two parameter sets are considered: (i.e., Ω = 1 and θ = 90°) and (i.e., and θ = 135°). These Coriolis parameters belong to the inertially stable regime in the nondiffusive limit Pe → ∞ according to the criterion (61). However, we see that it becomes unstable as Pe becomes finite, and the growth rate increases as Pe decreases. For a given parameter set of , we see that all growth-rate curves for any Pe reach the same maximum value σ0, 0 given by Eq. (80) derived in the limit Pe → 0 as kz → ∞. We also verified numerically that this asymptotic behavior is the same for other values of N at finite Pe. Therefore, the largest growth rate of the inertial instability is always σ0, 0, and we have analogously the same relation between the maximum growth rate of the inertial instability and the Coriolis parameters at finite Pe.

thumbnail Fig. 15.

Growth rate σr of the inertial instability as a function of kz for values of Pe at N = 1 and kx = 0 for (i.e., Ω = 1 and θ = 90°; solid lines) and for (i.e., Ω = 0.707 and θ = 135°; dashed lines). Dotted lines indicate asymptotes σ0, 0 of Eq. (80) predicted from the WKBJ analysis in the limit Pe → 0.

Open with DEXTER

In Table 2, we summarize the variation of the growth rate of the inertial instability with parameters and Pe, and the unstable regimes in the two limits Pe → ∞ and Pe → 0. The destabilization of the inertial instability by Pe was already reported previously (see e.g., Park et al. 2020), and we confirm in this paper that this destabilization holds beyond the traditional approximation.

Table 2.

Summary table for the variation of the inertial-instability growth rate and the inertially unstable regime

5.3. Effect of the viscosity at finite Re

Now, we investigate how the viscosity modifies the inflectional and inertial instabilities. Figure 16 shows the growth rate of the inflectional and inertial instabilities for various parameters at different Reynolds numbers Re. For viscous cases, we fix the Prandtl number Pr defined as Pr = Pe/Re. We clearly see that the viscosity has a stabilizing effect on both instabilities. For the inflectional instability in Fig. 16a, the inviscid growth rate curve is increased for (black solid line) compared to that of the traditional case at (gray dashed line), while the viscous growth rate at decreases as Re decreases. It is remarkable that the inflectional instability persists at low Reynolds number Re = 10, and that the curve at has the same order of magnitude as the inviscid growth rate at . Therefore, we can conclude the inflectional instability is stabilized by the viscosity while it is destabilized as increases. This trend is summarized in Table 1.

thumbnail Fig. 16.

(a) Growth rate σ of the inflectional instability as a function of the streamwise wavenumber kx for kz = 0, N = 1, Pr = 1, , f = 0 (i.e., Ω = 1 and θ = 90°) for Re = ∞ (black), Re = 50 (blue), and Re = 10 (red). The gray dash-dot line denotes the inviscid growth rate of the inflectional instability under the traditional approximation (i.e., ). (b,c) Growth rate σ of the inertial instability as a function of the vertical wavenumber kz for kx = 0, N = 1, , and f = 0.5 (i.e., Ω = 2.02 and θ = 82.9°) for (b) Pr = 1 and Re = ∞ (black), Re = 1000 (blue), and Re = 500 (red), and (c) Pr = 10−6 and Re = 104 (black), Re = 5000 (blue), and Re = 1000 (red). Solid lines are numerical results, and dashed lines in (b) and (c) are asymptotic predictions from Eqs. (86) and (89), respectively.

Open with DEXTER

Figures 16b and c show the growth rate of the inertial instability at kx = 0, N = 1, f = 0.5, and for two Prandtl numbers: Pr = 1 and Pr = 10−6. For both cases, the maximum growth rate is now attained at finite kz due to a stronger stabilization by the viscosity at large wavenumber kz. For large Reynolds numbers and a unity Prandtl number (Pr = 1), we can apply the multiple-scale analysis to perturbation equations and propose the viscous growth rate expressed in terms of the inviscid growth rate as

(85)

This expression of the growth rate has already been proposed in other works (see e.g., Arobone & Sarkar 2012; Yim et al. 2016). Using the asymptotic expression of the inviscid growth rate obtained from the WKBJ analysis in Sect. 4, we can express explicitly the viscous growth rate of the inertial instability at kx = 0 as

(86)

We see in Fig. 16b that this expression for the viscous growth rate (86) is in good agreement with numerical results at finite Reynolds numbers and Pr = 1.

From (86), we can further derive the optimal vertical wavenumber kz, c, at which the maximum growth rate is attained, by solving the equation ∂σ/∂kz = 0:

(87)

Applying the optimal wavenumber kz, c to the viscous growth rate (86), we find the critical Reynolds number Rec above which the inertial instability occurs.

(88)

In Fig. 17a, we display the optimal wavenumber kz, c (87) as a function of the Reynolds number Re and the corresponding critical Rec at 2Ω = 1, N = 4, and kx = 0 at various colatitudes θ. In Fig. 17b, we plot the critical Reynolds number Rec (88) as a function of 2Ω for different colatitudes. Asymptotic predictions from Eq. (88) are of the same order of magnitude as numerical results for Rec, and the small differences come from the small values of kz, c predicted from the WKBJ analysis, which works well for large kz. However, the interest of the asymptotic expression of Eq. (88) is that it can provide information about the critical Reynolds number Rec in wide parameter spaces without exhaustive numerical computations.

thumbnail Fig. 17.

(a) Optimal wavenumbers kz, c from (87) (solid lines) and kz, c0 from (90) (dashed lines) as a function of Re at 2Ω = 1, N = 4, and kx = 0 for colatitudes θ: 30° (blue), 45° (red), and 60° (green). Filled circles denote the critical Reynolds numbers. (b,c) Critical Reynolds number Rec as a function of the absolute Coriolis parameter 2Ω for N = 4, and (b) Pr = 1, (c) Pr = 10−6 at different colatitudes θ: 0° (black), 30° (blue), 45° (red), and 60° (green) from numerical results (crosses) and asymptotic results of Rec (solid lines in b) from (88) and Rec, 0 (dashed lines in c) from (91).

Open with DEXTER

For very small Prandtl numbers, which is the relevant regime for the radiative zone of stars (Lignières 1999), we can similarly propose the viscous growth rate using the asymptotic growth rate (79) in the limit Pe = 0:

(89)

In Fig. 16c, we see that this viscous growth rate agrees very well with numerical results at finite Reynolds numbers and Pr = 10−6. Moreover, we can further derive the optimal wavenumber kz, c0 and the critical Reynolds number Rec, 0 in the limit Pr → 0 as follows:

(90)

(91)

Figure 17a and c show the asymptotic predictions of kz, c0 and Rec, 0 and comparisons with numerical results at Pr = 10−6. The critical Reynolds number is roughly within the same order of magnitude as numerical results, but the difference is large since the predicted critical wavenumber kz, c0 is very small, of order O(1), which lies much below the validity range of the WKBJ approximation. However, we can see that the critical Reynolds number Rec, 0 increases with 2Ω except at the northern pole θ = 0° where it is unstable in the range 0 <  2Ω <  1. We also found both numerically and theoretically that Rec, 0 does not change significantly with colatitude θ at a given Ω.

6. Effective horizontal turbulent viscosity

6.1. Turbulent viscosity induced by the inertial instability

As small-amplitude perturbations grow due to the instability mechanism, they first follow the linear instability growth, then reach a saturated equilibrium state and undergo a transition to turbulence as the perturbation amplitude increases and nonlinear effects become crucial. The nonlinear saturation of instabilities has been investigated in different astrophysical contexts and has been used to predict the order of magnitude of angular momentum transport, mixing, or dynamo actions (Spruit 2002; Denissenkov & Pinsonneault 2007; Zahn et al. 2007; Fuller et al. 2019).

In this section, we use our results on horizontal shear instabilities to derive an effective turbulent viscosity. We consider perturbations that are developed by the horizontal shear instabilities and reach a saturated state due to nonlinear interactions, which we model here through a turbulent viscosity. We adopt the approach described by Spruit (2002), Fuller et al. (2019) that allows the saturated state when the turbulent damping rate γturb equals the growth rate of the instability σ. This leads to the following expression of the local effective turbulent viscosity in the horizontal direction νh, h:

(92)

where is the horizontal wavenumber in the y-direction required for instability. We use here the notation νh, h first introduced by Mathis et al. (2018) that corresponds to the turbulent transport triggered in the horizontal direction by the instability of a latitudinal horizontal shear. We thus make the choice to consider the horizontal wavenumber in the latitudinal y-direction to characterize the horizontal effective turbulent viscosity.

We first consider the case of the effective turbulent viscosity νh, h induced by the inertial instability, whose dispersion relation has been derived using the WKBJ approximation in Sect. 4. Since the horizontal shear flow U(y) is inhomogeneous in the y-direction, the associated wavenumber ky in the y-direction is not a constant but also a function of y. Moreover, the local wavenumber ky(y) only exists when the solution is wavelike in the regime between the two turning points. For the nondiffusive case at Pe = ∞, we have the WKBJ expression for as

(93)

Since the exponent changes with y, we define the horizontal wavenumber as the average of the absolute part of the exponent over the two turning points as follows:

(94)

Similarly, we compute νh, h using the asymptotic dispersion relation (79) in the highly diffusive limit Pe → 0:

(95)

Figure 18 shows an example of as a function of in the two limits Pe = ∞ (solid line) and Pe → 0 (dashed line). The asymptotic growth rates σ from Eq. (58) for Pe = ∞ and from Eq. (79) for Pe → 0 derived with the WKBJ analyses are used to compute νh, h. We see that the effective turbulent viscosity reaches its maximum at a finite and decreases monotonically with . The monotonic decrease of νh, h is due to the fact that the growth rate reaches asymptotically the maximum growth rate σmax as kz → ∞ while the horizontal wavenumber , proportional to kz, goes to infinity.

thumbnail Fig. 18.

Effective turbulent viscosity νh as a function of the horizontal wavenumber for f = 0.5, (i.e., Ω = 0.559 and θ = 63.4°), and kx = 0 for Pe = ∞ and N = 1 (solid line), and as Pe → 0 (dashed line).

Open with DEXTER

Similarly to Fuller et al. (2019), where the minimum of the wavenumber is considered to get νh, h, we consider the maximum peak of νh, h as the representative value of the effective horizontal turbulent viscosity for given values of the physical parameters such as f and . In Fig. 19a and b, we plot contours of max(νh, h) in the parameter space of (N/2Ω, θ) for two different stratified fluids with N = 1 and N = 2 without thermal diffusion (i.e., Pe = ∞). For both cases, the turbulent viscosity max(νh, h) reaches its maximum around the colatitude θ ≃ 80° in the northern hemisphere close to the equator around the value N/2Ω ∼ 1. The effective turbulent viscosity in the southern hemisphere θ >  90° is zero since there is no instability for negative f = 2Ω cos θ. This implies that a strong turbulence and a large effective turbulent viscosity are expected near the equator in the northern hemisphere due to strong instability. In Fig. 19c, we plot max(νh, h) in the parameter space of (2Ω, θ) for high-diffusivity fluids in the limit Pe → 0 using Eq. (95). In this case, the effective turbulent viscosity is not zero in the southern hemisphere due to the presence of the inertial instability for negative f. The effective turbulent viscosity max(νh, h) has a maximum around θ ≃ 70°. This emphasizes that nontraditional effects with the positive horizontal Coriolis parameter play a crucial role on the effective turbulent viscosity, especially near the equator.

thumbnail Fig. 19.

(a,b) Contours of the maximum of the effective turbulent viscosity max(νh, h) in the parameter space (N/2Ω, θ) for Pe = ∞ and (a) N = 1, (b) N = 2. (c) Contours of the maximum of the effective turbulent viscosity max(νh, h) in the parameter space (2Ω, θ) as Pe → 0.

Open with DEXTER

From νh, h, the latitudinally-averaged turbulent viscosity can be computed in two ways. On the one hand, we define the averaged turbulent viscosity following Zahn (1992) as

(96)

This average was used to compute the transport of angular momentum. On the other hand, we follow the definition by Mathis et al. (2018):

(97)

where is the averaged turbulent viscosity obtained in the way the mean transport of chemicals is computed. In Fig. 20, we plot these averaged viscosities. For the nondiffusive case (Pe = ∞), is larger than and they are reduced as N increases. The maxima of the two viscosities are reached around N/2Ω = 1 for both N = 1 and N = 2. For Pe → 0, the maximum is much lower than that for Pe = ∞, and it is reached around 2Ω = 1.

thumbnail Fig. 20.

Latitudinally-averaged turbulent viscosities (black) and (blue) computed from contours in Fig. 19.

Open with DEXTER

6.2. Turbulent viscosity induced by the inflectional instability

The inflectional instability can also induce a turbulent transport in stellar radiation zones. Therefore, we also compute the effective turbulent viscosity based on the growth rate σ of the inflectional instability. The issues with the inflectional instability are that we do not have an analytic expression for the growth rate σ, and it is difficult to define systematically the horizontal wavenumber . One possible candidate for is to choose based on the solution’s asymptotic behavior as |y| → ∞. From numerical computations of the maximum growth rate of the inflectional instability in the parameter space (kx, kz), we display in Fig. 21 the effective turbulent viscosity in the parameter space of at Pe = ∞ and f = 0. It is important to note that the maximum growth rate σmax of the inflectional instability is independent of f at Pe = ∞. We see that at a fixed N, max(νh, h) increases with for a weak stratification with N <  1 while it decreases with for a strong stratification with N >  1. The effective turbulent viscosity max(νh, h) is generally of order of the unity O(1) in the parameter space of and is smaller than the effective turbulent viscosity max(νh, h) induced by the inertial instability. We also verified numerically that the max(νh, h) of the inflectional instability decreases as the thermal diffusivity becomes finite and small. Therefore, we can conjecture that the effective turbulent viscosity induced by the inertial instability will be larger than that induced by the inflectional instability.

thumbnail Fig. 21.

Contours of the maximum of the effective turbulent viscosity max(νh, h) induced by the inflectional instability in the parameter space of for f = 0 (i.e., on the equator at θ = 90°) and Pe = ∞.

Open with DEXTER

7. Discussion with dimensional parameters

7.1. Inertial instability with latitudinal differential rotation

We discussed in the previous sections how the inertial instability is developed at a given latitude θ when the base shear of Eq. (4) is considered. In our analysis, we investigated the growth rate of the inertial instability in the dimensionless form by taking the velocity scale U0, length scale L0, and time scale t0 = L0/U0. We note that the length scale L0 is different from that of global large-scale shear flows in stellar interiors as we consider small-scale shear flows and associated perturbations on a local plane. An important thing to note is that the base flow with positive shear S0 = U0/L0 at y = 0 is used in the normalization. This implies that the velocity in the azimuthal direction always decreases as the colatitude θ increases. The latitudinal differential rotation in stars is, however, different from this case. For instance, stars can have a conical differential rotation in the simplest form as

(98)

where ϵ >  0 corresponds to the solar-like case, in which the equator rotates faster than the pole (e.g., ϵ ≃ 0.3 for the Sun), while ϵ <  0 corresponds to the antisolar case, where the equator rotates slower than the pole (Guenel et al. 2016). If instead of being constant, ϵ is proportional to r2, one obtains the simplest form of cylindrical differential rotation (Baruteau & Rieutord 2013). At a constant radius, the latitudinal differential rotation follows the same law as in the conical case.

In the local frame rotating with Ω0 at the radius r = R and the colatitude θ, the relative azimuthal velocity Uφ is

(99)

(see also Hypolite et al. 2018). By considering the increment in the latitudinal direction on the local frame:

(100)

we obtain the base shear

(101)

The latitudinal shear S is negative (resp. positive) in the northern hemisphere and positive (resp. negative) in the southern hemisphere if ϵ >  0 (resp. ϵ <  0). Moreover, the latitudinal shear S is zero at the equator since the differential rotation of Eq. (98) is either at the maximum when ϵ >  0 or the minimum when ϵ <  0.

The growth rate of Eq. (79) that we derived from the WKBJ analysis is local but the growth rate induced by the latitudinal differential rotation of Eq. (98) needs to be obtained by the global stability analysis (see e.g., Guenel et al. 2016). However, we can approximately predict the growth rate σθ in the stellar radiation zones using the latitudinal shear S and the maximum growth rate of Eq. (80), which can be expressed in a dimensional form as follows:

(102)

Equation (102) can further be expanded in terms of the stellar rotation Ω0, the colatitude θ, and ϵ as

(103)

Figure 22 shows the growth rate σθ versus the colatitude θ for various ϵ. The growth rate σθ increases with |ϵ| and reaches its maximum around θ ≃ 60° and 120°. It is zero at the poles and the equator. We see that, for the same |ϵ|, the growth rate of the antisolar case for ϵ <  0 is slightly larger than that of the solar-like case for ϵ >  0.

thumbnail Fig. 22.

Growth rate σθ as a function of the colatitude θ for various values of ϵ. Solid lines denote σθ with ϵ >  0 and dashed lines denote σθ with ϵ <  0.

Open with DEXTER

Furthermore, the latitudinally-averaged3 growth rate can be computed as

(104)

In Fig. 23, we display as a function of ϵ for the solar-like and antisolar cases. In the range |ϵ| < 1, almost does not depend on the sign of ϵ. The growth rate shows a linear relation with |ϵ| as

(105)

thumbnail Fig. 23.

Latitudinally-averaged growth rate as a function of |ϵ| for the solar-like case (ϵ >  0, solid line) and the antisolar case (ϵ <  0, dashed line).

Open with DEXTER

where aϵ ≃ 0.16. We can derive a similar scaling law for the growth rate in the limit Pe → ∞, more relevant to the solar tachocline case (Garaud 2020), if we use the growth rate (59) and a proper value of N.

7.2. Turbulent viscosities and characteristic time of turbulent transport

Using the averaged growth rate , we can define the turbulent viscosities induced by the horizontal shear due to the latitudinal differential rotation of Eq. (98). Consistently with the notations from Mathis et al. (2018), we define the turbulent viscosity in the latitudinal (horizontal) direction and the turbulent viscosity in the radial (vertical) direction . Following the assumption proposed by Spruit (2002), Fuller et al. (2019) that the shear instability grows and the momentum balances in nonlinear regime with the turbulent Reynolds stress, we find

(106)

(107)

where ∥ denotes the direction parallel to the stratification, ⊥ denotes the direction perpendicular to the stratification, k and k are the characteristic wavenumbers in the horizontal and vertical directions, and l and l are the length scales in the horizontal and vertical directions, respectively. We note here that we consider small length scales of the unstable modes that trigger turbulent transport, the length scales different from those of global large-scale shear flows. The eddy viscosities are directly proportional to the horizontal shear as in the prescription derived by Mathis et al. (2004) (see their Eq. (18)). We note that it would also be possible to define diffusion coefficients for chemicals and ; for this one should compute .

While the turbulent viscosities depend on the direction with the scaling

(108)

(see also, Mathis et al. 2018), the dynamical time scale τ that characterizes the turbulence induced by the horizontal shear, and the transport of momentum and chemicals it triggers both in the vertical and latitudinal directions, is

(109)

which is simply the inverse of the growth rate .

Using the relation of Eq. (105), we can express the characteristic time τ in terms of Ω0 as

(110)

where τ0 denotes the rotation period of the star at the pole. We note that the characteristic time scale τ is similar to that proposed by Mathis et al. (2018) in the form τ = 1/S where S(r, θ) = r sin θrΩ characterizes the radial (vertical) shear.

According to the scaling (110), we estimate the transport time τ for the Sun at the level of the tachocline as 44.5 days when we use Ω0 ≃ 433 nHz and ϵ = 0.3 (García et al. 2007). For solar-like stars within the mass range 0.9 − 1.1 M (M is the Solar mass), if we make the rough assumption that the latitudinal differential rotation with ϵ = 0.3 is maintained during the evolution4, we can roughly estimate that the transport time τ for the solar-like stars in the case of a slow initial rotation (Gallet & Bouvier 2015) is about 13.3 days, 9.7 days, and 16.8 days at the ages of 10 Myr, 100 Myr, and 1000 Myr, respectively. These time scales will be eventually used to compute the eddy viscosities in (109) upon the choice of the length scales.

For stars with a convective core, numerical simulations showed that the differential rotation in the core is mostly cylindrical (Browning et al. 2004; Augustson et al. 2016). For a 2 M star, Browning et al. (2004) found a rotation contrast within the core of between 30 and 60%. At the boundary between the convective core and the surrounding radiative zone, this is equivalent to ϵ ≃ 0.3 − 0.6. Using this estimate as well as rotation rates typical of A-type stars computed using observed surface velocities (see e.g., Zorec & Royer 2012) and MESA stellar structure and evolution models of a 2 M star with a Solar metallicity (Paxton et al. 2011) that provide us the stellar radius, the characteristic time τ ranges between 0.4 and 1 days during most of the main sequence and can go up to 1.5 days the end of the main sequence (around 1 Gyr).

If the turbulence triggered by the horizontal shear acts to damp its source, the horizontal differential rotation, as proposed by Zahn (1992) we thus predict a very efficient transport of angular momentum that can lead to the so-called shellular rotation where horizontal gradients of the angular velocity are weak. This also may be the origin of the observed very small radial extent of the Solar tachocline (Spiegel & Zahn 1992) and of an efficient mixing of chemicals in this region (Brun et al. 1999).

The behavior of the turbulence generated by horizontal shears has been recently explored in the nonrotating case by Cope et al. (2020) and Garaud (2020). They find when exploring the parameter space using Direct Numerical Simulations (DNS) a turbulent stratified regime where the turbulent transport can be modeled by an eddy diffusivity. This opens an interesting path to verify our predictions when such DNS will take rotation into account.

8. Conclusion

In this paper, we studied horizontal shear flow instabilities in stably-stratified and thermally-diffusive fluids on a rotating plane where the full Coriolis acceleration with both vertical and horizontal components of the rotation vector is taken into account. For the canonical shear flow in the hyperbolic tangent profile, there exist two types of shear instabilities: the inflectional instability due to the presence of an inflection point and the inertial instability due to an inertial imbalance in the presence of the Coriolis acceleration. In the presence of positive horizontal Coriolis parameter , we found that both the inflectional and inertial instabilities are strongly affected. For instance, the inflectional instability, whose maximum growth rate is known to be independent of N, f, and Pe in the traditional approximation (Deloncle et al. 2007; Park et al. 2020), has now a maximum growth rate that depends on N and Pe. The horizontal Coriolis parameter destabilizes the inflectional instability for strong stratification while it stabilizes the instability at a small N. The thermal diffusivity at finite Pe also plays a stabilizing role in the inflectional instability. On the other hand, the inertial instability is destabilized by the thermal diffusivity, and the unstable regime is widely extended. For instance, in the nondiffusive limit Pe = ∞, the inertially unstable regime is found to be (i.e., tan−1(N/2) < θ <  90°). More strikingly, it is inertially unstable for any f for high diffusivity fluids as Pe → 0 when (i.e., the inertial instability is active at any colatitude 0 <  θ <  180° for the absolute Coriolis parameter 2Ω >  0). These unstable regimes as well as the dispersion relations for the inertial instability are derived analytically using the WKBJ approximation for large vertical wavenumber kz, and this asymptotic analysis demonstrates a very good agreement with numerical results in the inviscid limit. Using the asymptotic expressions of the growth rate for the inertial instability, we also predicted the critical Reynolds number above which the growth rate becomes positive. Finally, we proposed prescriptions for the effective horizontal and vertical turbulent viscosities induced by the inertial and inflectional instabilities and found that the inertial instability plays an important role in the turbulent transport near the equator.

Observational and numerical studies suggest that stellar radiative zones have a mostly uniform rotation, whereas stellar convective zones are differentially rotating, but neutrally stratified. Therefore, one expects the present study to be relevant near the boundary between the convective and radiative zones. Because of the small values predicted for the time that characterizes the turbulence triggered by horizontal shear flows in such regions, this turbulent transport can have a strong impact on the structure and the evolution of stars, for instance by interacting with overshooting flows and by extracting angular momentum and chemical elements from the core to the envelope. In particular, the horizontal shear could be a crucial physical ingredient needed to explain the observed structure of the solar and stellar tachoclines (Spiegel & Zahn 1992; Brun et al. 1999). Other ingredients that are missing in the present work, such as magnetic fields (Gough & McIntyre 1998; Strugarek et al. 2011; Acevedo-Arreguin et al. 2013; Barnabé et al. 2017), may also play an important role.

Our predictions concerning the turbulent viscosities need to be confirmed by fully turbulent numerical simulations. In particular, global numerical simulations would allow us to validate the local approach used in this work, and the relevance of the latitudinally averaged quantities derived for stellar evolution codes. Besides, the local results could be used to build subgrid models for large-eddy simulations and stellar-evolution calculations to better capture small-scale transport processes.


1

The nondimensional parameter N2 is equivalent to the Richardson number Ri with our choice of normalization.

2

The configuration for the GSF instability is relevant to our case for the inertial instability where the rotation varies along the latitudinal y-direction due to the local horizontal shear we consider.

3

We here use the latitudinal average defined by Zahn (1992) for quantities related to angular momentum. For instance, he defined the shellular rotation as .

4

Using scaling laws computed in Brun et al. (2017) and grids of stellar models that take rotation into account using the STAREVOL code (Amard et al. 2019), Astoul et al. (2020) demonstrates that 0.15 <  |ϵ| < 0.4 that will not predict orders of magnitude for τ

Acknowledgments

The authors acknowledge support from the European Research Council through ERC grant SPIRE 647383 and from GOLF and PLATO CNES grants at the Department of Astrophysics at CEA Paris-Saclay. We thank the referee, Prof. A. J. Barker, for his constructive comments that allow us to improve the article.

References

  1. Abramowitz, M., & Stegun, I. A. 1972, Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables (Dover Publications) [Google Scholar]
  2. Acevedo-Arreguin, L. A., Garaud, P., & Wood, T. S. 2013, MNRAS, 434, 720 [Google Scholar]
  3. Aerts, C., Mathis, S., & Rogers, T. M. 2019, ARA&A, 57, 35 [Google Scholar]
  4. Amard, L., Palacios, A., Charbonnel, C., et al. 2019, A&A, 631, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Antkowiak, A. 2005, Ph.D. Thesis, Université Paul Sabatier de Toulouse [Google Scholar]
  6. Arobone, E., & Sarkar, S. 2012, J. Fluid Mech., 703, 29 [Google Scholar]
  7. Astoul, A., Park, J., Mathis, S., Baruteau, C., & Gallet, F. 2020, A&A, submitted [Google Scholar]
  8. Augustson, K. C., Brun, A. S., & Toomre, J. 2016, ApJ, 829, 92 [Google Scholar]
  9. Barker, A. J., Jones, C. A., & Tobias, S. M. 2019, MNRAS, 487, 1777 [Google Scholar]
  10. Barker, A. J., Jones, C. A., & Tobias, S. M. 2020, MNRAS, 495, 1468 [Google Scholar]
  11. Barnabé, R., Strugarek, A., Charbonneau, P., Brun, A. S., & Zahn, J.-P. 2017, A&A, 601, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Baruteau, C., & Rieutord, M. 2013, J. Fluid Mech., 719, 47 [Google Scholar]
  13. Belkacem, K., Marques, J. P., Goupil, M. J., et al. 2015a, A&A, 579, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  14. Belkacem, K., Marques, J. P., Goupil, M. J., et al. 2015b, A&A, 579, A31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Browning, M. K., Brun, A. S., & Toomre, J. 2004, ApJ, 601, 512 [Google Scholar]
  16. Brun, A. S., & Browning, M. K. 2017, Liv. Rev. Sol. Phys., 14, 4 [Google Scholar]
  17. Brun, A. S., Strugarek, A., Varela, J., et al. 2017, ApJ, 836, 192 [Google Scholar]
  18. Brun, A. S., Turck-Chièze, S., & Zahn, J. P. 1999, ApJ, 525, 1032 [Google Scholar]
  19. Charbonneau, P., & MacGregor, K. B. 1993, ApJ, 417, 762 [Google Scholar]
  20. Cope, L., Garaud, P., & Caulfield, C. P. 2020, J. Fluid Mech., 903, A1 [Google Scholar]
  21. Deloncle, A., Chomaz, J.-M., & Billant, P. 2007, J. Fluid Mech., 570, 297 [Google Scholar]
  22. Denissenkov, P. A., & Pinsonneault, M. 2007, ApJ, 655, 1157 [Google Scholar]
  23. Eckart, C. 1960, Hydrodynamics of Oceans and Atmospheres (Elsevier) [Google Scholar]
  24. Ekström, S., Georgy, C., Eggenberger, P., et al. 2012, A&A, 537, A146 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Fabre, D., & Jacquin, L. 2004, J. Fluid Mech., 500, 239 [Google Scholar]
  26. Fricke, K. 1968, Z. Astrophys., 68, 317 [NASA ADS] [Google Scholar]
  27. Fuller, J., Piro, A. L., & Jermyn, A. S. 2019, MNRAS, 485, 3661 [NASA ADS] [Google Scholar]
  28. Gagnier, D., & Garaud, P. 2018, ApJ, 862, 36 [Google Scholar]
  29. Gagnier, D., Rieutord, M., Charbonnel, C., Putigny, B., & Espinosa Lara, F. 2019, A&A, 625, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Gallet, F., & Bouvier, J. 2015, A&A, 577, A98 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Gallet, F., Bolmont, E., Mathis, S., Charbonnel, C., & Amard, L. 2017, A&A, 604, A112 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Garaud, P. 2020, ApJ, 901, 146 [Google Scholar]
  33. Garaud, P., Gagnier, D., & Verhoeven, J. 2017, ApJ, 837, 133 [Google Scholar]
  34. García, R. A., Turck-Chièze, S., Jiménez-Reyes, S. J., et al. 2007, Science, 316, 1591 [Google Scholar]
  35. Gerkema, T., & Shrira, V. I. 2005, J. Fluid Mech., 529, 195 [Google Scholar]
  36. Gerkema, T., Zimmerman, J. T. F., Maas, L. R. M., & van Haren, H. 2008, Rev. Geophys., 46, RG2004 [Google Scholar]
  37. Goldreich, P., & Schubert, G. 1967, ApJ, 150, 571 [Google Scholar]
  38. Gough, D. O., & McIntyre, M. E. 1998, Nature, 394, 755 [Google Scholar]
  39. Griffiths, S. D. 2008, J. Fluid Mech., 605, 115 [NASA ADS] [Google Scholar]
  40. Guenel, M., Baruteau, C., Mathis, S., & Rieutord, M. 2016, A&A, 589, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Hirschi, R., Meynet, G., & Maeder, A. 2004, A&A, 425, 649 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Høiland, E. 1941, in Avhandgliger Norske Videnskaps-Akademi i Oslo, I, Math.-Naturv. Klasse, 11, 1 [Google Scholar]
  43. Hypolite, D., Mathis, S., & Rieutord, M. 2018, A&A, 610, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Knobloch, E., & Spruit, H. C. 1982, A&A, 113, 261 [NASA ADS] [Google Scholar]
  45. Kulenthirarajah, L., & Garaud, P. 2018, ApJ, 864, 107 [Google Scholar]
  46. Lignières, F. 1999, A&A, 348, 933 [NASA ADS] [Google Scholar]
  47. Maeder, A. 2003, A&A, 399, 263 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Maeder, A. 2009, Physics, Formation and Evolution of Rotating Stars. [Google Scholar]
  49. Maeder, A., & Zahn, J.-P. 1998, A&A, 334, 1000 [NASA ADS] [Google Scholar]
  50. Maeder, A., Meynet, G., Lagarde, N., & Charbonnel, C. 2013, A&A, 553, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Marques, J. P., Goupil, M. J., Lebreton, Y., et al. 2013, A&A, 549, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Mathis, S., & Zahn, J. P. 2004, A&A, 425, 229 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Mathis, S., Palacios, A., & Zahn, J.-P. 2004, A&A, 425, 243 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Mathis, S., Neiner, C., & Tran Minh, N. 2014, A&A, 565, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Mathis, S., Prat, V., Amard, L., et al. 2018, A&A, 620, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Matt, S. P., Brun, A. S., Baraffe, I., Bouvier, J., & Chabrier, G. 2015, ApJ, 799, L23 [Google Scholar]
  57. Moss, D. 1992, MNRAS, 257, 593 [NASA ADS] [Google Scholar]
  58. Park, J. 2012, Ph.D. Thesis, Ecole Polytechnique [Google Scholar]
  59. Park, J., & Billant, P. 2012, J. Fluid Mech., 707, 381 [Google Scholar]
  60. Park, J., & Billant, P. 2013, J. Fluid Mech., 725, 262 [Google Scholar]
  61. Park, J., Billant, P., & Baik, J.-J. 2017, J. Fluid Mech., 822, 80 [Google Scholar]
  62. Park, J., Prat, V., & Mathis, S. 2020, A&A, 635, A133 [EDP Sciences] [Google Scholar]
  63. Paxton, B., Bildsten, L., Dotter, A., et al. 2011, ApJS, 192, 1 [Google Scholar]
  64. Pinçon, C., Belkacem, K., Goupil, M. J., & Marques, J. P. 2017, A&A, 605, A31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Prat, V., & Lignières, F. 2013, A&A, 551, L3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Prat, V., & Lignières, F. 2014, A&A, 566, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Prat, V., Guilet, J., Viallet, M., & Müller, E. 2016, A&A, 592, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Richard, D., & Zahn, J.-P. 1999, A&A, 347, 734 [Google Scholar]
  69. Rogers, T. M. 2015, ApJ, 815, L30 [Google Scholar]
  70. Schmid, P., & Henningson, D. S. 2001, Stability and Transition in Shear Flows (New York: Springer-Verlag) [Google Scholar]
  71. Solberg, H. 1936, in Procès Verbaux Ass. Météor., UGGI, 6e Assemblée Générale, Edinburgh, Mém. et Disc., 2, 66 [Google Scholar]
  72. Spiegel, E. A., & Zahn, J. P. 1992, A&A, 265, 106 [NASA ADS] [Google Scholar]
  73. Spruit, H. C. 1999, A&A, 349, 189 [NASA ADS] [Google Scholar]
  74. Spruit, H. C. 2002, A&A, 381, 923 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Strugarek, A., Brun, A. S., & Zahn, J. P. 2011, A&A, 532, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  76. Strugarek, A., Bolmont, E., Mathis, S., et al. 2017, ApJ, 847, L16 [Google Scholar]
  77. Talon, S., & Charbonnel, C. 2005, A&A, 440, 981 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  78. Talon, S., & Zahn, J. P. 1997, A&A, 317, 749 [NASA ADS] [Google Scholar]
  79. Ud-Doula, A., Owocki, S. P., & Townsend, R. H. D. 2009, MNRAS, 392, 1022 [Google Scholar]
  80. Wang, P., McWilliams, J. C., & Ménesguen, C. 2014, J. Fluid Mech., 755, 397 [Google Scholar]
  81. Yim, E., Billant, P., & Ménesguen, C. 2016, J. Fluid Mech., 801, 508 [Google Scholar]
  82. Zahn, J. 1983, in Saas-Fee Advanced Course 13: Astrophysical Processes in Upper Main Sequence Stars, eds. A. N. Cox, S. Vauclair, & J. P. Zahn, 253 [Google Scholar]
  83. Zahn, J.-P. 1992, A&A, 265, 115 [NASA ADS] [Google Scholar]
  84. Zahn, J.-P., Brun, A. S., & Mathis, S. 2007, A&A, 474, 145 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  85. Zeitlin, V. 2018, Phys. Fluids, 30, 061701 [Google Scholar]
  86. Zorec, J., & Royer, F. 2012, A&A, 537, A120 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

All Tables

Table 1.

Summary table for the inflectional-instability growth rate by its variation with parameters , Pe, and Re

Table 2.

Summary table for the variation of the inertial-instability growth rate and the inertially unstable regime

All Figures

thumbnail Fig. 1.

(a) Schematic of the radiative and convective zones, colored as yellow and orange, respectively, for the configuration of low-mass stars rotating with angular speed Ω0. (b) Horizontal shear flow U(y) in a local nontraditional f-plane at a colatitude θ in the radiative zone; zc denotes the transition altitude between the radiative and convective zones.

Open with DEXTER
In the text
thumbnail Fig. 2.

Eigenvalue spectra in the space (σi, σr) at Re = ∞,  N = 1, Pe = 1, f = 0.5, and (i.e., Ω = 1.031 and θ = 76°) for (a) (kx, kz) = (0.25, 0), (b) (kx, kz) = (0, 6), and (c) (kx, kz) = (0.25, 6). Triangles denote the maximum growth rates.

Open with DEXTER
In the text
thumbnail Fig. 3.

Eigenfunctions of the most unstable modes in Fig. 2 at Re = ∞,  N = 1, Pe = 1, f = 0.5, and : (a) the mode shape and (b) velocity v in the physical space (x, y) at z = 0 for (kx, kz) = (0.25, 0), and σ = 0.0712, (c) the mode shape and (d) velocity v in the physical space (y, z) at x = 0 for (kx, kz) = (0, 6), and σ = 0.250. In (a) and (c), black, blue, and red lines denote the absolute value, real part, and imaginary part of the mode shape , respectively.

Open with DEXTER
In the text
thumbnail Fig. 4.

Contours of the maximum growth rate max(σr) in the parameter space of (kx, kz) (solid lines) for Re = ∞ and N = 1 for (a) , (b) , and (c) . For (a) and (b), this corresponds to (Ω, θ) = (1, 90° ), and (Ω, θ) = (1.031, 76°) for (c). Black and blue solid lines are contours of the maximum growth rate of the inflectional and inertial instabilities, respectively, and gray dashed lines in the background are contours of the maximum growth rate for the same parameters but in the traditional f-plane approximation where .

Open with DEXTER
In the text
thumbnail Fig. 5.

Inviscid growth rate of the inflectional instability for various parameter sets of and N at f = 0 (i.e., θ = 90°), Pe = ∞, and kz = 0.

Open with DEXTER
In the text
thumbnail Fig. 6.

Inviscid growth rate of the inertial instability for various Péclet numbers at N = 1 and kx = 0 for f = 0 and , i.e., (Ω, θ) = (1, 90° ).

Open with DEXTER
In the text
thumbnail Fig. 7.

Function σt(y) (solid lines) at N = 1 for (i.e., Ω = 1.031 and θ = 76°). Gray and white areas denote the regimes where is positive and negative, respectively. The dashed line represents an example of the growth rate σ = 0.250.

Open with DEXTER
In the text
thumbnail Fig. 8.

(a) Contours of the maximum growth rate σ0 from Eq. (59) in the parameter space for N = 2. Black dashed lines represent the upper limit of the instability from Eq. (61) for different N, and the white dashed line represents fmax from (64) where the maximum growth rate is attained. (b) Contours of σ0 from Eq. (59) in the parameter space at f = 1. (c) Contours of the first-order term σ1 in the parameter space for N = 2 and the first branch at m = 1.

Open with DEXTER
In the text
thumbnail Fig. 9.

Numerical results (lines) and asymptotic results (circles) from (58) for various values of and N at f = 1, kx = 0, Re = ∞, and Pe = ∞.

Open with DEXTER
In the text
thumbnail Fig. 10.

Numerical results (lines) of the growth rate σ as a function of kz for various parameter sets of : (1,1) (black, Ω = 0.707 and θ = 45°), (1,2) (blue, Ω = 1.118 and θ = 63.4°), (1,3) (red, Ω = 1.581 and θ = 71.6°), (−0.5, 1) (green, Ω = 0.559 and θ = 116.6°), (−1, 1) (yellow, Ω = 0.707 and θ = 135°), at Pe = 0.1, kx = 0, N = 1, and Re = ∞, and asymptotic results (circles) from Eq. (79) in the limit Pe → 0.

Open with DEXTER
In the text
thumbnail Fig. 11.

Contours of σ0, 0 in the parameter space of .

Open with DEXTER
In the text
thumbnail Fig. 12.

Maximum growth rate σmax = σ0, 0 of Eq. (80) as a function of the absolute Coriolis parameter 2Ω at various colatitudes θ. Solid lines and dashed lines denote the results in the northern and southern hemispheres, respectively.

Open with DEXTER
In the text
thumbnail Fig. 13.

Maximum growth rate σmax of the inflectional instability as a function of the Brunt-Väisälä frequency N for different values of at Re = Pe = ∞.

Open with DEXTER
In the text
thumbnail Fig. 14.

Maximum of σmax of the inflectional instability over N as a function of for different values of Pe (solid lines) at kz = 0 and f = 0.

Open with DEXTER
In the text
thumbnail Fig. 15.

Growth rate σr of the inertial instability as a function of kz for values of Pe at N = 1 and kx = 0 for (i.e., Ω = 1 and θ = 90°; solid lines) and for (i.e., Ω = 0.707 and θ = 135°; dashed lines). Dotted lines indicate asymptotes σ0, 0 of Eq. (80) predicted from the WKBJ analysis in the limit Pe → 0.

Open with DEXTER
In the text
thumbnail Fig. 16.

(a) Growth rate σ of the inflectional instability as a function of the streamwise wavenumber kx for kz = 0, N = 1, Pr = 1, , f = 0 (i.e., Ω = 1 and θ = 90°) for Re = ∞ (black), Re = 50 (blue), and Re = 10 (red). The gray dash-dot line denotes the inviscid growth rate of the inflectional instability under the traditional approximation (i.e., ). (b,c) Growth rate σ of the inertial instability as a function of the vertical wavenumber kz for kx = 0, N = 1, , and f = 0.5 (i.e., Ω = 2.02 and θ = 82.9°) for (b) Pr = 1 and Re = ∞ (black), Re = 1000 (blue), and Re = 500 (red), and (c) Pr = 10−6 and Re = 104 (black), Re = 5000 (blue), and Re = 1000 (red). Solid lines are numerical results, and dashed lines in (b) and (c) are asymptotic predictions from Eqs. (86) and (89), respectively.

Open with DEXTER
In the text
thumbnail Fig. 17.

(a) Optimal wavenumbers kz, c from (87) (solid lines) and kz, c0 from (90) (dashed lines) as a function of Re at 2Ω = 1, N = 4, and kx = 0 for colatitudes θ: 30° (blue), 45° (red), and 60° (green). Filled circles denote the critical Reynolds numbers. (b,c) Critical Reynolds number Rec as a function of the absolute Coriolis parameter 2Ω for N = 4, and (b) Pr = 1, (c) Pr = 10−6 at different colatitudes θ: 0° (black), 30° (blue), 45° (red), and 60° (green) from numerical results (crosses) and asymptotic results of Rec (solid lines in b) from (88) and Rec, 0 (dashed lines in c) from (91).

Open with DEXTER
In the text
thumbnail Fig. 18.

Effective turbulent viscosity νh as a function of the horizontal wavenumber for f = 0.5, (i.e., Ω = 0.559 and θ = 63.4°), and kx = 0 for Pe = ∞ and N = 1 (solid line), and as Pe → 0 (dashed line).

Open with DEXTER
In the text
thumbnail Fig. 19.

(a,b) Contours of the maximum of the effective turbulent viscosity max(νh, h) in the parameter space (N/2Ω, θ) for Pe = ∞ and (a) N = 1, (b) N = 2. (c) Contours of the maximum of the effective turbulent viscosity max(νh, h) in the parameter space (2Ω, θ) as Pe → 0.

Open with DEXTER
In the text
thumbnail Fig. 20.

Latitudinally-averaged turbulent viscosities (black) and (blue) computed from contours in Fig. 19.

Open with DEXTER
In the text
thumbnail Fig. 21.

Contours of the maximum of the effective turbulent viscosity max(νh, h) induced by the inflectional instability in the parameter space of for f = 0 (i.e., on the equator at θ = 90°) and Pe = ∞.

Open with DEXTER
In the text
thumbnail Fig. 22.

Growth rate σθ as a function of the colatitude θ for various values of ϵ. Solid lines denote σθ with ϵ >  0 and dashed lines denote σθ with ϵ <  0.

Open with DEXTER
In the text
thumbnail Fig. 23.

Latitudinally-averaged growth rate as a function of |ϵ| for the solar-like case (ϵ >  0, solid line) and the antisolar case (ϵ <  0, dashed line).

Open with DEXTER
In the text

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

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

Initial download of the metrics may take a while.