Coupling between internal waves and shearinduced turbulence in stellar radiation zones: the critical layers
^{1} Laboratoire AIM, CEA/DSM – CNRS – Université Paris Diderot, IRFU/SAp Centre de Saclay, 91191 GifsurYvette, France
email: lucie.alvan@cea.fr; stephane.mathis@cea.fr
^{2} Geneva Observatory, University of Geneva, chemin des Maillettes 51, 1290 Sauverny, Switzerland
email: thibaut.decressin@unige.ch
Received: 31 January 2013
Accepted: 5 March 2013
Context. Internal gravity waves (IGW) are known as one of the candidates for explaining the angular velocity profile in the Sun and in solartype mainsequence and evolved stars due to their role in the transport of angular momentum. Our contribution deals with critical layers, which are defined as the locations where the Dopplershifted frequency of the wave approaches zero (i.e., they correspond to corotation resonances).
Aims. The IGW propagate through stably stratified radiative regions, where they extract or deposit angular momentum through two processes: radiative and viscous dampings and critical layers. Our goal is to obtain a complete picture of the effects of these processes.
Methods. First, we expose a mathematical resolution of the equation of propagation for IGW in adiabatic and nonadiabatic cases near critical layers. Then, the use of a dynamical stellar evolution code, which treats the secular transport of angular momentum, allows us to apply these results to the case of a solarlike star.
Results. The analysis reveals two cases depending on the value of the Richardson number at critical layers: a stable one, where IGW are attenuated as they pass through a critical level, and an unstable turbulent case, where they can be reflected/transmitted by the critical level with a coefficient larger than one. Such overreflection/transmission can have strong implications on our vision of angular momentum transport in stellar interiors.
Conclusions. This paper highlights the existence of two regimes defining the interaction between an IGW and a critical layer. An application exposes the effect of the first regime, showing a strengthening of the damping of the wave. Moreover, this work opens up new ways concerning the coupling between IGW and shear instabilities in stellar interiors.
Key words: hydrodynamics / waves / turbulence / stars: rotation / stars: evolution
© ESO, 2013
1. Introduction
Thanks to helio and asteroseismology, we are able to extract a huge amount of information about solar and stellar structures and compositions (e.g. TurckChièze & Couvidat 2011; Aerts et al. 2010) and their internal differential rotation profile (Garcia et al. 2007; Beck et al. 2012; Deheuvels et al. 2012). We know that internal rotation modifies the stellar structure since it generates flows, instabilities, and chemical elements mixing, which modify stars’ evolution, for example, their lifetime and their nucleosynthetic properties (e.g. Maeder 2009, and references therein). Moreover, in order to understand the obtained solar/stellar rotation profiles and the related rotational history, it is essential to develop a complete theory incorporating the different angular momentum transport mechanisms occuring in stellar interiors (e.g. Mathis 2010, and references therein). In this context, one legitimately asks the question about the origin of the internal rotation profiles. Four main processes are responsible, in different ways, for the secular angular momentum transport in radiative interiors. First, a largescale meridional circulation is driven by structural adjustments of stars, externally applied torques, and internal stresses (e.g. Zahn 1992; Mathis & Zahn 2004; Decressin et al. 2009). Next, rotation profiles may be subject to different hydrodynamical shear and baroclinic instabilities and turbulence (e.g. Knobloch & Spruit 1982; Talon & Zahn 1997; Maeder 2003; Mathis et al. 2004). Then, fossil magnetic fields, trapped during early phases of stellar evolution once radiation zones have been formed (Braithwaite & Spruit 2004; Duez & Mathis 2010), can transport angular momentum through large scale torques and Maxwell stresses (e.g. Gough & McIntyre 1998; Mathis & Zahn 2005; Garaud & Garaud 2008; Strugarek et al. 2011). Finally, internal gravity waves (IGW) that are excited at the convection/radiation boundaries constitute the fourth mechanism able to transport angular momentum over large distances in stellar radiation zones (e.g. Press 1981; Goldreich & Nicholson 1989; Schatzman 1993; Zahn et al. 1997; Talon & Charbonnel 2005; Mathis & de Brye 2012). We note that all these processes are not necessarily present at the same time everywhere in the HertzsprungRussel diagram. Moreover, they act on various characteristic timescales in stars of different masses and ages.
The object of this paper is the propagation of IGW and the way they interact with the shear (i.e. the differential rotation) of the surrounding fluid. They are common in the terrestrial atmosphere and oceans (Eckart 1961; Chapman & Lindzen 1970), which is why they are fairly well known in geophysics. We will use this advantage in the stellar case. We here draw attention to the mechanism whereby IGW exchange energy with the mean flow, independently of other dissipative processes such as thermal and viscous diffusion. Indeed, when the frequency of excited waves is of the same order as the angular velocity of the fluid (we present the accurate definition later), a phenomenon of resonance occurs, which affects the properties of both the wave and the shear of the surrounding fluid. This phenomenon is called a critical layer. Under the assumption of a perfect fluid (neither a heat conductor nor viscous), Booker & Bretherton (1967) and Lindzen & Barker (1985) have provided first results about critical layers in the geophysical case. They have shown that, depending on the value of the Richardson number of the fluid, which compares the relative strength of the shear and the stable stratification (see Eq. (3)), the waves might be either attenuated or reflected by the critical layer. This reflection may even be an overreflection together with an overtransmission (see also Sutherland & Yewchuk (2004) for a laboratory evidence). At the same time, Koppel (1964), Hazel (1967), Baldwin & Roberts (1970), Van Duin & Kelder (1986) have completed this work by taking into account the conduction of heat and the viscosity of the fluid. Surprisingly, their conclusions about the role of critical layers are identical. In the astrophysical case, Rogava et al. (2007) have studied critical layers for gravitoalfven waves. However, all these authors have produced their studies in Cartesian coordinates, assuming that the stiffness of the domain was thin enough to neglect its curvature. In the case of stellar radiation zones where critical layers may play an important role (e.g. Barker & Ogilvie 2010; Barker 2011; Rogers et al. 2012), these equations should thus be generalized to the case of spherical coordinates to be able to treat deep spherical shells.
Therefore, after exposing our notations and assumptions (Sect. 2), we present a complete mathematical study of critical layers in the case of a perfect fluid (Sect. 3) and of a nonperfect fluid (Sect. 4) in spherical geometry. Then, we determine the related mean vertical flux of angular momentum transported by IGW (Sect. 5). Next, we implement our theoretical results in the dynamical stellar evolution code STAREVOL and apply our formalism to the evolution of a one solarmass star (Sect. 6). Finally, we present the conclusion and perspectives of this work (Sect. 7).
2. Definition, notations, and hypotheses
The star we consider is composed of (at least) a convective and a radiative region. IGW are excited at the boundary between these two regions. They propagate in the radiative zone and are evanescent in the convective zone (e.g. Press 1981). As explained in the Introductions, the study of IGW in geophysics is usually made in Cartesian coordinates, justified by a local approach. In the stellar case, however, a global approach using the spherical coordinates (r,θ,ϕ) is necessary.
In the frame of Zahn (1992), we choose a shelular angular velocity for the studied star’s radiation zone, , considering that, because of the strong stable stratification, the shear instability decreases the horizontal gradient of the angular velocity (e.g. Talon & Zahn 1997; Maeder 2003; Mathis et al. 2004). This consequently can be considered as only dependent on radius. For the moment, we neglect the action of Coriolis and centrifugal accelerations while the Doppler shift due to differential rotation is retained. We also neglect the action of a potential magnetic field.
Now, we need to introduce some quantities to describe the properties of the fluid in the studied radiative zone. Each scalar field X is written as (1)where we have introduced its horizontal average on an isobar, , and its associated fluctuation X′. The thermodynamic variables employed are the density ρ, the pressure p, the temperature T, and the specific entropy S. Next, the stratification is described in terms of the BruntVaïsälä frequency (2)where is the mean gravity in the Cowling approximation, where fluctuations of the gravific potential are neglected (see Cowling 1941), and the adiabatic exponent. The relative importance of the stable stratification restoring force and the shear destabilizing effects is quantified by the Richardson number (3)When Ri is small, the velocity shear overcomes the stabilizing buoyancy and turbulence and mixing occur (e.g. Talon & Zahn 1997). In contrast, when Ri is large, the fluid remains stable. Finally, for the case of a nonperfect fluid (Sect. 4), we introduce the viscosity ν of the fluid and the coefficient of thermal conductivity κ. These notations will be recaled at the proper moment.
IGW themselves are characterized by their relative frequency (4)where σ_{w} is their excitation frequency (from the base of the convective zone in lowmass stars or the top of the convective core in intermediate and highmass stars), m corresponds to a Fourier expansion along the longitudinal direction (Eq. (10)), and is the difference between the angular velocity at the level r and at the border with the convective zone. The introduction of m leads to the definition of two classes of waves. Prograde (respectively retrograde) waves correspond to negative (respectively positive) values of m.
Finally, all these notations allow us to define properly a critical layer. It arises for a wave whose relative frequency σ(r) becomes zero. The coherence with the qualitative definition given in the introduction is respected since σ(r_{c}) = 0 means that the excitation frequency of the wave equals − m times the angular velocity of the fluid. It corresponds to a corotation resonance. It is important to highlight that the position of the critical layer depends on both wave characteristics (m and σ_{w}) and shear properties (ΔΩ(r)). In our approximation, σ(r_{c}) = 0 means that is constant, so the critical levels are isobaric surfaces of the star.
3. Case of the perfect fluid
In this first approach, we treat the hydrodynamic equations assuming that the fluid is neither viscous nor a heat conductor.
3.1. Equation of propagation of IGW near a critical layer
Our aim is to calculate the Eulerian velocity field of a fluid particle. We introduce the time t and the unit vectors (ê_{r},ê_{θ},ê_{ϕ}), which are associated with the classical spherical coordinates (r, θ,ϕ). The velocity field is (5)where u is the velocity associated with the wave and the azimuthal velocity field of the differential rotation. We note that we have neglected all other largescale velocities such as meridionnal circulation. Then, we give the linearized equations of hydrodynamics governing IGW dynamics in an inertial frame, i.e. the momentum, continuity, and energy equations: (6)where D_{t} = ∂_{t} + ΔΩ∂_{ϕ}. The equation of energy conservation is obtained using the linearized equation of state (7)and assuming the ideal gaz law (8)where ℛ is the gaz constant and c_{P} the specific heatperunit mass at constant pressure.
We define the Lagrangien displacement ξ as u = D_{t}ξ and following Rieutord (1986), we expand it on the vectorial spherical harmonics basis () defined by (9)where are the spherical harmonics with l ∈ N and m ∈ ˜ − l,l¨. Thus, (10)Then, we decompose ρ′ and p′ using spherical harmonics: and we obtain a new system made up of radial equations of momentum (13)of mass conservation (14)and of energy in the adiabatic limit (15)The combination of these relations leads to the system presented by Press (1981): (16)where .
In stellar radiative regions, the transport of angular momentum is dominated by lowfrequency IGW with σ ≪ N, where N is the BruntVaïsälä frequency defined in Eq. (2). It allows us to apply the anelastic approximation (Press 1981), where acoustic waves are filtered out. Introducing (17)we obtain the following equation of propagation: (18)As the righthand side of Eq. (18) is of order with H_{P} the characteristic pressure or density height scale, it can be neglected if (19)which is the case here. Finally, we obtain the equation of propagation of IGW in a perfect fluid: (20)with (21)In Cartesian coordinates, this equation is called the TaylorGoldsteinSynge (TGS) equation. We observe that the value r = r_{c}, where σ(r_{c}) = 0, is a singular point for this equation. Thus, we now focus on the study of the behaviour of this equation around such a critical point. Then for r ≈ r_{c}, the equation of propagation becomes (22)where (23)is the value of the Richardson number at the critical level, and (24)the horizontal wavenumber at the critical layer. We saw in the introduction (Eq. (3)) that the Richardson number is relevant to distinguish between the relative importance of shear and stratification. We describe this divergence of cases from a quantitave point of view in the following part.
3.2. Mathematical resolution with the method of Frobenius
The Frobenius method offers an infinite series solution for a secondorder ordinary differential equation of the form in the vicinity of the singular point z_{0} = 0. For more details, a mathematical description of this method can be found in Teschl (2011).
In our case, p(r) = 0, , and the theory brings out two cases depending on the value of l, m, and Ri_{c}:

Case 1: The solution is (25)where A_{1} and B_{1} are constants representing the wave amplitude.

Case 2: We define the complex parameter (26)and the solution is given by (27)where A_{2} and B_{2} are constants.
We exclude this special value of Ri_{c} and choose to consider only two different cases for the rest of the paper: and . The results can then be extended to the case without recourse to the logarithmic solution (Van Duin & Kelder 1986).
We now discuss the hydrodynamical behaviour corresponding to the situations where and Applying the classical method exposed in Drazin & Reid (2004) to the generalized spherical TGS equation (Eq. (22)), we can identify that the first regime corresponds to the case where the fluid stays stable with respect to the vertical shear instability at the critical layer. In the other case, this instability and thus turbulence can develop. This clearly shows how it is necessary to go beyond the currently used formalisms for stellar evolution, where IGW and vertical shear instability are considered as uncoupled. Indeed, if a fluid becomes shear unstable, mixing occurs that modifies the local stratification and thus IGW propagation (see e.g. Brown & Sutherland 2007; Nault & Sutherland 2007). Thus, we now distinguish the case of critical layers when the fluid is stable from the one when it is unstable.
3.3. The case of stable critical layers: Ri
It is now time to understand the physical behaviour of the solutions given above. First of all, we propose to cut the physical domain into two parts: above and below the critical layer. As , η_{l,m} defined by (26)is a purely imaginary number. In order to clearly distinguish between real and imaginary parts, we introduce (28)The two solutions can be written as (29)where (resp. ) is available when r > r_{c} (resp. r < r_{c}). Booker & Bretherton (1967) and Ringot (1998) both explain a way to connect these solutions, considering that the term (r − r_{c})^{1/2 + iαl,m} can be compared to an upward propagative wave of the form e^{ikr} and respectively that (r − r_{c})^{1/2 − iαl,m} can be compared to a downward propagative wave of the form e^{−ikr}. We obtain the following identification: (30)In order to connect the solutions, we observe the behaviour of (r − r_{c}) above and below the critical level. As r − r_{c} decreases from positive to negative values, its complex argument changes continuously from 0 to − π (Ringot 1998). Mathematically, we get It follows that the solutions above and below the critical layer can be written as (33)Physically, these equations can be explained this way. Starting from above the critical layer, the downward propagating wave passes through the critical layer and is attenuated by a factor equal to e^{ − αl,mπ}. At the same time, starting from below the critical layer, the upward propagating wave is attenuated by the same factor and its amplitude becomes equal to A. We also underline that both waves take a phase difference when they cross the critical layer. In Fig. 1, we represent the attenuation rate of different waves defined by the numbers l and m passing through a critical layer. The greater the ratio , the stronger is the attenuation Att for the same value of Ri_{c}. The axis scale depends on l and m because the condition of validity of this result is . We so deduce that waves of high ratio (which do not necessarily correspond to high order) are strongly attenuated if they reach their critical layer.
Fig. 1 Attenuation rate Att = of the wave passing through a critical layer as a function of the Richardson number. We observe that Att increases with 

Open with DEXTER 
We have not yet discussed a latter point: the choice of the method of resolution. Most of the publications concerning IGW use another process to solve the equation of propagation (Press 1981; Zahn et al. 1997; Mathis 2009). In fact, the WentzelKramersBrillouinJeffreys (WKBJ) theory is particularly adapted to the resolution of this equation. However, it is not convenient in our case because it imposes a condition on the value of the Richardson number as demonstrated in appendix A. It is shown that the WKBJ approximation is available only if . Despite this restriction, we write the solution: by separating the domain into two parts, we obtain (34)As (35)it becomes: (36)It is comforting to see that both methods (Frobenius and WKBJ) give the same solution to a multiplicative constant when the Richardson number becomes high: (37)
Fig. 2 Assumed neighbourhood of an unstable critical layer for the calculation of the IGW’s reflection and transmission coefficients. We assume that the unstable region around the critical layer (in red) has a thickness given by 2δ (zone II). The surrounding regions where IGW are propagative are in beige (zones I and III). 

Open with DEXTER 
3.4. The unstable case: Ri
In the unstable regime, the Frobenius method also gives a solution but we are not able to identify upward and downward propagating waves because of shearinduced instability and turbulence. As a consequence, we can not connect the solutions at the critical layer. In order to avoid this difficulty, we here propose to solve Eq. (22) following the method developped by Lindzen & Barker (1985). They applied it in the case of Cartesian coordinates, and our work is to generalize it to spherical coordinates. The parameter η_{l,m} defined in Eq. (26) is real in this case, and we introduce (38)Equation (22) becomes (39)We are seeking solutions of the form , where Φ_{l,m} is the solution to the Bessel equation (40)Consequently, Φ_{l,m} is a combination of the Bessel’s modified functions I_{ηl,m}(X) and I_{ − ηl,m}(X): (41)The final solution is given by: (42)We would like to calculate the reflection and transmission coefficients of the wave passing through the unstable critical layer. We assume that the fluid has the profile described in Fig. 2. We decompose the studied region into three zones defined by the value of the quantity . In zones I and III, the Richardson number is high enough to allow us to apply the WKBJ method described in the previous section. However, around the studied critical layer, in unstable zone II, we must use the solution with the modified Bessel functions. Moreover, we assume here that the thickness 2δ of this zone is small in comparison with the characteristic length of the problem. Then, we can consider that the wavenumber k_{V} is constant as (43)We consider a wave coming from the overside of the critical layer (zone I). It is partly transmitted toward zone III and partly reflected backward to zone I. We thus write the solutions corresponding to the three zones: (44)where Ψ_{I}(r) is available if (r − r_{c}) > δ, Ψ_{II}(r) if − δ < (r − r_{c}) < δ and Ψ_{III}(r) if (r − r_{c}) > − δ. The coefficients A, B, T, and R are calculated by to the following four continuity equations:(45)They correspond to the continuity of the solution and of its first derivative, which physically means that both displacement and mechanical stresses are continuous. After some algebra, we obtain the coefficients R and T of reflection and transmission, which depend on the stiffness δ, the vertical and horizontal wavenumbers k_{Vc} and k_{Hc} (see Eqs. (43) and (24)) , and the variable η_{l,m} (Eq. (26)). In order to simplify the formula, we note I_{ ± ηl,m} instead of I_{ ± ηl,m}(k_{Hc}δ). Then, we obtain (46)with and (47)with We have now calculated the transmission and reflection coefficients of an IGW, through an unstable region around a given critical layer, where . We represent the level lines of R and T in Fig. 3 for l = 4 and m = ± 3. They are plotted as a function of the Richardson number at the critical layer Ri_{c}, growing from 0 to its maximum value defined by . For instance, in Fig. 3, . The other variable is the halfthickness δ of the unstable layer (zone II), arbitrarily chosen, that points the nonlocal character of unstable turbulent layers. The main result is that both coefficients are greater than 1 when is small enough. Consequently, for a low Richardson number at the critical layer, the wave can be overreflected and overtransmitted at the same time. This means that, in contrast to the first stable case, the wave takes potential energy from the unstable fluid and converts it into kinetic energy. In other words, the turbulent layer acts as an excitation region. If R < 1 and T < 1, we speak about IGW tunneling (Sutherland & Yewchuk 2004; Brown & Sutherland 2007; Nault & Sutherland 2007). Another remark concerns the dependency of R and T on m. We note that the sign of m does not matter since only its square appears in the expressions. Physically, it shows that the critical layer’s action is the same on prograde and retrograde waves. This point is of importance because we know that other dissipative processes occuring during the propagation of IGW discriminate between both types of waves. In order to visualize the action of the critical layer on different waves, Fig. 2 shows the level lines R = 1 for 1 ≤ l ≤ 5 and 1 ≤ m ≤ l. We previously said that it is useless to consider negative values of m since R depends only on m^{2}. The pair (l,m) is indicated on each line, followed by the value of . Lines in the same colour correspond to the same value of l. These lines mark out the limit to observe overreflection for a chosen wave. We observe that the higher the value of , the stronger is the condition on Ri_{c} to observe an overreflection.
Fig. 3 Level lines of reflection R and transmission T coefficients of an IGW at a critical layer as a function of the Richardson number Ri_{c} and the thickness δ of the critical layer (zone II in Fig. 2). The top panels represent the level lines for an arbitrarily chosen value of (l,m) = (4, ± 3), while the bottom panels show level lines R = 1 and T = 1 for different couples (l,m). 

Open with DEXTER 
3.5. Choice of the method
We decided to apply different methods to solve the stable and unstable cases. However, it could be legitimate to wonder if both methods are equivalent from a mathematical point of view. In this section, we present a short comparison between the solutions obtained with the method of Frobenius and the one with Bessel functions. The modified Bessel function I_{ηl,m}(X) can be computed using (48)At the neighbourhood of the critical layer, X tends to 0 and the firstorder expression is (49)Close to the critical layer, the global solution given in Eq. (42) is then (50)Now, we recall the expression of the solution given by the method of Frobenius (Eq. (27)), rewritten with the previous notations (51)In conclusion, for a fixed couple (l,m), and vary in the same way as function of r − r_{c}.
4. Case of the nonperfect fluid
Until now, we have studied the role of the critical layers assuming that the fluid was perfect. In order to make the problem more realistic, we include in the second part the viscosity ν of the fluid and the coefficient of thermal conductivity κ (e.g. Koppel 1964; Hazel 1967; Baldwin & Roberts 1970; Van Duin & Kelder 1986).
4.1. Equation of propagation of IGW near a critical layer
The linearized equations of hydrodynamics in Eq. (6) become (52)The following method for forming the propagation equation is adapted from the work of Baldwin & Roberts (1970). We assume that the mean density of the fluid takes nearly a constant value, that is to say that is small compared with the characteristic lengths of the problem. First, we project Eq. (52) onto ê_{r} and apply the operator ∇^{2}. Then we apply D_{t} − κ∇^{2}. After combination with the two other equations, we obtain (53)As in the previous section, we decompose the radial velocity on the basis of spherical harmonics (54)where σ = σ_{w} + mΔΩ(r). Moreover, it is easier to work with dimensionless numbers. For this reason, we introduce the notations detailed in Table 1.
Then, Eq. (53) becomes for each pair (l, m) such as l ∈ N and m ∈ [[− l,l]]: (55)where Δ_{l} is the scalar spherical Laplacian operator: (56)Lastly, we introduce and due to a development close to the critical layer obtain the sixthorder equation: (57)
where now . Equation (57) can be compared with the one obtained by Press (1981), who has neglected the viscosity ν. Moreover, we can see that if we take ν = κ = 0, and without forgetting that P_{r} also depends on ν, Eq. (57) is identical to Eq. (22) for the perfect fluid.
4.2. Mathematical resolution
The resolution of Eq. (57) requests several substitutions and is quite complex. The detailed calculation can be found in Appendix B, and we give only the main steps here. We drew our inspiration from Hazel (1967), Baldwin & Roberts (1970), Koppel (1964), Van Duin & Kelder (1986), who solved the same equation in Cartesian coordinates. However, the resolution in spherical coordinates has some differences. The aim is to rewrite the equation here under a known form: the Whittaker differential equation; then, after some algebra, we can show that Eq. (57) can be written in the following form: (58)The solutions of Eq. (58) are thus the Whittaker functions (Abramowitz & Stegun 1965): (59)where is the confluent hypergeometric function of Kummer (60)with and (k)_{0} = 1, being the usual gamma function, (61)and (62)Because of this solution, we obtain the expression of the radial displacement of the wave: (63)We will clarify the function f(κ,P_{r}) of the thermal diffusivity coefficient later.
There is still the last step to carry out. We must find a curve along which we integrate the solution, that is to say, we determine a and b. During the calculation detailed in the Appendix B, Eq. (B.3), (64)has not been used yet. A sufficient condition to make this equality true is . Considering that t^{3} = t^{3}e^{3iθt} (θ_{t} being the complex argument of t), we have (65)There are several possibilities for the choice of curves. Koppel (1964) proposed building a basis of six solutions using the curves (,,) represented in Fig. 6 (left). Therefore, the solutions of the sixthorder equation (Eq. (57)) are linear combinations of (66)and V_{i}(η) with i ∈ { 1,2,3 } , where V_{i}(η) corresponds to U_{i}(η) with the opposite sign for M_{l,m}.
Fig. 4 Left: curves (,,) defining the basis (U_{1},U_{2},U_{3},V_{1},V_{2},V_{3}). Right: curves (,,) defining the basis (u_{1},u_{2},u_{3},v_{1},v_{2},v_{3}). 

Open with DEXTER 
4.3. Application to IGW
It is now time to apply these mathematical results to IGW. We recall the solution for the perfect fluid obtained in Sect. 3. For the moment, it is not necessary to distinguish between the stable regime, i.e. , and the unstable one, i.e. . The Frobenius solution is available everywhere if we take as a complex number. Remembering that , the radial Lagrangian displacement is (67)For easier reading, we have removed the indices r;l,m. Subscript P designates the solution for the perfect fluid.
Concerning the solution for nonperfect fluids, Hazel (1967) proposed using another basis for the solutions of Eq. (57). The curves (, , ) are represented in Fig. 6. The new basis (u_{1},u_{2},u_{3},v_{1},v_{2},v_{3}) is defined by

u_{1} = U_{1} + U_{2} + U_{3} and v_{1} = V_{1} + V_{2} + V_{3},

and ,

and ,
where X^{∗} is the complex conjugate of X. Consequently, the solution for the nonperfect fluid (subscript NP) can be written as (68)A relation between both solutions and exists if we consider that the Reynolds number is great. In the case of the Sun (e.g. Brun & Zahn 2006), the microscopic viscosity in the radiative zone is weak. For this reason, it is appropriate to consider that the Reynolds number is much greater than 1. This assumption leads to the relation (69)Baldwin & Roberts (1970) give tables for the asymptotic behaviour of u_{i} and v_{i} (i ∈ 1,2,3). The solutions u_{1}, u_{2}, v_{1} and v_{2} diverge when R_{e} → + ∞. They are therefore physically unacceptable and we deduce that (70)Tables 2 and 3 give the asymptotic expressions of u_{3} and v_{3} as a function of the sign of r − r_{c}.
Expressions of u_{3} and v_{3} when R_{e} → + ∞, below the critical layer.
Expressions of u_{3} and v_{3} when R_{e} → + ∞, above the critical layer.
In the stable case where , is a real number. The same findings as in the previous section can be made: after the passage through the critical layer, the wave is attenuated by a factor e^{ − αl,mπ}. In the unstable case, where , is a complex number, so we cannot interpret the solution as upward and downward propagating waves. But the expressions given in Tables 2 and 3 remain comparable to those in the first section, and we deduce that the calculation of R and T will lead to the same result: the possibility of overreflection and overtransmission.
4.4. Radiative and viscous dampings
4.4.1. General equations
Summary of the different cases.
We volontarily left aside the factor f(κ,P_{r}) in Eq. (66). In order to establish its expression, Zahn et al. (1997) used the equation of the propagation of IGW, taking into account heat diffusion but with a viscosity coefficient (ν) equal to zero. Here, we generalize their result for ν ≠ 0 (i.e. P_{r} ≠ 0) and obtain (71)where (72)We introduced the general expression for N^{2}, the BrüntVäisälä frequency, to be able to take chemical stratification into account. Then, we have (73)with and where is the pressure heighscale, the temperature gradient, and the mean molecular weight (μ) gradient. Moreover, we introduced the generalized equation of state (EOS) given in Kippenhahn & Weigert (1990): (74)where and . Next, σ = σ_{w} + mΔΩ is the Dopplershifted frequency of the wave relative to the fluid rotation with an excitation frequency σ_{w}, r_{c} and r_{CZ} are respectively the positions of the critical layer and the boundary between the studied radiative zone and the convection region where IGW are initially excited. This damping is another source of attenuation independent of the presence of a critical layer. Moreover, as shown in Sects. 3.4. and 4., we have to consider IGW reflected and transmitted by unstable critical layers in addition to those initially excited by convection. Then, we introduce the notation (75)
where r_{1} and r_{2} are respectively the emission point of the IGW and the position with r_{1} > r_{2} (in the opposite case, where r_{1} < r_{2}, limits in the integral have to be reversed). This will enable us to describe radiative and viscous dampings in every configuration. We note that in stellar radiation zones P_{r} ≪ 1 (Brun & Zahn 2006) inducing that the damping is mostly radiative. We now compare it with the effects of critical layers.
4.4.2. Prograde and retrograde waves
For a same environment, prograde waves have a frequency lower than retrograde ones (e.g. Eq. (4)). We chose, for example, a couple of IGW with the same excitation frequency, σ_{0}, the same number l and opposite azimuthal degrees m. Thus, we compared a prograde wave of frequency σ_{1}(r) = σ_{0} − m ΔΩ and a retrograde one of frequency σ_{2}(r) = σ_{0} + m ΔΩ. We obtained σ_{1}(r) − σ_{2}(r) = − 2 m (Ω(r) − Ω_{CZ}) < 0 in the presence of negative gradient of Ω and σ_{1}(r) − σ_{2}(r) > 0 if the gradient is positive. We stress that τ given in Eq. (71) varies globally as . As a consequence, assuming that a negative Ωgradient is present near the excitation layer (see Sect. 6), the radiative damping is stronger for prograde waves than retrograde waves. Therefore, prograde IGW are absorbed by the fluid much closer to their region of excitation, while the retrograde waves are damped in a deeper region. This process is responsible for the net transport of angular momentum by IGW in stars. On the other hand, critical layers do not introduce such net bias between prograde and retrgrade IGW because their effects scale with m^{2}.
4.4.3. Dependency in l and m
The second comment concerns the variation of radiative and viscous dampings as a function of l and m. On one hand, looking at the multiplicative factor, which is in front of the integral in Eq. (72), we can roughly write that . On the other hand, the expression of the attenuation factor due to stable critical layers is proportional to , which is always greater than unity since m ≤ l. The comparison between radiative and viscous dampings and the one due to stable critical layers for an IGW with given will be examined in Sect. 6.2.3.
4.4.4. Location
Lastly, radiative and viscous dampings occur throughout the whole propagation of IGW. In contrast, critical layers are localised. Moreover, there is a condition for a wave to reach a critical layer: the rotation speed of the fluid must be of the same order as the wave frequency for a chance of observing σ = 0. Finally, all IGW are concerned by radiative and viscous dampings, which increase around critical layers since τ ∝ σ^{4}.
In Table 4, we sum up the different cases studied in these work. Depending both on the properties of the fluid and the studied wave, the wave is submitted to different phenomena.
5. Transport of angular momentum
As stated in the introduction, our goal is to study the transport of angular momentum in stellar radiation zones and to unravel the role of critical layers. Therefore, the first step is to recall the flux of angular momentum transported by propagative IGW and by the shearinduced turbulence. To illustrate our purpose, we focus here on the case of a lowmass star, where IGW are initially excited at the border of the upper convective envelope (see Fig. 5) by turbulent convection (Garcia Lopez & Spruit 1991; Dintrans et al. 2005; Rogers & Glatzmaier 2005; Belkacem et al. 2009; Brun et al. 2011; Lecoanet & Quataert 2013) and by tides if there is a close companion (Zahn 1975; Ogilvie & Lin 2007) with an amplitude A (corresponding results for massive stars with an internal convective core can be easily deduced by reversing signs).
Fig. 5 Two configurations in a lowmass star with an external convective envelope. Left: the case of a stable critical layer (CL) where IGW are damped. Right: the case of an unstable critical layer where IGW can be overreflected/transmitted. 

Open with DEXTER 
5.1. Angular momentum fluxes
5.1.1. Angular momentum flux transported by propagative IGW
First, we have to calculate the angular momentum flux transported by a propagative monochromatic wave over a spherical surface. It is given by the horizontal average of the Reynolds stresses associated with the wave (e.g. Zahn et al. 1997): (76)where . Besides, Mathis (2009) shows that (77)where ℱ_{E;l,m,σ} is the horizontal average of the energy flux in the vertical direction expressed by Lighthill (1986) as (78)So finally, the angular momentum flux is given by (79)To calculate this angular momentum flux, we need the expressions of u_{r;l,m} and . The expression of u_{r;l,m} is immediately accessible from Eq. (33), solution of Eq. (22). It is a little more complicated for because we need to go back in the calculation leading to Eq. (22). The expression of results from the first equation of the system given in Eq. (16), which can be reduced to (80)applying the anelastic approximation and neglecting terms of order 1/L^{2}.
The total vertical flux of angular momentum transported by propagative IGW is then given by (81)and we define the associated socalled action of angular momentum as (82)
5.1.2. Angular momentum flux transported by shearinduced turbulence
In the case of shearunstable regions, such as region II in the regime where , IGW are unstable (e.g. Drazin & Reid 2004). After the nonlinear saturation of the instability, a steady turbulent state is reached. Then, as established, for example, in Zahn (1992); Talon & Zahn (1997), and now confirmed by numerical simulations (see Prat & Lignières 2013), the vertical flux of angular momentum transported by shearinduced turbulence is given by (83)where (84)The adopted value for the critical Richardson number is Ri^{c} = 1/6, and ν_{H} is the horizontal turbulent viscosity for which we assume the prescription derived by Zahn (1992).
5.1.3. Equation of transport of angular momentum
We now refocus these results in the wider frame of the complete angular momentum transport theory. Considering the other transport mechanisms we presented in the introduction, the angular momentum transport equation that takes into account meridional flows, shearinduced turbulence, and IGW (e.g. Talon & Charbonnel 2005; Mathis 2009) becomes (85)The first term on the righthand side corresponds to the advection of angular momentum by the meridional circulation, where U_{r} = U_{2} P_{2}(cosθ) is its radial component. The Lagrangian derivative takes into account the radial contractions and dilatations of the star during its evolution, which are described by ṙê_{r}. According to our hypothesis, we do not take into account the transport by the Lorentz force, associated with magnetic fields in stellar radiative zones (Mathis & Zahn 2005).
The major difference with previously published equations is that IGW and shearinduced turbulence transports of angular momentum are not summed linearly since they are, as we demonstrated before, intrinsically coupled. Therefore, for stable regions, one must consider the Reynolds stresses of IGW (Eq. (82)), while for unstable regions, only the vertical turbulent flux given in Eq. (90) must be taken into account.
5.2. Stable critical layer
5.2.1. Case of the perfect fluid
This is the simplest case of a stable critical layer in a perfect fluid. We apply Eq. (79) to obtain the expressions for the transported fluxes by propagative IGW below and above the critical layer: (86)where A is the initial amplitude of the IGW at r = r_{CZ}, (see Eq. (28)) and (87)where are the associated Legendre polynomials. Booker & Bretherton (1967) obtained similar results in Cartesian coordinates. We see the expected attenuation of the flux by a factor e^{ − 2παl,m} when the wave passes through the critical layer. Moreover, ℱ_{J;l,m,σ ± } depends on m (and not on m^{2}). As a consequence, we recover the classical result that prograde waves (m < 0) and retrograde ones (m > 0) have opposite angular momentum flux (respectively a deposit and an extraction). Finally, the monochromatic action of angular momentum ℒ_{J;l,m,σ}(r) = 4πr^{2}ℱ_{J;l,m,σ} is constant in each region because of the absence of dissipation.
5.2.2. Case of the nonperfect fluid
We saw in the previous section that the solution of the equation of propagation in the case of a nonperfect fluid is similar to the one obtained for a perfect fluid. For this reason, we are able to apply the same method for the calculation of ℱ_{J;l,m,σ}. We obtain (88)The difference to Eq. (86) comes from the introduction of radiative and viscous dampings (Eq. (75)).
The conclusion is that, in the stable case, the attenuation due to the passage through a critical layer is added to dampings due to dissipation. This observation leads us to simply implement the role of stable critical layers as an additional term in the damping coefficient.
5.3. Unstable critical layer
5.3.1. Region I: r_{c} + δ ≤ r ≤ r_{CZ}
Using the summary given in Fig. 5 for the unstable case and the obtained results for regions where IGW are propagative (Sect. 5.2.2), we obtain: (89)where we identify the transport induced by the incident wave, which propagates downward, and the one induced by the reflected one, which propagates upward. We see that the unstable critical layer is a second excitation source for IGW propagating in region I and that the angular momentum transport will be modified in such situation, particularly when overreflection (R > 1) occurs.
5.3.2. Region II: r_{c} − δ < r < r_{c} + δ
In this unstable region, we directly use results obtained in Sect. 5.1.2. to describe the flux of angular momentum transported by the shearinduced turbulence, i.e.: (90)
5.3.3. Region III:
Using the same methodology as for region I, we get in a straightfoward way (91)where we indentify the transport induced by the transmitted wave, which propagates downward. As in region I, we can see that the unstable critical layer constitutes a secondary excitation source for IGW propagating in region III and that the angular momentum transport will be modified, particularly when overtransmission (T > 1) occurs.
Since the general theoretical framework has been given, we now explore the possibility of the existence of the two different regimes (stable and unstable) along the evolution of a given star. As a first application, we chose to study the case of a solartype star which has already been studied without critical layers by Talon & Charbonnel (2005).
6. A first application: the evolution of a solartype star
6.1. The STAREVOL code
We used the onedimensional hydrodynamical Lagrangian stellar evolution code STAREVOL V3.10 (see Lagarde et al. (2012) and references therein for a detailed description of the input physics). We recall the main characteristics and parameters used for the modelling that are directly relevant for the present work. We use the Schwarzschild criterion to determine the position of convective zones, and compute their temperature gradient according to the mixing length theory with a α_{MLT} = 1.75. The solar composition is taken from Asplund et al. (2005) with the Ne abundance from (Cunha et al. 2006). We generated the opacity table for temperatures higher than 8000 K following Iglesias & Rogers (1996) by using their website^{1}. The opacity table at lower temperatures follows Ferguson et al. (2005)^{2}. The mass loss rate is determined following Reimers (1975) with a parameter η_{R} = 0.5. The increase of mass loss due to rotation is taken into account following Maeder & Meynet (2001). However due to the small mass loss and velocity of our models this effect remains weak.
In radiative regions, we follow Mathis & Zahn (2004) formalism for the transport of angular momentum and chemicals as well as the prescription from Talon & Zahn (1997) for the vertical turbulent transport (Eq. (84)). We assume that convective regions rotate as a solid body. The treatment of IGW follows Talon & Charbonnel (2005, 2008), with the difference that the volumetric excitation by Reynolds stresses in the bulk of convective zones (e.g. Goldreich & Kumar 1990; Goldreich et al. 1994; Belkacem et al. 2009) is consistently computed at each timestep as a function of their physical properties.
We start with an initial model of 1.0 M_{⊙} at solar metallicity, with an initial surface velocity of 70 km s^{1}. The rotation profile is initially flat. We add magnetic braking through the following law: with a constant K = 3 × 10^{30}. This value has been calibrated to reproduce the surface velocity determined in the Hyades by Gaige (1993).
6.2. The effects of critical layers
6.2.1. Location of critical layers
We studied theoretically the impact of the critical layer for a given IGW, thus assuming that there are some radii where the relation σ(r) = 0 is satisfied. As a consequence, the first question we may answer due to the simulation concerns the existence of such critical layers and their location in the radiative zone. Figure 6 shows that critical layers do exist in the radiative core of the solarlike star studied here. In the three panels, we superimpose the rotational velocity of the star’s interior as a function of the normalized radius with the position of potential critical layers, which are indicated by colored squares that correspond to positions where σ = σ_{w} + mΔΩ = 0. Each panel corresponds to a given value of the excitation frequency σ_{w}. As expected, the positions of the critical layers only depend on the azimuthal number m of the wave, and not on the degree l. Through this plot, we confirm that critical layers exist in the solarlike star. However, we did take into account the fact that not all waves can reach these positions.
Fig. 6 Position of the critical layers for retrograde IGW for three different frequencies, superimposed with the rotational velocity profile in the star. 

Open with DEXTER 
First of all, we only consider retrograde waves because the prograde ones are damped immediately after their excitation (see Sect. 4.4.2). Moreover, in the lefthand panel of Fig. 7, we represente in logarithmic scale the luminosity of a given IGW at the moment of its initial stochastic excitation by the turbulent convection as a function of its degrees l and m, following the spectrum adopted in Talon & Charbonnel (2005). The difference between the three plots is the value of the excitation frequency σ_{w}. Through this representation, we can see that the maximum of excitation lies in a domain where l and m are close and quite small. Moreover, it shows that for each excitation frequency, the amplitude of the excited wave depends on l and m and may be close to zero. As a consequence, some critical layers represented in Fig. 6 may belong to a nonexcited wave or to a wave completely damped at this depth. Fortunately, results obtained show that some waves really meet their critical layer.
Fig. 7 Left: luminosity of the waves at the location of their excitation as function of l and m. Right: ratio between τ_{CL} and τ_{tot} = τ_{CL} + τ_{rad} (see Eq. (92)). 

Open with DEXTER 
6.2.2. Interaction between waves and critical layers
Concerning the way IGW interact with the surrounding fluid, the theory predicts two possible regimes depending on the value of the Richardson number Ri_{c} at the critical level. In our simulation, it appears that for every detected critical layer, the relation , which correponds to the stable regime, is verified. As a consequence, we establish that the second unstable regime (with the associated possible tunneling or overreflection and transmission) does not occur for the solarlike star of our simulation; forthcoming studies may explore other types of stars at different evolutionary stages, to see if this regime can occur. Therefore, for the solartype star studied here, we only implement in STAREVOL the terms related to stable critical layers: each time a wave passes through a critical layer, it is damped with a coefficient (see Eq. (88)).
6.2.3. Effect of critical layers
Since all interactions between waves and critical layers are of the same kind in the star studied here, we can concentrate on the quantitative importance of their effect on the transport of angular momentum. We know that in the stable regime (Sect. 5.2.2.), the wave passing through its critical layer is damped by the factor given in Eq. (88) that is added to radiative damping (here P_{r} ≪ 1 and the viscous damping is thus negligible), which has already been taken into account in previous works (e.g. Talon & Charbonnel 2005). In the righthand panel of Fig. 7, we thus choose to represent the ratio between τ_{CL}, the rate of attenuation due to the passage of the wave through its stable critical layer, and the sum τ_{tot} of this rate and that of the radiative damping. The explicit formula is (92)In contrast to the lefthand panel, only waves are represented here which meet a critical layer. That is the reason why several white zones are seen. They correspond to the waves which have been attenuated before reaching the depth of their critical layer. In the case of low σ_{w} (upper panels in Fig. 7), the highdegree waves (l > 45) are simply not excited, as we can see on the right. In red zones, the role of critical layers is important in comparison with the radiative damping, while dark blue regions are those where the radiative damping dominates. Therefore, this figure shows that critical layers should be taken into account.
6.2.4. Evolution of the rotation profile
Fig. 8 Evolution of the rotation profile, where the role of the critical layers in the transport of angular momentum is taken into account. The curves are labelled according to the corresponding ages in Gyr. Parameters are indicated above where K is the braking constant. Right: comparison between the evolution of the rotation with (red) and without (black) taking into account the effect of the critical layers. 

Open with DEXTER 
We now concentrate on the evolution of the rotation profile, when momentum deposition due to the interaction between waves and critical layers is taken into account. The other transport mechanism (IGW’s radiative damping effects, shearinduced turbulence, and meridional circulation) have been previously implemented in the code (Talon & Charbonnel 2005).
In the centre of the star, the rotation velocity is lower and increases the radiative damping of retrograde waves (see Eq. (72)). This forms an angular momentum extraction front, which propagates from the core to the surface to damp the differential rotation. Three fronts are seen in the top panel of Fig. 8 and have the same form as those already obtained in Talon & Charbonnel (2005). To isolate the action of critical layers on the evolution of the rotation profile, we superimposed in the bottom panel of Fig. 8 the surface velocity as function of the evolution time with (black line) and without (red line) critical layer effects between 2.8 × 10^{8} and 3.5 × 10^{8} years. Both curves are nearly identical, which shows that despite their local action as represented in Fig. 6, IGW’s stable critical layers do not disturb the dynamical evolution of surface velocity in the case of the studied star. This can be easily understood since the radiative damping becomes mostly efficient around the positions of critical layers because of its dependance on σ^{4}. Moreover, if in Talon & Charbonnel (2005) the IGW’s action and the shearinduced turbulence have been added as uncoupled physical mechanisms, it has been demonstrated that in the upper region (r > r_{c}) where IGW are propagative, the coefficient ν_{V} is negligible, while in the inner one (r < r_{c}) the differential rotation has been damped, leading to the same result with a transport dominated by the meridional circulation.
This clearly indicates that unstable critical layers will lead to major modification of the transport of angular momentum in stellar interiors. A systematic exploration of different types of stars for different evolutionary stages will be undertaken in the near future. Moreover, in order to give quantitative information, it is necessary to improve the way waves are excited in this model. Future work will implement a prescription about penetrative convection processes. The best way to do this is to use a realistic numerical simulation of such a mechanism to obtain the excitation spectrum at the base of the convective zone (Rogers & Glatzmaier 2005; Brun et al. 2011; Alvan et al. 2012).
7. Conclusion
In this paper, we study in detail a new mechanism of interaction between IGW and the shear of the mean flow that occurs at corotation layers in stably stratified stellar radiation zones. Taking advantage of the work in the literature about atmospheric and oceanic fluids, we highlight the similarities with such stellar regions and propose an analytical approach adapted to the related case of deep spherical shells. The use of spherical coordinates brings differences in the equations and make their resolution more complicated, but the final results are comparable. We then demonstrate the intrinsic couplings between IGW and the shearinduced instabilities and turbulence that can thus not be added linearly, as done previously in stellar evolution literature. We highlight the existence of two regimes where the interactions between IGW and the shear at critical layers are strongly different:

In the first case, the fluid is stable and the amplitude of IGW isoverdamped by the critical layer compared to the classical case,where only radiative and viscous dampings are taken intoaccount.

In the second case, the fluid is unstable and turbulent and the critical layer acts as a secondary excitation region. Indeed, through overreflection/transmission (when R > 1 and T > 1), energy is taken from the unstable shear that increases the amplitude of an incident IGW. Moreover, even in the case of simple reflection and transmission where R < 1 and T < 1, this demonstrates the existence of IGW tunneling through unstable regions as identified by Sutherland & Yewchuk (2004) in laboratory experiments and by Brown & Sutherland (2007), Nault & Sutherland (2007) in geophysics.
Therefore, these mechanisms opens up a new field of investigations concerning angular momentum tranport processes by IGW in stellar interiors.
Indeed, even if according to our first evolutionary calculations with STAREVOL, only the first stable regime exists in solartype stars, we expect to find stars where the unstable regime and possible tunneling or overreflection/transmission take place. Moreover, while the formalism presented in this work is general, several uncertainties remain. First, concerning the stochastic excitation of IGW by convection, the model used in our evolutionary code certainly underestimates the wave flux since it considers only the volumetric excitation in the bulk of the convective envelope, while convective penetration should also be taken into account. This will influence the measured action of critical layers since it is proportional to the initial amplitude of the IGW. Then, only retrograde waves are simulated here, considering that prograde ones are immediately damped and do not penetrate deeply in the radiation zone. This should normally not affect our results because the critical layers we detected are located in the deep radiation zone, but formal equations take both types of waves into account.
The last point to bear in mind is that no latitudinal dependence on the angular velocity is considered here. We explained the reason for this choice in the introduction, but one must not forget this approximation. Finally, since our goal is to obtain a complete and coherent picture of the transport of angular momentum in stellar radiation zones for every stellar type or evolutionary stage, it will be important to extend this work to cases of gravitoinertial waves, where the action of Coriolis and centrifugal accelerations is considered (e.g. Lee & Saio 1997; Dintrans & Rieutord 2000; Mathis 2009; Ballot et al. 2010) and magnetogravitoinertial waves, when radiation zones are magnetized (e.g. Rudraiah & Venkatachalappa 1972; Kim & MacGregor 2003; MacGregor & Rogers 2011; Mathis & de Brye 2012).
Acknowledgments
The authors thank the referee for her/his comments that allowed them to improve the paper. This work was supported by the French Programme National de Physique Stellaire (PNPS) of CNRS/INSU, by the CNESSOHO/GOLF grant and asteroseismology support in CEASaclay, by the Campus Spatial of the University ParisDiderot, and by the TOUPIES project funded by the French National Agency for Research (ANR). S.M. and L.A. are grateful to Geneva Observatory, where part of this work has been achieved. T.D. acknowledges financial support from the Swiss National Science Foundation (FNS) and from ESFEurogenesis.
References
 Abramowitz, M., & Stegun, I. 1965, Handbook of mathematical functions [Google Scholar]
 Aerts, C., ChristensenDalsgaard, J., & Kurtz, D. W. 2010, Asteroseismology [Google Scholar]
 Alvan, L., Brun, A. S., & Mathis, S. 2012, in SF2A2012: Proc. of the Annual Meet. of the French Society of Astronomy and Astrophysics, ed. S. Boissier, 289 [Google Scholar]
 Asplund, M., Grevesse, N., & Sauval, A. J. 2005, in Cosmic Abundances as Records of Stellar Evolution and Nucleosynthesis, eds. T. G. Barnes III, & F. N. Bash, ASP Conf. Ser., 336, 25 [Google Scholar]
 Baldwin, P., & Roberts, P. H. 1970, Mathematika, 17, 102 [CrossRef] [Google Scholar]
 Ballot, J., Lignières, F., Reese, D. R., & Rieutord, M. 2010, A&A, 518, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Barker, A. J. 2011, MNRAS, 414, 1365 [NASA ADS] [CrossRef] [Google Scholar]
 Barker, A. J., & Ogilvie, G. I. 2010, MNRAS, 404, 1849 [NASA ADS] [Google Scholar]
 Beck, P. G., Montalban, J., Kallinger, T., et al. 2012, Nature, 481, 55 [NASA ADS] [CrossRef] [Google Scholar]
 Belkacem, K., Samadi, R., Goupil, M. J., et al. 2009, A&A, 494, 191 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Booker, J., & Bretherton, F. 1967, J. Fluid Mech., 27, 513 [NASA ADS] [CrossRef] [Google Scholar]
 Braithwaite, J., & Spruit, H. C. 2004, Nature, 431, 819 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Brown, G. L., & Sutherland, B. R. 2007, AtmosphereOcean, 45, 47 [CrossRef] [Google Scholar]
 Brun, A. S., & Zahn, J.P. 2006, A&A, 457, 665 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brun, A. S., Miesch, M. S., & Toomre, J. 2011, ApJ, 742, 79 [NASA ADS] [CrossRef] [Google Scholar]
 Chapman, S., & Lindzen, R. 1970, Atmospheric tides, Thermal and gravitational (London: Hollis & Carter) [Google Scholar]
 Cowling, T. G. 1941, MNRAS, 101, 367 [NASA ADS] [Google Scholar]
 Cunha, K., Hubeny, I., & Lanz, T. 2006, ApJ, 647, L143 [NASA ADS] [CrossRef] [Google Scholar]
 Decressin, T., Mathis, S., Palacios, A., et al. 2009, A&A, 495, 271 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Deheuvels, S., Garcia, R., Chaplin, W. J., et al. 2012, ApJ, 756, 19 [NASA ADS] [CrossRef] [Google Scholar]
 Dintrans, B., & Rieutord, M. 2000, A&A, 354, 86 [NASA ADS] [Google Scholar]
 Dintrans, B., Brandenburg, A., Nordlund, Å., & Stein, R. F. 2005, A&A, 438, 365 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Drazin, P. G., & Reid, W. H. 2004, Hydrodynamic Stability (Cambridge Mathematical Library) [Google Scholar]
 Duez, V., & Mathis, S. 2010, A&A, 517, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Eckart, C. 1961, Phys. Fluids, 4, 791 [NASA ADS] [CrossRef] [Google Scholar]
 Ferguson, J. W., Alexander, D. R., Allard, F., et al. 2005, ApJ, 623, 585 [NASA ADS] [CrossRef] [Google Scholar]
 Gaige, Y. 1993, A&A, 269, 267 [NASA ADS] [Google Scholar]
 Garaud, P., & Garaud, J. D. 2008, MNRAS, 391, 1239 [NASA ADS] [CrossRef] [Google Scholar]
 Garcia, R. A., TurckChièze, S., JiménezReyes, S. J., et al. 2007, Science, 316, 1591 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Garcia Lopez, R. J., & Spruit, H. C. 1991, ApJ, 377, 268 [NASA ADS] [CrossRef] [Google Scholar]
 Goldreich, P., & Kumar, P. 1990, ApJ, 363, 694 [NASA ADS] [CrossRef] [Google Scholar]
 Goldreich, P., & Nicholson, P. D. 1989, ApJ, 342, 1079 [NASA ADS] [CrossRef] [Google Scholar]
 Goldreich, P., Murray, N., & Kumar, P. 1994, ApJ, 424, 466 [NASA ADS] [CrossRef] [Google Scholar]
 Gough, D., & McIntyre, M. 1998, Nature, 394, 755 [NASA ADS] [CrossRef] [Google Scholar]
 Hazel, P. 1967, J. Fluid Mech., 30, 775 [NASA ADS] [CrossRef] [Google Scholar]
 Iglesias, C. A., & Rogers, F. J. 1996, ApJ, 464, 943 [NASA ADS] [CrossRef] [Google Scholar]
 Kim, E.j., & MacGregor, K. B. 2003, ApJ, 588, 645 [NASA ADS] [CrossRef] [Google Scholar]
 Kippenhahn, R., & Weigert, A. 1990, Sky Tele, 80, 504 [NASA ADS] [Google Scholar]
 Knobloch, E., & Spruit, H. C. 1982, A&A, 113, 261 [NASA ADS] [Google Scholar]
 Koppel, D. 1964, JMP, 5 [Google Scholar]
 Lagarde, N., Decressin, T., Charbonnel, C., et al. 2012, A&A, 543, A108 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lecoanet, D., & Quataert, E. 2013, MNRAS, 430, 2363 [NASA ADS] [CrossRef] [Google Scholar]
 Lee, U., & Saio, H. 1997, ApJ, 491, 839 [NASA ADS] [CrossRef] [Google Scholar]
 Lighthill, J. 1986, Provided by the SAO/NASA Astrophysics Data System [Google Scholar]
 Lindzen, R., & Barker, A. 1985, J. Fluid Mech., 151, 189 [NASA ADS] [CrossRef] [Google Scholar]
 MacGregor, K. B., & Rogers, T. M. 2011, Sol. Phys., 270, 417 [NASA ADS] [CrossRef] [Google Scholar]
 Maeder, A. 2003, A&A, 399, 263 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Maeder, A. 2009, Physics, Formation and Evolution of Rotating Stars (Steven Stahler, Physics Today) [Google Scholar]
 Maeder, A., & Meynet, G. 2001, A&A, 373, 555 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mathis, S. 2009, A&A, 506, 811 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mathis, S. 2010, Astron. Nachr., 331, 883 [NASA ADS] [CrossRef] [Google Scholar]
 Mathis, S., & de Brye, N. 2012, A&A, 540, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mathis, S., & Zahn, J.P. 2004, A&A, 425, 229 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mathis, S., & Zahn, J.P. 2005, A&A, 440, 653 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mathis, S., Palacios, A., & Zahn, J.P. 2004, A&A, 425, 243 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Nault, J. T., & Sutherland, B. R. 2007, Phys. Fluids, 19, 6601 [NASA ADS] [CrossRef] [Google Scholar]
 Ogilvie, G. I., & Lin, D. N. C. 2007, ApJ, 661, 1180 [NASA ADS] [CrossRef] [Google Scholar]
 Olver, F. W. J. 1974, Asymptotics and Special Functions, ed. A. Press (New York: Marcel Dekker, Inc.) [Google Scholar]
 Prat, V., & Lignières, F. 2013, A&A, 551, L3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Press, W. 1981, ApJ, 245, 286 [NASA ADS] [CrossRef] [Google Scholar]
 Reimers, D. 1975, Problems in Stellar Atmospheres and Envelopes, eds. B. Baschek, W. H. Kegel, & G. Traving (Berlin: Springer), 229 [Google Scholar]
 Rieutord, M. 1986, Geophys. Astrophys. Fluid Dynamics, 39, 163 [NASA ADS] [CrossRef] [Google Scholar]
 Ringot, O. 1998, Ph.D. Thesis, École doctorale d’astronomie d’Ile de France [Google Scholar]
 Rogava, A., Gogoberidze, G., & Poedts, S. 2007, ApJ, 664, 1221 [NASA ADS] [CrossRef] [Google Scholar]
 Rogers, T. M., & Glatzmaier, G. A. 2005, MNRAS, 364, 1135 [NASA ADS] [CrossRef] [Google Scholar]
 Rogers, T. M., Lin, D. N. C., & Lau, H. H. B. 2012, ApJ, 758, L6 [NASA ADS] [CrossRef] [Google Scholar]
 Rudraiah, N., & Venkatachalappa, M. 1972, J. Fluid Mech., 52, 193 [NASA ADS] [CrossRef] [Google Scholar]
 Schatzman, E. 1993, A&A, 279, 431 [NASA ADS] [Google Scholar]
 Strugarek, A., Brun, A. S., & Zahn, J.P. 2011, A&A, 532, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Sutherland, B. R., & Yewchuk, K. 2004, J. Fluid Mech., 511, 125 [NASA ADS] [CrossRef] [Google Scholar]
 Talon, S., & Charbonnel, C. 2005, A&A, 440, 981 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Talon, S., & Charbonnel, C. 2008, A&A, 482, 597 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Talon, S., & Zahn, J.P. 1997, A&A, 317, 749 [Google Scholar]
 Teschl, G. 2011, Ordinary Differential Equations and Dynamical Systems (Providence, Rhode Island: American Mathematical Society) [Google Scholar]
 TurckChièze, S., & Couvidat, S. 2011, Rep. Progr. Phys., 74, 086901 [NASA ADS] [CrossRef] [Google Scholar]
 Van Duin, C., & Kelder, H. 1986, J. Fluid. Mech., 169, 293 [NASA ADS] [CrossRef] [Google Scholar]
 Zahn, J. 1992, A&A, 265, 115 [NASA ADS] [Google Scholar]
 Zahn, J.P. 1975, A&A, 41, 329 [NASA ADS] [Google Scholar]
 Zahn, J.P., Talon, S., & Matias, J. 1997, A&A, 322, 320 [Google Scholar]
Appendix A: Validity of the WKBJ approximation
The form of the equation to solve is (A.1)The first step is to introduce the Liouville transformation (e.g. Olver 1974), i.e., (A.2)We deduce and Eq. (A.1) becomes (A.5)where
In the present case we have (A.6)The WKBJ approximation is available when N^{2} ≫ σ^{2}. Consequently, we get
and (A.10)where Ri is the Richardson number of the fluid defined in Eq. (3).
Finally, we determine that the condition for applying the WKBJ approximation is Φ ≪ 1, which leads to .
Appendix B: Mathematical treatment for the nonperfect fluid case
We introduce (B.1)where a and b are the limits of a domain which will be defined later. The equation of propagation can then be written as (B.2)and (B.3)The new variable u, defined by , transforms the original equation into (B.4)Then, we introduce s = Dt^{3} where D ∈ C^{ ∗ }: (B.5)Finally, leads to (B.6)To obtain the final equation, we define and we get the following Whittaker equation: (B.10)
All Tables
All Figures
Fig. 1 Attenuation rate Att = of the wave passing through a critical layer as a function of the Richardson number. We observe that Att increases with 

Open with DEXTER  
In the text 
Fig. 2 Assumed neighbourhood of an unstable critical layer for the calculation of the IGW’s reflection and transmission coefficients. We assume that the unstable region around the critical layer (in red) has a thickness given by 2δ (zone II). The surrounding regions where IGW are propagative are in beige (zones I and III). 

Open with DEXTER  
In the text 
Fig. 3 Level lines of reflection R and transmission T coefficients of an IGW at a critical layer as a function of the Richardson number Ri_{c} and the thickness δ of the critical layer (zone II in Fig. 2). The top panels represent the level lines for an arbitrarily chosen value of (l,m) = (4, ± 3), while the bottom panels show level lines R = 1 and T = 1 for different couples (l,m). 

Open with DEXTER  
In the text 
Fig. 4 Left: curves (,,) defining the basis (U_{1},U_{2},U_{3},V_{1},V_{2},V_{3}). Right: curves (,,) defining the basis (u_{1},u_{2},u_{3},v_{1},v_{2},v_{3}). 

Open with DEXTER  
In the text 
Fig. 5 Two configurations in a lowmass star with an external convective envelope. Left: the case of a stable critical layer (CL) where IGW are damped. Right: the case of an unstable critical layer where IGW can be overreflected/transmitted. 

Open with DEXTER  
In the text 
Fig. 6 Position of the critical layers for retrograde IGW for three different frequencies, superimposed with the rotational velocity profile in the star. 

Open with DEXTER  
In the text 
Fig. 7 Left: luminosity of the waves at the location of their excitation as function of l and m. Right: ratio between τ_{CL} and τ_{tot} = τ_{CL} + τ_{rad} (see Eq. (92)). 

Open with DEXTER  
In the text 
Fig. 8 Evolution of the rotation profile, where the role of the critical layers in the transport of angular momentum is taken into account. The curves are labelled according to the corresponding ages in Gyr. Parameters are indicated above where K is the braking constant. Right: comparison between the evolution of the rotation with (red) and without (black) taking into account the effect of the critical layers. 

Open with DEXTER  
In the text 