Internal gravity waves in geophysical fluids

In this second Chapter, we recall some basic notions on stably-stratified flows before focusing on internal gravity wave dynamics. In Section 1, we illustrate the occurence of stably-stratified flows in nature and, in Section 2, we derive the Boussinesq approximation. Except for a very brief presentation of the Kelvin-Helmholtz instability in Section 3, we start discussing about internal gravity waves from this section on. The linear properties of the wave are discussed in this Section. The mechanisms that can lead the wave field to breaking are addressed in Sections 4 to 6: parametric and buoyancy-induced instabilities (Section 4), interaction with a shear flow (Section 5), interaction with a sloping boundary (Section 6). We briefly discuss about the statistical properties of the breaking wave field in Section 7 and introduce some notions on mixing in Section 8.


Stably-stratiied luids are ubiquitous
An incompressible luid is stably-stratiied if it displays a vertical density gradient with negative sign (or a vertical temperature gradient with positive sign): the luid becomes less dense, or warmer, as the altitude increases. Such luids are ubiquitous in nature: the oceanic water masses, the stratosphere (which is the part of the atmosphere comprised between:: : 10 and 50 kms) are stably-stratiied. This is exempliied in the oceanic case in Figure 1: the temperature proile measured in the Mindanao trench during the [1929][1930] Snellius expedition in eastern Indonesia is displayed as a function of depth. The gradient is positive down to 3500m but compressibility efects become important below this depth leading to a slight warming of the water masses. Such a change of sign in the temperature gradient therefore does not imply that the luid is unstably stratiied but simply results from compression of the water columns. The concept of potential temperature, which is discussed in detail in Chapter 4 of the lecture notes, is introduced to account for compressibility efects. Indeed, the vertical gradient of the potential temperature ield is now everywhere positive.
In this example, the huge depth of the temperature record, close to 10 km, along with the sparsity of the sampling, implies that the overall (large scale) behaviour of the temperature proile is captured. If one considers much smaller scales however, the vertical potential temperature gradient is not always positive. An example is provided in the stratosphere in Figure 2, which displays the potential temperature as a function of altitude. While the gradient is positive in overall, negative gradients exist locally, which are due to turbulent motions. In the stratosphere, the temperature and velocity luctuations are mostly due to internal gravity waves so that this turbulence results from the breaking of the waves, as we shall further discuss it .
Stably-stratiied media are also found in the interior of stars. In the Sun for instance, a stably stratiied region exists below the convective zone, in the so-called radiative zone. As we shall show it, internal gravity waves can propagate in a stably-stratiied medium. It has been proposed by Schatzman (1996) that the superposition of these waves results in a net transport of mass toward the interior of the Sun, which would account for the under-abundance of Lithium observed in this star: the Lithium would be transported  High resolution (20cm) temperature proile, just above the tropopause. Boxes marked with letters (a ... g) contain examples of sheets. Two close-ups of sheets "b" and "e" are displayed with the three proiles measured on the main gondola. Ongoing mixing is clearly present in the sheet "b" (from Dalaudier et al., 1994). text. The solar pond ( Figure 3) is a nice illustration. The idea is to capture heat by a stable stratiication. The principle is very simple: the bottom of a pond is made dark and very salty, in order for the bottom water, when heated by the sun (more precisely, when heated by the dark bottom) to remain heavier than the layer above. The temperature of the bottom water may reach 70 to 100 degrees and is used, for instance, to warm water inside a pipe that would cross the pond in its lowest part.
Stably-stratiied luids are also fo und in the cooling circuits of nuclear reactors.

Stably-stratiied luids give rise to phenomena with environmental im plications
The solar pond just discussed provides one example that stably-stratiied luids may give rise to phenomena with environmental implications. Another common phenomenon, which has a negative impact on the environment, is the well-known situation of a thermal inversion in the atmosphere of a valley in winter. The process is again simple (f.i. Stull, 1988): the ground cools when no longer heated by the sun and, when its temperature becomes lower than that of the neighbouring air, the air cools in turn in transfering heat to the ground. A very cold layer of air thus exists at the bottom of the valley, which is surmounted by warmer air, leading to a stable (potential) temperature gradient. This stable situation may persist for the whole day in winter. If cities are located at the bottom of the valley, the pollutants emitted by road traic and factories remain trapped in the stably-stratified layer of air, leading to high concentrations of pollutants. Two examples are provided in Figure 4, in the Grenoble valley and the Chamonix valley.

Quantiication of the stability: the Brunt-VaisaUi frequency
The Brunt-VaisaEi frequency is the frequency of a luid particle displaced adiabatically from its equilibrium position along the vertical. The demonstration is very standard. Let us do a thought experiment and consider one luid particle of unit volume, at equilibrium position z in a stably-stratiied luid with density proile p(z). When the particle is displaced adiabatically by a vertical distance (, the local density at the new position is p(z + (). The fluid particle is thus subjected to two forces, its weight p(z) and the Archimedean force -p(z + (). Applying Newton's law to the luid particle results in g � the equation ( + N 2 ( = 0, where N 2 =p(z) dz is the Brunt-Vaisala, or buoyancy, frequency. Stability implies that N 2 is positive.
A standard approximation to describe stably-stratiied incompressible luid motions is the Boussinesq approximation, described in Section 2.2; in this case, the density and the temperature proiles are linearly related. Therefore, the Brunt-Vaisala frequency may also be deined as N 2 = _Y dT . When the luid is compressible, the potential tempera-T(z) dz ture is used instead of the temperature and the expression of the Brunt-Vaisala frequency is the same as the previous one with the temperature profile being simply replaced by the potential temperature proile. As announced above, the generic luid motions in a stably-stratiied luid are internal gravity waves. Indeed, the sum of the weight and of the Archimedean force acting on a luid particle displaced fr om equilibrium position is (at irst order in() -?(z, which is a restoring fo rce. This is the restoring buoyancy force. A restoring force in a continuous medium generates wave motions. We shall show in Section 3 that linear waves are indeed a solution of the equations of motions. Before doing so, we provide two illustrations of internal gravity waves in the atmosphere ( Figure 5) and in the ocean ( Figure 6).
In the atmosphere, the main mechanims for gravity wave generation is the wind blowing over topography. The generated waves are called lee waves and are made visible by clouds. Figure 5 thus displays the cloud cover over the sea when the blowing wind interacts with two islands. The island located on the left of the picture is too high for the air advected by the wind to go over the island and the advected air simply lows around the island. In this case, the luid particles are not lifted and no wave is generated. Instead, a quasi-two-dimensional pattern (in a horizontal plane) of vortices are produced, of the Karmann street type. By contrast, the island located on the right of the picture is low enough for the advected air to go over the island and a quite diferent pattern is observed: in this case, waves are generated, with horizontal wave length of the order of the island size.
In the ocean, the two main mechanisms for the generation of internal waves are the wind and the tide. The wind forcing is indirect: the wind mixes the upper layer of the ocean (over a depth of order lOOm) and motions at the bottom of this layer perturbs the stably-stratiied luid below which generates internal gravity waves. By contrast, the tide generates internal gravity waves directly, through its interaction with topography (seamounts, continental slopes, etc.) Internal gravity waves generated in the neighbourhood of a continental slope by the tide are displayed in Figure 6. The waves are in the interior of the water masses but, because their amplitude is large and they propagate not far from the surface, the velocity ield they induce perturbs the fr ee surface which makes them visible using a radar on board of a satellite. Internal waves (propagating below the free surface) generated by the interac tion of the tide with a continental slope. This is a radar image acquired by the Synthetic Radar Aperture aboard of a satellite.

What about unstably stratiied lows?
Unstably stratiied lows are also very common in natural media. In this situation, no waves are generated since the buoyancy force is not a restoring one. On the opposite, this force lifts the luid particle further away from the position it has been brought. The commonest example occurs in the atmospheric boundary layer : the ground is heated by the Sun during the day, which heats the air (the direct heating of air by the Sun is much less important, especially if the air is dry) . The vertical gradient of the potential tem perature proile is therefore negative in this case. Such a situation is highly unstable and gives rise to large scale turbulent motions, which mix the luid. Therefore, if pollutants are emitted at the ground level, the local concentrations are generally smaller than in a stable situation, but may be found very far from their emission region. For instance, in summer in the Chamonix valley, high concentrations of heavy metals (such as Mercury) have been fo und at 4300m altitude, at Dome du Gouter (Veysseyre et al., 2000) . Let (x, y, z ) be a cartesian coordinate system in an absolute reference frame where the luid velocity is 1 = (u, v, w). The luid motions in this reference fr ame are governed by the compressible Navier-Stokes equations. These equations are composed of • the momentum equation, which is, when expressed per unit mass: (2.1) D where Dt : t + l. ' refers to the material derivative, p is the pressure ield, -2 0 \ 1 is the Coriolis acceleration due to earth rotation and l is the dynamical viscosity; • the mass conservation equation: where p is the density ield; • and a dynamical equation for the internal energy e, stemming from the irst law of thermodynamics applied to a system made of a unit mass of volume v = 11 p: Q is the amount of heat added to the system v and -p � is the work of pressure forces on the system.

The Boussinesq approximation
Water is a very weakly compressible luid so that oceanic luid motions may be as sumed to be incompressible over the typical vertical scale of the motions. The fact that density, that is the mass of the luid system, still varies in space and time may be in compatible with the assumption that the luid velocity is incompressible. A consistent, and very well-known, approximation that satisies the former requirement while assum ing that the luid is close to incompressiblity is the Boussinesq approximation. We now describe in detail this approximation (see also Cushman-Raisin, 1994, p. 37).
In the following, we set p = p0 + p where p0 is some constant reference density and p refers to the deviation ield about p0. In the ocean for instance, p0 is the mean value of the density, equal to 1028kglm3 (at atmospheric pressure) and PI Po � w-3. We also introduce reference values for the temperature and the salinity: T o = l0°C and So = 34glkg.
How does the mass conservation equation become '.l = 0? Let us ind in which sense the mass conservation equation can be approximated by '.l = 0.
We introduce the small parameter E such that JI Po= E p(ll I p o , with p(l) I Po = 0(1). We also decompose 1 into J(O) + EJ(l); hence, J(O) is the velocity ield in the limit E-+ 0, that is, in the limit of ininitely small density luctuations. Introducing this decomposition of the ields into the mass conservation equation, one gets (2.4) At order zero in E, the mass conservation reduces to '.J (O). This implies that the mo mentum equations for an incompressible luid are also valid in this limit (that is, with 1 � J( 0 l) . This also implies that the gravity acceleration g should be large enough for the product gpl p0 to stay inite.
How does the energy equation become the usual heat equation? For any constant volume transformation, the internal energy e is equal to CvT where C v is the heat capacity per unit mass and T is the temperature. Equation (2.3) can thus be rewritten as Since the state variable T has been introduced, we need an equation of state. In the limit of very small density deviation about a constant reference density, the luid may be assumed to be incompressible, as just seen. Hence, p only weakly depends on pressure. Density luctuations arise from temperature and salinity changes about the background state so that, at irst approximation, a linear equation of state may be assumed where a is the thermal expansion coeicient and 3 is the haline contraction coeicient.
Assuming that the heat added to the luid volume v results from difusion process, Q can be related to T using Fourier law: pQ = "V 2 T. Replacing Q by this expression in equation (2.5) and using � = -pl.l, one gets (k is the thermal conductivity of the luid) (2.7) In the limit p « p0, that is, 1 � ](O) one eventually recover the usual heat equation The coeicient -0 is the thermal difusity in the luid, usually denoted " · If the density Po v luctuation p mainly arises from temperature luctuations, the above equations become, using (2.6): And what about the momentum equation? We use the assumptions p « p0 and p � Po· In the momentum equations, p comes into play in two terms of quite diferent physical meaning: (i) in the acceleration and Coriolis terms, pDlj Dt and -2p0 \ l, and (ii) in the gravity term gp. While p can be replaced by Po in the former two terms, it certainly cannot be replaced also in the latter term because the luid would be homogeneous otherwise, with constant density p0. Therefore • pDlj Dt -* poDlj Dt ; -2p0 \ u --* -2p00 \ u.
• -lp+ p §-* -lp+ p §, decomposing the dynamical pressure p into Po (hydrostatic pressure) plus p and using the hydrostatic balance of the rest state dPo / dz = -p0g. The equations of motions become The Boussinesq equations are equations (2.9) and (2.10) above along with the incom pressibility condition l.u = 0. We recall that this approximation is strictly valid in the limit u = l( o ) and p «Po·

Non-dimensional parameters
Models of geophysical (or astrophysical) lows are often designed in the laboratory. The idea is to model not the whole system but only a part of it, whereof dynamics are to be studied. For instance, the formation and dynamics of the Great Red Spot of Jupiter has been modelled in a laboratory experiment using two counter-rotating disks, in order to reproduce the zonal shear that prevails on this planet (as in any strongly rotating planet) and controls the dynamics. The scales are very diferent in the natural and laboratory systems but comparison can be performed by considering the non-dimensional equations of motions. These non-dimensional equations involve a few non-dimensional parameters and the application of the laboratory results to the natural system can be done if the non-dimensional parameters are the same. The other interest of working with non-dimensional equations, especially when natural systems are to be studied, is that one deals with values of order 1. Such a non-dimensionalization is actually necessary when the equations of motions are solved on a computer.
We now list the main non-dimensional parameters that govern the dynamics of a Boussinesq low (that is, a low whose dynamics can be described by the Navier-Stokes equations in the Boussinesq approximation) subjected to the Coriolis force.
Let U be the scale of the velocity ield, L be a typical length scale of the luid motions and A few other non-dimensional parameters can be inferred from these four ones. These are for instance the Ekman number, which compares viscous and Coriolis efects: Ek = v 20 £2 = ReRo; and the ratio of rotation to stratiication efects (also referred to as the Prandtl ratio): � = Ro .

Linearization about a basic state
Since we are interested in motions of small amplitude relative to a basic state, it is useful to decompose the ields into a part associated with the basic state and a part associated with the small amplitude motions. Thus we write: u = J + u', p = R + p ' , p = P + p' where J and R may depend upon the spatial coordinates but (this is our assumption) do not depend upon time, R and P are in hydrostatic balance, and the amplitude of u', p' and p' is much smaller than that of U, R and P respectively. Introducing this decomposition into the Boussinesq equations and linearizing about the basic state yield the linearized Boussinesq equations. We shall consider two speciic cases.
• The basic state is a luid at rest in hydrostatic balance. In this case, U = 0 and, therefore, the amplitude of u' is much smaller than 1. The linearized Boussinesq equations become in this case (v is the kinematic viscosity) : (2.11) (2.12) (2.13) • The basic state is a parallel shear low U(z) with density ield R(z) in hydrostatic balance.
The linearized Boussinesq equations of motions become in this case A linear wave in a homogeneous and steady medium 1 is characterized by a frequency n and a wave vector k. For a plane wave in an ininite medium, the motion is of the where A is the wave amplitude. The properties of the wave are completely set by the dispersion relation, namely a relation between n and k: where the function F also involves the luid parameters upon which the wave propagation depends (N, g, etc.) The dispersion relation does not depend upon space and time because of the homogeneity and steadiness assumptions.
Since the dispersion relation characterizes the wave ield, its derivation is an essential step in the wave study. One can distinguish two methods for this purpose. In the classical method, one irst needs to derive the linearized equations of motion. Since the medium is homogeneous, the coeicients of these equations are constant. The general form of the wave solution is then substituted in these equations, from where the dispersion relation is obtained.
A general method that gives access directly to the dispersion relation, without using the equations of motion, can be used instead (Lighthill, 1978). Because it is simpler than the classical method and very powerful, we shall rely on this second method to compute the dispersion relation.

A general method to compute the dispersion relation
The frequency D of an oscillating system, whether consisting of a inite number of degrees of freedom (like the pendulum and or chain of springs) or evolving in a continuous medium, can be computed directly in three diferent ways (Lighthill, 1978): 1. From the restoring force. The square of D is then obtained as fo llows: 12 = restoring force . displacement x mass ' 2. From energetic considerations. 12 is then deined as: where the generalized stifness is the coeicient of the displacement squared in the expression of the potential energy and the generalized inertia is the coeicient of the temporal derivative of the displacement squared in the expression of the kinetic energy; 3. From dimensional analysis. Before applying this general method to internal gravity waves, let us show how it works for a simple discrete system. For this purpose, we consider the harmonic oscillator, whose frequency is D = ;, where K is the stifness of the spring and m the mass of the object attached to the spring.
In case 1., the restoring force is -K.
x where x is the displacement of the mass from equilibrium position. Hence, 12 should be equal to IK.xl divided by lxl times m, which yields the actual expression for D. In case 2., the potential energy is 0.5Kx2 so that the "generalized" stifness is K and the kinetic energy is 0.5mx 2 so that the generalized inertia is m. Finally, in case 3., let us assume that 1 is a function of K and m, of the form K:m3. By writing that the unit of n, m and K are rad.s-1, kg (also Newton.m-1.s2) and Newton.m-1 respectively, one gets a= -3 = 1/2.

Internal gravity waves in an ininite medium
The dispersion relation from general principles. In the following, we assume that the stably-stratiied medium is homogeneous, which implies that N is constant. We shall assume that the medium is steady throughout the paper.
To compute n, we shall use the second expression provided above (case 2.). We irst need to compute the potential energy to infer the generalized stifness. The po tential energy is the work done by the restoring force from the current position to the equilibrium position. We showed in section 1.4 that the restoring force is gp ( so that dz Ep = jo gzzdz = -N 2 p0 jo zdz = 0.5p0N 2 ( 2 . The generalized stifness is therefore p0N2• Let us now compute the kinetic energy to get the generalized inertia. As we show it in the next section, the motion induced by a monochromatic internal gravity wave occurs in a two-dimensional vertical plane, made by k and § (this is the propagation plane of the wave). We assume that this plane coincides with the (x, z) plane. Hence, the kinetic energy is Ek = 0.5p0(u2 + w2). rom continuity equation: kxu + kzw = 0, implying that u = -tan} w, where } is the angle of the wave vector with the horizontal.
2 . 2 , Since, w = (, one gets Ek = 0.5p0(1 +tan })( = 0.5po� } . The generalized inertia is cos Po therefore --2 -. cos } Hence, the frequency (and dispersion relation) of a monochromatic internal gravity waves in a homogeneous medium is (3.5) It follows that internal gravity waves are dispersive waves, because the phase velocity c = 1/lkl depends upon k, with an anisotropic dispersion relation, because n depends upon the angle that k makes with the horizontal.
Main properties of internal gravity waves. The properties of internal gravity waves are described in several textbooks (Lighthill, 1978;Leblond and Mysak, 1978;Cushman-Roisin, 1994) and we refer to these books for a detailed account. We make a brief summary of the important properties which we illustrate with laboratory experi ments.
We cannot avoid discussing the beautiful experiment of Mowbray and Rarity (1967), which was designed to check the dispersion relation (3.5). A horizontal cylinder oscillates vertically at a constant frequency 1 in a constant N luid. The vertical motion of the cylinder perturbs the low and relation (3.5) show that, if .12 :; N2, waves will be produced. Linear waves in a non dissipative medium carry the momentum and the  Mowbray and Rarity (1967). Visual ization is made in a vertical plane and the dark vertical bar is the cylinder support. energy of their source without material transport (f.i. Andrews et al., 1987). The energy is transported at the group velocity 8D/8 k . For internal gravity waves, the group velocity is perpendicular to k . Hence the energy should be carried away by the waves along directions making an angle e with the vertical. This peculiar pattern is displayed in Figure 7 and is referred to as the St-Andrews cross.
Perhaps the most important property of internal gravity waves in a homogeneous stratiied medium is that only a time scale is imposed, namely N-1 . There is no scale selection by the stratiied medium 2 . It follows that in laboratory experiments of internal gravity waves forced by a paddle, the size of the paddle (or, more generally, of the forcing mechanism) selects the scale (f.i. McEwan, 1971). In the experiment by Mowbray and Rarity discussed above, the diameter of the cylinder selects the width of the energetic region. It may thus be expected that in experiments where the forcing mechanism at frequency n has no length scale, all wave vectors which make an angle cos-1 (D/N) with the horizontal will be excited, acording to the dispersion relation (3.5). In the labora tory experiments of McEwan and Robinson (1975) and Benielli and Sommeria (1998), internal gravity waves were forced by parametric instability: a tank stably-stratiied with salted water oscillates vertically. In a frame of reference attached to the tank, the gravity appears to oscillate which triggers a parametric instability which generates in turn in ternal gravity waves. Though the forcing mechanim does not impose any length scale, it eventually appears that the size of the tank provides the scale selection: only the largest wavelength was excited.
When Coriolis efects are important , the waves are named "inertia-gravity waves", with dispersion relation D 2 = N 2 cos 2 e + f 2 sin 2 B where f is twice the projection of the angular vector of the Earth onto the local vertical. In a laboratory experiment, f is Note that ) denotes here the angle that the wave vector makes with the vertical and that the intrinsic wave frequency is denoted as w (from Peacock and Tabaei, 2005).
twice the angular frequency of the rotating tank. It follows that f2 now becomes double bounded: P :; f2 2 :; N 2 . The dispersion relation can also be written as so that, again, one angle should be selected in a rotating stably-stratiied medium by a forcing mechanism with frequency l. The laboratory experiment by Peacock and Tabaei (2005) nicely illustrates this relation ( Figure 8).

The Kelvin-Helmholtz instability in a stably-stratiied luid
Evidence for the Kelvin-Helmholtz instability in natural media. The Kelvin Helmholtz instability is also described by the linearized Boussinesq equations, when the basic low is made of a velocity proile of the hyperbolic (or error function) type in a stably-stratiied medium. The Kelvin-Helmholtz instability develops at the interface between two luids lowing horizontally at diferent velocities. This coniguration is strongly unstable when the two luids have the same density: a perturbation ampliies as soon as the Reynolds number associated with the shear low is of order 1 or larger. As the perturbation ampliies, nonlinear efects come into play which make the instability saturate. The nonlinear development of the instability is manifested as Kelvin-Helmholtz vortices. Examples for this instability are ubiquitous in geophysical lows, whether in the atmosphere (Figure 9) or in the ocean ( Figure 10). In the former case, the instability results from a vertically sheared wind; in latter case, the instability is due to tidal motion over a sill. The Richardson number. A shear low may exist in a constant N stably-stratiied medium but a coniguration most encountered in real lows is the one when the two luids have also diferent densities (the lower luid being of course heavier), leading to a vertical density proile of the hyperbolic tangent type, like the velocity proile. In this case, the density diference between the two luids comes into play and the Froude number is no longer the most appropriate parameter to account for stratiication efects. A bulk (or global or overall) Richardson number is rather introduced, whose deinition comes naturally when the Boussinesq equations are made non-dimensional. The Richardson number is a measure of the stabilizing efect of buoyancy relative to the destabilizing efect of the shear. Let us start from the Boussinesq equations linearized about a basic state consisting of the stratiied shear low (U(z), R(z)) (see the inal set of equations at the end of section 2). We nondimensionalize these equations by choosing the following scales: U, the velocity scale, is half the velocity diference across the shear low; 8, the length scale, is the thickness of the shear low; .p, the density scale, is half the density diference across the shear low. T, the time scale, is equal to 8/ U. It follows that the pressure scale is p0U2.
The non-dimensional equations are (using the same notation for dimensional and dimensionless variables): These equations involve three non-dimensional parameters: the Reynolds number, m �L Re =-, the Prandtl number and the bulk Richardson number J = guz. v

Po
The Richardson number criterion. When the two luids have diferent densities such that the lower luid is heavier, the instability cannot develop whatever the density diference: as luid particles raise because of the instability, kinetic energy is converted into potential energy so that, intuitively, a too strong density diference prevents the instability from developing. Chandrashekar (1961) derived an heuristic criterion for local instability based upon energetic considerations, which was demonstrated rigorously by Howard (1961) (see also Drazin and Reid, 1981, p. 326). The energetic considerations compare kinetic and potential energies, that is nonlinear and stratiication efects, which relates to the Richardson number. Since the bulk Richardson number does not provide any information about the local stability of the shear low, a local Richardson number was introduced for this purpose: The stability criterion is the following: a necessary condition for instability is Ri(z) < 1/4 somewhere in the luid. It has been shown by Hazel (1972) that this condition is necessary and suicient for a velocity and a density proile of the hyperbolic tangent type.

Resonant interaction theory
Resonant interactions in the luid dynamical context were discovered independently by Phillips (1960) and Hasselmann (1962) for surface gravity waves and progressively extended to internal gravity waves by considering, in particular, the interaction of two surface waves with an internal gravity wave (Ball, 1964); the case of an internal wave triad was eventually considered by Phillips (1966) and Hasselmann (1967). Extensive work was next performed in the oceanographic context (see the reviews by Phillips, 1981;Muller et al., 1986;Staquet and Sommeria, 2002). It was shown in particular that energy transfers in the spectrum of large scale, low frequency internal gravity waves in the ocean (the Garrett-Munk spectrum) mostly consist of resonant interactions (McComas and Bretherton, 1977). In the atmosphere, the decrease of density with altitude ampliies the waves, so that resonant interactions are only rarely expected (a signiicant fraction of the wave energy is still backscattered downward, see f.i. Barat and Cot, 1992) . Perhaps one exception is the Antartic atmospheric boundary layer, which is surmonted by a very cold layer of air in which the waves may be trapped and interact resonantly (Rees et al., 2000).
We briely summarize this theory in the following, which we illustrate by results from Koudella and Staquet (2005).
Wave steepness. In the following, we consider a monochromatic internal gravity wave (k, 0). We show below that luid motions may be assumed to occur in the propagation plane (k, §) of the wave if the wave amplitude is ininitely small. Henceforth we assume two-dimensional luid dynamics in the vertical ( x, z ) propagation plane. Let ( 4.1) denote the velocity ield induced by the wave, where 0 = F(k) obeys the dispersion relation (3.5).
The normalized amplitude, or steepness, of the wave is deined as: where ex is the phase velocity in the x-direction. As illustrated in Figure 11, this deinition implies that isopycnals (i.e. constant density surfaces) are nowhere overturned when s < 1 (the wave is also said statically stable) and that the isopycnals are locally overturned when s > 1 (statically unstable wave). The case s = 1 implies that vertical isopycnals have locally formed. Since Umax = Au and Cx = k n = � , the steepness is also deined as The magnitude of the advective term l· 7 [ or of the vortex stretching/tilting term {Vi relative to the magnitude of the baroclinic term 7 x L pez scales as the square of the wave Po steepness. This shows that, when s 2 « 1, the dynamics are not only weakly nonlinear but the luid motions are mostly conined in a vertical plane if the vertical vorticity is initially zero. The linear stability analysis of a monochromatic internal gravity wave to three-dimensional perturbations by Klostermeyer (1991), later discussed in this paper, actually showed that this result holds as long ass S 0.7.
Resonance relations. The resonant interaction theory describes the interaction among a triad of internal gravity waves of ininitely small steepness ( s « 1). Consider a triad of waves (k0, 00), (k1 , 01) and (k2, 02) with Oi = F(ki), 0; i ; 2. Since the waves are nonlinearly interacting, This relation is a spatial resonance relation. If, moreover, (4.5) with Oi of any sign, eicient interactions can occur among the wave triad. This relation is a temporal resonance relation. If the temporal resonance relation is not satisied, the triad interaction produces harmonic perturbation, at frequency Oo + 01 + 02, without cumulative efects on wave amplitudes. These selective interactions are named 'resonant interactions'. They are selective because they imply both a spatial and a temporal resonance relation. It was shown by Hasselmann (1967) that "sum interactions" I01I + I02I = IOol are unstable. In the limit of small-scale secondary waves, that is lk1l, lk2l » lko l , the instability is of the parametric type: 0 1 :' 02 :' -Oo/2. Hasselmann also showed that "diference interactions" 101 1-1021 = I Oo l are stable.

E· �
The energy and the quasi-momentum of the waves (Pi = o: ki for wave i) are conserved within the triad: Eo + E1 + E2 = constant and Po + P1 + P2 = constant. But the Eo Er E2 wave action, deined for each wave by n �-, is not conserved: -+ -+ -= constant Ho no nl n2 (Muller et al., 1986).
Resonant interaction equations. Resonant interactions among a triad occur over a slow time scale, inversely proportional to the wave steepness s (we recall that the theory is valid for s « 1). A slow time scale his therefore introduced, in addition to the (rapid) t time scale t ( := N-1 ) over which the waves evolve: t1 = -. To describe the interaction, s the amplitude and the phase of each wave is assumed to vary over the slow time scale t1. A two-time scale asymptotic expansion of the velocity and density ields is performed and, at irst order in s, the evolution equations of the amplitudes Ai are obtained (f.i. McEwan and Plumb, 1977): where Ti is the viscous damping rate, S is the interaction coeicient, which depends upon the wave vectors of the triad, and the bar denotes the complex conjugate.

Stability of a primary wave of ininitely small steepness
Theoretical results. The amplitude equations ( 4.6)-( 4.8) can be used to investigate the stability of a monochromatic internal gravity wave, termed the primary wave. Let (k0, n0) be the primary wave and assume that A1, A2 « A0. Equations (4.6)-(4.8) may then be linearized about the basic state of the primary wave. An evolution equation for A1 and A2 is obtained (neglecting viscous efects for simplicity): (4.9) It follows immediately that if D1 D2 > 0, exponential growth of the secondary waves is possible, with growth rate equal to S[Ao[�. This condition was shown by Hassel mann (1967) to be equivalent to sum interactions, that is, the triad is unstable. The growth rate of the secondary waves is plotted as a function of k1 in Figure 12 for the run we consider below, which we shall refer to as run R1 for simplicity. For very large wavenumbers, resonance occurs through parametric instability and the growth rate satu rates: there is no scale selection in the inviscid case, arbitrary small scales being excited by the instability. When viscous efects are taken into account, A0 must be larger than a threshold amplitude for ampliication to occur. The smallest scales are damped by viscosity, whatever its (inite) value, so that the viscosity provides the scale selection. In the case displayed in Figure 12, a wavenumber close to 5 is selected. Figure 12. Theoretical predictions of the inviscid and viscous growth rates for an un stable triad as a function of wavenumber I k11, for run R1 (from Koudella and Staquet, 2005) .
Constant contours of the density and vorticity ields of the primary wave are displayed in Figure 13a. The primary wave vector is equal to (1, 1) in a (2r) 2 periodic domain and the wave steepness is equal to 0.36. Since a monochromatic wave is an exact solution of the inviscid Boussinesq equations, a small amplitude noise is superimposed upon the primary wave at initial time to trigger the instability. The Brunt-Vaisala frequency has a constant value equal to 1. This is run Rl.
For simplicity, we assume that frequencies are all positive. Since we shall focus on unstable interactions, the temporal resonance relation should be written as no = nl + 0 2 (sum interactions) and the associated spatial resonance condition is k0 = k1 + k 2 .
The latter values also show that the secondary waves ( k1, Or) and ( k 2 , 0 2 ) are res onantly excited by the primary wave through parametric instability. Indeed, one has 0 1 :: 0 2 :: 00/2. One would have 0 1 much closer to 0 2 if the wavenumbers of the secondary waves were larger but such small scale waves are damped by viscosity. The modulus of the secondary wave vectors in triad T is of order 5 instead, consistently with the scale selection observed in Figure 12.
How do these parametrically excited secondary waves manifest themselves in physical space? Constant contours of the perturbation vorticity are plotted at a given time in Figure 14a. (The dark band on the Figure is the region where the primary wave vorticity is negative and may ignored in this section.) We recall that the perturbation is made at initial time of the ininitely small amplitude noise, from where a few resonant triads emerge as time elapses. The perturbation vorticity appears to be organized as bands of alternate sign, which make an angle 81 such that N cosB1 :: w0j2. Also, the sign of a given band changes sign after one primary wave period and recovers its sign after two periods. These observations are consistent with the perturbation being parametrically excited by the primary wave.
The fact that the bands are associated with parametrically excited secondary waves is conirmed in Figure 14b. This calculation difers from the one displayed in fr ame (a) by the value of the viscosity, which has been decreased by a factor 100. The same bands are visible, with the same inclination (because the excited waves have a frequency close to f10/2), but the wavelength is much smaller because of the smaller scale selection by the lower viscosity.
Energetics of the parametric instability. In this section, we analyse the mechanism through which the perturbation grows. As usual in this type of study, we consider the energetics of the perturbation. The novel point here is that the spatio-temporal dependency of the forcing mechanisms of the perturbation is described (Koudella and Staquet, 2005). We start with a rotation of the reference fr ame (f.i. Mied, 1976) so that Z is aligned with the primary wave vector ko in the new (0, X, Z) frame. We decompose the velocity, density and pressure fields into a part associated with the primary wave and a part associated with the perturbation. The Boussinesq equations linearized about the primary wave state are, for the kinetic and potential energy of the perturbation (Lombard   Koudella and Staquet, 2005) .
and Riley, 1996): (4.10) ( 4.11) where (U', W', P', p') are the perturbation ields in the rotated (0, X, Z) fr ame and BF is the buoyancy lux. In this fr ame, the primary wave velocity is (Uo, 0) ; we denote its density ield by Ro. These equations immediately show that the kinetic energy of the perturbation is forced by the term : -U'W'� , which involves the shear of the primary wave (as for . . 1 dR o any unstable shear low) and the potential energy 1s forced by the term: N 2 p'W' dZ , which involves the vertical gradient of the primary wave density ield.
Let us consider the fo rcing of the kinetic energy. If we assume that the perturbation is made of a superposition of linear intenal gravity waves, the correlation U'W' is one signed. It is not diicult to show here that U'W' :: 0 in the rotated reference frame: the velocity ield of the perturbation is aligned with the vorticity bands displayed in Figure  14 so that U' has the same sign as W'. On the other hand, dU0jdZ oscillates in time and space. At a given position, this implies that the forcing term - the perturbation is forced, when dUo/ dZ is negative; but the fo rcing would be negative when dUo/ dZ becomes positive at that position, implying that the perturbation would give back the acquired energy. How is a non zero net forcing possible? The point is that the perturbation will not give its kinetic energy back if the correlation term U'W' is zero, or very small, during the half primary wave period where dU0jdZ is positive. This will be the case if the perturbation energy is of potential form during this half-period.
It fo llows that optimum energy fo rcing occurs when the perturbation energy is of kinetic fo rm during the half-period when dUo/ dZ < 0 and of potential fo rm during the other half-period when dU0/dZ > 0.
Since dU0/dZ is of the form exp(IP0), the energy of the perturbation should be of this form too. It follows that the velocity and density ields of the perturbation should behave as (U1, W1, '1) 'e xp(I� 0 ). In other words, the instability must be of the parametric type for optimal forcing to occur. What about the potential energy of the perturbation? It can be shown that its fo rcing term �2 p1W1 d d � o is always positive. Since dR0 / dZ oscillates at the primary wave period, this implies that the perturbation potential energy is forced twice per primary wave period: at the locations where dRo/dZ is maximum, where the vertical density gradient experienced by the perturbation is decreased, and at locations where dR0 / dZ is minimum, where the vertical density gradient experienced by the perturbation is increased. The former region is referred to as RSS (reduced static stability regions) and the latter as ISS (increased static stability regions). By contrast, as shown above, the kinetic energy is forced once, in regions where the primary wave shear dUo/ dZ is negative, and this is where vorticity bands are visible in Figure 14. Due to the continuous conversion of kinetic energy into potential energy, it can be shown that potential energy is converted into kinetic energy in a ISS region, the reverse occuring in a RSS region. As a consequence, the perturbation potential energy always increases in RSS regions. A summary of the spatial occurence of these forcing mechanisms is displayed in Figure  13b.
As mentionned above, the perturbation growth is controlled by a few parametrically resonant triads, which may be assumed to be independent as they grow. Hence, the energetics of the perturbation may be reduced to that of a single resonant triad and approximate expressions of the perturbation energy and of energy transfer terms between the primary wave and the perturbation may be computed from this assumption. It can be shown in particular that the phase average over 2r of the kinetic energy transfer term ( 1 1 dUo ) .
( ) -U W dZ to the potentml energy transfer term N 2 p W dZ 1 s equal to cos fh -80 , that is, � 0.91 for run Rl. Hence, potential energy is transfe rred at a higher rate than kinetic energy in an unstable primary wave of small amplitude. This is actually a very general result in stratiied turbulence, where it is found that potential energy transfer toward small scales is higher than kinetic energy transfers. The analogy can be pursued in the present case since parametric instability excites small scale waves. Two arguments can be provided to account for this behaviour : the incompressibility constraint, which applies to the velocity but not to the density ield (Lesieur, 1997) and the subset of two-dimensional interactions which prevent kinetic energy from being transferred toward small scales (Holloway and Ramsden, 1988).
What does happen as the perturbation ampliies? The ampliying perturbation vorticity ield in dUo/ dz < 0 regions may become unstable to a (large scale) Kelvin Helmholtz instability. On the other hand, the perturbation potential energy continuously increases in RSS regions, eventually leading to unstable luid layers; hence, a (small-scale) buoyancy-induced instability may also occur there. However, the unstable luid layers are stabilized by the shear of the primary wave so that, in a two-dimensional computa tion, a Kelvin-Helmhotz instability eventually occurs. In a three-dimensional medium by contrast, the buoyancy-induced instability can develop in the plane perpendicular to the primary wave propagation plane because there is no shear there. Hence, in a three-dimensional fluid, both instabilities compete. Koudella and Staquet (2005) showed that the buoyancy-induced instability grows faster than the Kelvin-Helmholtz instabil ity, unless viscous efects are strong enough, thereby conirming earlier heuristic analysis (Munk, 1981) . It follows that an unstable primary wave of small steepness eventually breaks down through a three-dimensional buoyancy-induced instability (if viscous efects are low enough).

Stability of a primary wave of inite steepness
Theoretical results. When the primary wave steepness is smaller than 1 but not ininitesimal, the use of a linear stability analysis against two-dimensional perturbations allows one to go beyond the resonant interaction theory. The early analyses neglected Coriolis efects. Drazin (1977) and Mied (1976) independently considered the case of a perfect luid while Klostermeyer (1983) extended their results to the case of a real luid. The fundamental conclusion of this two-dimensional analysis is that a primary wave is always unstable, whatever the associated Richardson number (that is whatever the intensity of stratiication N relative to the wave-induced shear) . The instability is again of the parametric type in the case of small scale perturbation. In the limit of a vanishing primary wave steepness, Drazin (1977) showed that the results of the linear stability analysis coincide with those of the resonant interaction theory. By contrast, when the primary wave amplitude is so strong that isopycnals are over turned (s > 1) , the analysis must be conducted in a three-dimensional fr amework. A linear stability analysis of this basic state to two-and three-dimensional perturbations was conducted by Klostermeyer (1991) ; only one orientation of the primary wave vec tor and of the plane of the three-dimensional perturbation was considered (this plane was perpendicular to the primary wave propagation plane (k, g)). Klostermeyer (1991) found that the most ampliied perturbation is three-dimensional and grows at small scale through a parametric instability mechanism. This pioneering work was complemented by Lombard and Riley (1996) and Sonmor and Klaassen (1996) . Largest growth rates were found for perturbations at an angle oblique to the primary wave propagation plane but, as Lombard and Riley (1996) irst realized, the features of the three-dimensional instabil ity (wave number and amplitude) strongly depend upon the primary wave propagation angle (see their Figure 6) .
Numerical results. These results are illustrated by three computations, referred to as runs R2, R3 and R4 (Koudella, 1999) . The irst one (run R2) is analogous to run R1 above, except that the initial primary wave now evolves in three dimensions and is perturbed by a three-dimensional noise of small amplitude (Figure 15a) . Such a small steepness wave is unstable to parametric instability in its propagation plane implying that the low remains two-dimensional during the instability growth (Figure 15b) . The secondary resonant waves amplify most in RSS regions, as we showed it above. This leads to overturned isopycnals (s is locally greater than 1), which are unstable. A buoyancy induced instability develops in the plane perpendicular to the propagation plane of the primary wave and the low becomes three-dimensional (Figure 15c). Like the parametric instability, the buoyancy-induced instability ampliies arbitrary small scale waves so that the viscosity provides the scale selection mechanism. In Figure 15c, a wavelength equal to := 1/7 of the domain size along the y-direction is visible. Because the instability has already reached a nonlinear stage in this frame, pairs of counter-rotating vortices (one pair per wavelength) are actually visible in the Figure. These small energetic scales are eventually dissipated by viscosity and the low relaminarizes, because the primary wave is not forced.
It is remarkable that the validity of the resonant interaction theory extends much beyond its theoretical domain of validity, as we shall show it. We present in Figure 16 a run analogous to run R2 except that the initial primary wave steepness is twice larger (s = 0.724, run R3). An in run R2, the isopycnals locally amplify up to overturning while remaining two-dimensional: the primary wave becomes unstable through parametric instability during this stage (irst two frames of the Figure). The overturned isopycnals become unstable through a buoyancy-induced instability (third frame). Because the viscosity is twice higher in this run than in run R2, a larger wavelength is selected along the y-direction, equal to 1/3 of the domain size. The whole low eventually becomes turbulent and dissipates.
Details of the formation of the buoyancy-induced instability in this run are presented in Figure 17. Snapshots of the density ield are displayed in a plane permendicular to the propagation plane of the primary wave . Three wavelengths are thus visible in the second frame (instead of := 7 for run R2) . The nonlinear development of the instability gives rise to coherent structures with a mushroom fo rm (ifth fr ame), which yield a pair of vortices of alternate sign. A zoom of the formation of the mushroom structures is displayed on the third row.
The formation of mushroom structures as a result of a buoyancy-induced instability has also been observed in laboratory experiments. Figure 18a stems from the experi ment of Benielli and Sonmeria (1998), already referred to in Section 3.3: it displays a visualization in a plane perpendicular to the propagation plane of a large scale internal gravity wave excited through parametric instability; two mushroom structures are visi ble. Figure 18b stems from the laboratory experiment of a stably-stratified shear layer performed by Schowalter et a!. (1994).
The above results show that overturned isopycnals are spontaneously unstable to a three-dimensional buoyancy-induced instability. This is illustrated in Figure 19 through the temporal evolution of three isopycnal surfaces fo r an internal gravity wave of initial steepness larger than 1 (run R4). In this case, the viscosity is 50% larger than in run R3 and two wavelengths only appear only the y-direction. As previously noted, the whole wave ield eventually breaks down through this instability. Such large amplitude internal gravity waves are encountered in the atmosphere, due to the rarefaction of air (lowering of density) as the wave propagates upward. Other mechanisms leading to a local ampliication of the wave are discussed in the next section. secondary instability develops, leading to local overturnings of the isopycnals; in frame (c) , a three-dimensional instability develops upon overturned isopycnals, and in (d), local breaking occurs. BVP is the Brunt-Vaisila period 21rjN (from Koudella, 1999).
Inluence of Coriolis efects on the instability of a primary wave. When Coriolis efects are taken into account, the linear stability analysis of Dunkerton (1997) (see also references therein) and numerical experiments of Lelong and Dunkerton (1998) show that the instability mechanism depends upon the intrinsic frequency of the primary wave. When 1 is close to f, the velocity ield induced by the wave is nearly steady (relative to the perturbation growth rate) and nearly horizontal (with a sinusoidal vertical dependency) so that this velocity ield is close to a parallel steady shear low: it may thus bear a Kelvin-Helmholtz instability. Such an instability is found to occur whatever the primary wave amplitude. By contrast, when n is larger than j, the behaviour depends upon the primary wave amplitude. Thus, resonant interactions occur for a vanishingly small amplitude (Miyazaki and Adachi, 1998). When s has a inite value, even larger than 1, a shear instability prevails as long as 1 is not too large. When s > 1, this instability is more and more modiied by convection as 1 increases and, in the limit of no rotation efect (! » f), it becomes a three-dimensional convective instability .   Figure 16, at successive times. The last two fr ames are a zoom of the buoyancy induced instability. Visualizations are made in the plane perpendicular to the primary wave propagation plane (from Koudella, 1999) .
5 Interaction of internal gravity waves with an ambient shear low

Geophysical context
In geophysical lows, internal gravity waves generally interact with an ambient shear low: for instance, the shear low is the wind in the atmosphere (Bretherton, 1966) or a current in the ocean (Badulin et a!., 1990). The former case is very common : the blowing of a wind over a mountain range generates lee waves which next interact with the wind. Recent research in the Sun has revealed that the diferential rotation that exists in this star creates a very strong shear at the basis of the convective zone, which must interact with the waves that are produced there (Talon et al., 2002). How does this shear low interact with the waves is unknown. An analogous interaction occurs when the waves encounter a temperature (or a density) front. Actually, when Corio lis efects come into play, such a front creates a vertical shear via the thermal wind balance so that the wave interact with both a density fr ont and a shear. And when there is no ambient (large scale) shear or density front , as in the deep ocean, the wave fi eld interacts with the shear low it induces.
Why is the wave-vortex interaction an important geophysical problem? A very simple example is provided by lee waves. Since lee waves are produced by the wind, energy is extracted fr om the wind to generate these waves. As a result, the wind speed decreases over a mountain range because of this generation process, as if the mountain were exerting a drag fo rce on the wind. This drag force has been introduced (as a subgrid-scale model) fo r the alpine range in the ECMWF3 weather fo recast model, thereby improving the weather predictions fo r this region (Lott, 1995). Another example may be provided by the permeability of the polar vortex. The polar vortex forms over the pole in the winter hemisphere and is a quasi-two-dimensional stratospheric vortex. It extends from the bottom of the stratosphere (� 8 kms at the pole) up to an altitude of 70 kms or so, while its horizontal extent is of a few thousand kms. Such a quasi-two-dimensional vortex is known to be a barrier to transport: laboratory experiments by Heijst et al. (1991) fo r instance show that dye injected in an isolated quasi-two-dimensional vortex remains capped inside the vortex and do not mix with the luid outside the vortex. Atmospheric measurements at mid-latitudes still led to the suspiscion that air possibly leaks outside the polar vortex and a mechanism was seeked for to explain this leakage. Mcintyre (1995) suggested that the interaction of inertia-gravity waves (which are internal gravity waves subjected to rotation) with the vortex edge could result in an irreversible transport of mass across the vortex edge. We fi nally report about a subtle mechanism of wave-vortex interactions in the ocean. Klein et al. (2003) showed that meso-scale eddies (which are Figure 19. Breaking of an internal gravity wave with overturned isopycnals at t = 0 (s = 1.086, resolution 256 3 ). Three constant density surfaces are plotted at successive times. First row, from left to right, in unit of the Brunt-Viisili period (BVP): t = 9.7, 10.6, 11.5. Second row, from left to right: t = 11.8, 12.3, 14.1 (from Koudella, 1999).
responsible for the transport of heat toward the poles), through their efects on near inertial (0 � f) oscillations, can induce a spatially heterogeneous difusion that depends on the vorticity sign. In turn, the heterogeneous difusion breaks the cyclone-anticyclone symmetry on short time scales (compared to the difusive time scales) and favors the emergence of cyclonic structures.

Ray theory
In the following, we focus upon the interaction in an ininite medium of internal gravity waves with an ambient large scale horizontal shear low J, and, when rotation is present, with a density front as well.
We consider a monochromatic wave of intrinsic frequency 0 and wave vector k (or a wave packet with main intrinsic frequency and wave vector n and k). The intrinsic frequency is the frequency measured in a frame of reference attached to the ambient low. The simplest, and main, efect of the interaction on the internal wave is the change of the intrinsic wave frequency through the Doppler efect: wo =fl+k.U, (5.1) w0 being the frequency of the source that emits the waves. This source is supposed to be ixed in the frame of reference relative to which 0 is measured. , changes as the wave propagates in the changing velocity ield, and may approach its lower or upper bound. In this case, further propagation is no longer possible and, in the linear limit, the wave may either be trapped or relected. Two academic situations are usually considered. Within an atmospheric context, the mean low is a horizontal wind with a vertical shear U(z)ix; within an oceanic context, the mean low is a horizontal current with a horizontal shear U(y)ix.
The change in . as the wave packet propagates into the current is most easily pre dicted within the WKB approximation (see Olbers, 1981, for a very clear presentation of this approximation). This approximation relies upon the assumption that the properties of the luid medium that afect the wave propagation (U and N, in the present case) vary slowly in time and space relative to the wave intrinsic frequency and wavelength respectively. Hence, the medium may be assumed to be uniform and steady over a length scale of order \k\-1 and over a time of order n-1. Under this assumption, the evolution of the wave vector is known along a ray (deined as dx / dt = c 9 + U) and is driven by the gradients of the ambient velocity and buoyancy ields. In this paper, we shall consider that the luid medium is steady (i.e. oUfot = 0, oN jot = 0). In this case, the absolute frequency is constant along a ray: dw0jdt = 0, where djdt = ojot + (c 9 + U).1 denotes the material derivative following a ray. The equations governing the refraction of the wave vector along a ray, known as the ray equations, are: that is, in the present case, with U = (U(y, z), O,O)  Changes of the wave amplitude are inferred from the conservation of wave action.
For any slowly varying background, the action A = E j., where E is the wave energy, (5.7) Equation (5.7) implies that the action contained in a small volume oV moving with the absolute group velocity is conserved, that is The form of the WKB theory we use is the approximation of geometrical optics but, for simplicity, the terminology WK B approximation will be employed hereafter.

Horizontal mean low with a vertical shear U(z)ix
This situation was irst investigated without rotation efect, when a wave packet propagates upwards in a vertical shear low such that its intrinsic fr equency decreases. If there exists an altitude at which n vanishes, the wave cannot propagate beyond this altitude, known s the critical level. Hence the critical level acts as a wave ilter. This level is also deined as the altitude at which the component of the phase velocity of the wave along the wind direction equals the wind velocity ( Figure 20).
The irst theoretical approach to this problem was performed by Bretherton (1966), who considered the Boussinesq equations in a non-rotating fr ame, linearized about a basic state deined by the wind and a linear stratiication, under the WKB approximation.
Note that this approximation implies that the Richardson number (Ri = N2/(dU/dz)2) is much larger than 1 everywhere (Olbers, 1981). For a wave packet approaching the critical level, Bretherton (1966) showed that the vertical component of the group velocity Cgz decreases as j2, where J is the distance of the wave packet to the critical level, while kz increases as J -1 . From the dispersion relation (3.5), n � 0 (since kx and ky are constant) . Even if these results are correct as the critical level is approched, the WKB approximation diverges in the immediate vicinity of that level. The theory predicts that this level is reached in an ininite time (so that there should be no transmitted component) and that the wave-induced energy increases without bound, as 7-1 .
Relaxing the WKB approximation, but still considering the linearized equations of motions and assuming that the Richardson number is everywhere greater than 0.25, Booker and Bretherton (1967) showed that, for a monochromatic wave, momentum is actually transferred from the wave to the mean low at the critical level, except fo r a weak transmitted component whose energy is the incident wave energy reduced by the factor exp[-2r(Ric -1/4)11 2 ], where Ric is the value of the Richardson number at the critical level. The transfer of momentum to the mean low is manifested as a discontinuity in the vertical lux of wave-induced horizontal momentum pu'w'(z) at the critical level (where the bar denotes an average over a horizontal wavelength) ; this lux is otherwise constant with altitude in the absence of critical level (Eliassen and Palm, 1961 ) . No component is relected in this linear limit.
For a high enough wave amplitude, the incoming wave-induced energy is higher than that absorbed by the shear low and, consequently, the energy density increases in the neighbourhood of the trapping level. Breaking eventually occurs, thereby partly dissi pating the accumulated energy (since momentum is conserved however, momentum will be locally deposited by the breaking waves) . This process is well-known to occur in the atmosphere. However, as noted by Mcintyre (2000), there is a subtle efect which implies that, unless viscous efects are too strong, the linear theory always fails after some time: as the wave approaches the critical level, its intrinsic phase velocity decreases faster than the wave-induced velocity along the wind direction does because of absorption so that their ratio, which is one measure of the wave steepness, increases. Hence, breaking may occur as well.
This situation was studied in detail through several experimental works (f.i. Koop, 1981;Koop and McGee, 1986) and is exempliied in Figure 21a for a weak amplitude wave and in Figure 21b for a large amplitude wave. Koop and McGee (1986) found that when breaking occurs, it does so through a two-dimensional Kelvin-Helmholtz instability when no critical level exists and through a three-dimensional buoyancy-induced instability when there is a critical level. The former behaviour can be explained as fo llows: as the waves propagate in the shear low, the intrinsic frequency decreases so that wave induced motions are close to horizontal, slow and of small steepness; such motions are therefore prone to Kelvin-Helmholtz instability. The latter behaviour is consistent with the general analysis of Lelong and Dunkerton (1998), who show that a gravity wave of low frequency and high steepness breaks through a buoyancy-induced instability. In numerical experiments, the breaking problem was irst addressed in a two-dimensional vertical plane. The idea was that, as kz -+ 0, a strong wave-induced vertical shear develops (which actually dominates the mean low shear) , so that destabilisation should occur via a Kelvin-Helmholtz instability. The two-dimensional simulations realized by Winters and D 'Asaro ( 1989) actually revealed an unexpected behaviour once the strong wave shear has built up, consisting of a long regime of overturned isopycnals without breaking. This led Winters and Riley (1992) to design a model for a gravity wave packet in the neighbourhood of the critical level, as a basic state for a linear stability analysis. They found that the most unstable modes are both of the Kelvin-Helmholtz type ( occuring in the propagation plane of the wave packet) and of the convective type ( occuring in a (a) The wave amplitude is small enough for complete absorption of the wave ield by the shear low to occur; (b) for a larger wave amplitude, the wave breaks below the critical level (from Koop and McGee, 1986).
plane perpendicular to it) . Three-dimensional numerical simulations (Fritts et al., 1994;Winters and D'Asaro, 1994) later conirmed that the dominant instability is a three dimensional convective instability, which the two-dimensional simulations were unable to reproduce. Only Winters and D'Asaro (1994) studied the inluence of the breaking of a wave packet upon the luid medium, for Ri = 25 and for three values of the wave packet initial amplitude. The main result is as follows: the energy carried by the wave packet is irst (in time) transferred to the mean low, the more so when the initial wave amplitude is smaller, in agreement with linear theory. A relected component is next produced, whose energy content grows with the initial wave amplitude (no component seems to be transmitted, whatever the amplitude). Consequently, as noted by Winters and D'Asaro (1994), very little energy is left for luid mixing and for potential vorticity changes.
The work of Booker and Bretherton (1967) was extended to rotating lows by Jones (1967) and by Wurtele et al. (1996) . Using a linearized approach, Jones (1967) showed that, in addition to the classical critical level (characterized by 0 = 0), two additional singular levels exist, corresponding to 0 = ±j. His paper contains another important result: in a rotating luid, the vertical lux of wave-induced horizontal momentum is no longer conserved (away from the critical level) . It should be replaced by the vertical lux of wave-induced angular momentum for this conservation property to hold again (away from the singular levels). Wurtele et al. (1996) further showed that the wave becomes evanescent in the neighbourhood of the critical level (since the wave can only propagate for 02 > j2) so that the only efective singular levels are those where 0 2 = j2 . Wurtele et al. (1996) also investigated numerically the situation of a time developing wave, being emitted from a source at t = 0. In this case, the singularity develops in time as well : at early times, the propagating wave crosses the 0 = f (for instance) level and decays in the evanescent region. As time elapses and the wave reaches a steady state (i. e. its amplitude becomes constant) , the singularity develops but non linear efects develop as well. As a result, the wave breaks in the neighbourhood of the singular level and a relected component is emitted; the wave would be absorbed at the singular level in a linear regime. These results are thus analogous to those found by Booker and Bretherton (1967) at a critical level in a non rotating luid. Wurtele et al. (1996) also showed that the behaviour is notably diferent for a continuous spectrum of fr equencies (such as a lee wave generated by a low over an arbitrary topography) , in which case no singular behaviour is encountered.
In the discussion above, the monochromatic wave propagates in the shear low so that its intrinsic fr equency decreases. When the wave propagates in a shear low so that its intrinsic frequency increases, it relects onto the horizontal plane where 1 = N, in the linear limit. Jones (1968) showed that overrelection occurs when the local Richardson number is smaller than 0.25 somewhere, that is, the shear low is potentially unstable; the wave extracts energy from that low when relecting. When nonlinear efects come into play (due to the high initial amplitude of the wave for instance) and the background shear low is stable, Sutherland (2001) showed that a horizontally periodic wave packet permanently deposits momentum to the mean low at altitudes close to and below the relecting level; when the wave packet is horizontally compact, a substantial part of the inite amplitude wave packet energy may be transmitted.

Horizontal mean low with a horizontal shear U(y)ix
Theoretical results. The behaviour of an internal gravity wave packet in a hori zontally sheared current U (y )ix was studied thirty years ago within an oceanographic context, using the linear theory in the WKB approximation (Ivanov and Morozov, 1974;Olbers, 1981;Basovich and Tsimring, 1984;Badulin et al., 1985;Badulin and Shrira, 1993) . Coriolis efects are ignored in all studies except in Olbers' work. When the wave packet enters into the current and propagates against it so that its intrinsic frequency increases (Eq. (5.6) with x; = y), it cannot propagate beyond the position Yt at which 1 = N. The mean low being barotropic, this position actually is a vertical plane if N is constant; this plane is hereafter referred to as the trapping plane ( Figure 22). Since the properties of the medium in which the wave propagates vary only with y, kx and kz remain unchanged (as well as w0 as already noted, because the medium is steady) . From the linear dispersion relation (3.5), one easily infers that ky goes to ininity as 1 -+ N. More precisely, it can be shown from WKB theory that ky ' y-112 as t -+ 0, where t refers again to the distance of the wave packet to the trapping plane; moreover, Cgy ' y312 and c9x , Cgz -+ 0 as well (Staquet and Huerre, 2002). The latter property implies that the wave packet slows down in the neighbourhood of the trapping plane so that its energy density locally increases (in other words, the wave-induced energy accumulates in the neighbourhood of this plane) . However, the WKB theory also predicts that this energy tends to ininity (E ' y-312) and that the trapping plane is reached in an ininite time. As in the critical level situation, the two latter results are unphysical and stem from the failure of the WKB theory in the immediate vicinity of the trapping plane. Olbers (1981) actually noted that, if E' y-1', the asymptotic behaviour of the wave is Trapping plane: Figure 22. Sketch of a wave packet approaching a trapping plane, within the WKB ap proximation, and propagating against the current (k.U = kx.U < 0). Two rays are shown in the horizontal (x, y) plane, with Cg h + U(y)ix being the component of the absolute group velocity in that plane. The intrinsic frequency of the wave packet, w0 -U(yt)kx, from the Doppler relation, is equal to the Brunt-Viisili frequency at the trapping plane (adapted from Staquet and Huerre, 2002).
regular or singular depending upon whether L < 1 or L : 1. In the former case, the wave reaches the plane in a inite time, and this plane is for instance a relexion plane. In the latter case, the wave reaches the plane in an ininite time; very strong gradients of the wave-induced Reynolds stress form in the close neighbourhood of the plane, which yield momentum exchange with the mean low. Information upon the actual wave behaviour when L : 1 can be obtained by solving the linearized equations of motions. Ivanov and Morozov (1974) thus found that the total wave-induced energy may indeed increase, as opposed to the critical level situation studied by Booker and Bretherton (1967). In this situation, momentum is transfe rred fro m the shear low to the wave, so that the potential for wave breaking exists. The three-dimensional numerical study of Staquet and Huerre (2002), whereof a summary is provided below, shows that an inertia-gravity wave packet may indeed break in the neighbourhood of a trapping plane. When the wave propagates along the current such that its intrinsic frequency de creases, ky decreases as well. If ky decays down to zero, n reaches a minimum value Omin, obtained by setting ky to 0 in the dispersion relation. The behaviour of the wave as 0 --+ Om in may be guessed by using WKB theory. Note that the theory becomes less and less valid as the plane gets closer since the wavelength along the y-direction increases. Assuming that the theory remains valid, ky v y112' Cgy v y112 and E v y -1/2 as y--+ 0. Here, according to Olbers' (1981) argument , L = 1/2 so that the 0 = Omin plane is a relexion plane for the wave . Oilers et al. (2003) solved numerically the equation for the amplitude of a hydrostatic internal gravity wave emitted away from the shear low and propagating toward it. The wave behaviour (transmitted, relected, over-relected) depends upon the stability of the shear low. When the latter low is inertially stable, the wave is always relected, with a possible transmitted component. When the shear low is inertially unstable by contrast, over-relection is possible, with a relection coeicient up to 3.25 (for the case considered in the paper) . One may wonder whether a high enough relexion coeicient would not lead to the instability of the relected wave. Note the analogy between this behaviour and that found by Jones (1968) for an unstable vertical shear low.
Numerical results. The behaviour of a wave packet in a rotating, constant-N medium propagating into a barotropic shear low U(y)ix has been investigated numerically by Staquet and Huerre (2002). The shear low consists of a horizontal shear layer (with a tanh proile) while the wave packet is a plane wave whose amplitude is modulated by a gaussian function along the y-direction. The parameters of the wave (its wave vector) , of the shear low (its maximum amplitude) and of the medium (the Brunt-ViisaJi and Coriolis frequencies) are chosen so that (i) the wave intrinsic frequency increases as the wave propagates in the shear low and (ii) a trapping plane exists.
The wave behaviour is illustrated in Figure 23 for a cyclonic shear low (that is, the vorticity of the shear low is of the same sign as the Corio lis frequency) . Constant contours of the density ield are displayed at successive times in a vertical (y, z) plane. At t = 0 (Figure 23a) , the wave packet is hardly visible due to its small steepness (s = 0.26). Because the horizontal shear low does not displace the isopycnals, it is not visible either. Since N is constant , the trapping plane is a vertical plane and its intersection with the (y, z ) plane is marked with a vertical line in the Figure. The wave packet exhibits two major changes as it propagates toward the shear low: the isopycnals steepen and the wave amplitude increases. The fo rmer efect is accounted fo r by noting that, as ky increases, the incompressibility condition reduces to kyuy + kzUz � 0. Hence, the vector (ky , kz) is perpendicular to the phase lines in the (y, z ) plane. Since ky _, +o while kz remains constant, the phase lines steepen. The second efect results from the trapping of the wave . The local increase of the wave amplitude makes the wave packet break ( Figure  23e) and small-scale motions are produced. The latter motions are quickly dissipated however because the primary wave packet is not fo rced. As well, the shear low is hardly modiied by the momentum deposit that occurs during wave breaking.
The stage of the low that fo llows breaking dramatically changes when the shear low is anticyclonic ( Figure 24) . Indeed, in this situation, the shear low is subjected to an inertial instability, through which small-scale motions are most ampliied. The point is that the small-scale motions resulting from wave breaking act as a perturbation to the shear low, which triggers the instability. The medium is therefore considerably modiied by the breaking of the wave in this situation, because it initiates the inertial instability of the (very energetic) shear low. The latter instability results in momentum and mass transport: the shear of the background low is weakened and a passive scalar is transported across the trapping plane, namely across the shear low.

Horizontal mean low with both a horizontal and a vertical shear
What does happen if a wave packet in a stably-stratiied rotating medium interacts with a (thermal wind) balanced shear low U(y, z) involving both a horizontal and a vertical shear? As discussed above, the situation is not simple. For instance, when the intrinsic frequency 2 increases and approaches N, the wave packet should be relected by the vertical shear aU I az but trapped by the horizontal shear aU I ay. Also, since the shear low satisies the thermal wind balance, a buoyancy ield B(y, z) exists as well such au aB · that f-= --. az ay We have investigated this situation in a simple context, when a wave packet conined both in the y and z directions interacts with a horizontal shear layer with a sinusoidal vertical dependency (Edwards and Staquet, 2005). We chose the parameters such that l increases because of the horizontal shear of the background low. We irst explored the parameter range by solving the ray equations and performed three-dimensional di rect numerical simulations (DNS) to investigate the inluence of nonlinear efects on the behaviour of the wave packet .
A general behaviour is observed in the WKB theory, which is illustrated in Figure 25. The ray trajectories are displayed in frame (a) and the intrinsic fr equency l is plotted versus time in frame (b). In fr ame (a), the wave packet at initial time is represented by a set of twenty rays aligned along two perpendicular segments. The l = N surface, plotted with a thick dashed line, displays two important locations when the interaction with the wave packet is considered: (i) where the local radius of curvature is minimum, corresponding to a maximum value of aU I ay (location 1) and where the local radius of curvature is maximum, corresponding to a minimum value of aUiay (location 2). In the former case, the l = N surface is nearly vertical and is a trapping surface; in the latter case, the surface is nearly horizontal and is relecting. All rays propagate toward the l = N surface and reach it, either in the neighbourhood of location 1 or of location 2. In the former case (location 1), the rays are trapped at the surface and propagate along it downward, with a nearly vertical group velocity, toward location 2. Note that the group velocity has strongly decreased when location 1 is reached. All rays sooner or later reach location 2 and relect there. The rays then propagate in the interior of the shear low within a wave guide made by the l = N surface. rame (a) also displays grey points, at which the steepness of the wave packet exceeds 1. This suggests that breaking may occur there, resulting in irreversible mass and momentum transport. Note however that WKB theory is no longer valid at the trapping plane (and at a relecting surface) and that, most importantly as we shall see, molecular efects have been ignored in this analysis.
DNS results are displayed in Figure 26 through constant contours of the luctuating buoyancy ield b'. It should be stressed that the steepness of the wave packet is twice smaller than in the baroclinic cases and that the varying horizontal shear (aU I ay) is not larger. Hence, the wave-shear interaction (which scales like s2 aU I ay) is smaller by a factor 4 at least. The interaction between the shear low and the wave packet is therefore weak and one may consider that the buoyancy luctuations displayed in Figure 26 solely belong to the wave ield. The DNS behaviour is close to the WKB prediction up to the time the wave packet reaches the trapping surface at location 2 (at t :' 176 :' 28 BVP). This is attested in fr ame b) where the intrinsic fr equency predicted by WKB theory and the numerical simulation are compared. Molecular efects deeply change the subsequent wave packet behaviour in the DNS, for two reasons. First, ky increases as the packet approches the trapping surface, implying that small scales along the y-direction are produced. These small scales are of course very sensitive to molecular efects. The Results are compared over the duration of the DNS (t=240, that is 38.2 BVP) (from Edwards and Staquet, 2005).
second reason, connected to the irst, is that the packet slows down as it approaches the trapping plane, which makes it also prone to molecular efects. As a consequence, the nearly steady, small-scale packet is dissipated locally and does not penetrate into the wave guide. The general conclusion that can be drawn from this study is that a single wave packet is unlikely to modify its environment. In geophysical lows however, waves are most often generated by a permanent source (like the interaction of the tide with the topography in the ocean) or, at least, are emitted during a long time with respect to their intrinsic period. Since, in our study, the wave packet slows down as it approches the trapping surface, a continuous emission of such packets would result in their superposition in the neighbourhood of the relecting surface, and possibly in breaking. Higher resolution simulations should also be conducted to reduce the inluence of molecular efects on the wave packet behaviour. turbulence. Coriolis efects are ignored in this paragraph. When an internal gravity wave relects onto a sloping boundary, the frequency is con served upon relection, as fo r any relecting wave in a steady medium. The dispersion relation (3.5) thus implies a peculiar property (Phillips, 1966): the angle that the wave vector makes with the horizontal is preserved as well upon relexion4 . This geometrical property implies that a ray tube (made of curves everywhere tangent to the group ve locity) relecting on a sloping boundary may be fo cused (see Figure 27). Since the lux of energy across the tube section is preserved during relexion, the energy lux increases after relexion. A particularly interesting situation occurs when the angle that the inci- Figure 27. Relexion of an internal wave on a sloping wall near the critical angle. The angle ) of propagation with respect to the vertical is conserved, so the relexion angle Er is much smaller than the incidence angle Ei , resulting in wave fo cusing (from Staquet and Sommeria, 2002). dent ray makes with the horizontal is equal to the angle of the slope ("critical angle" ). In this case, the linear theory predicts a relected wave of ininite amplitude, ininitesimal wavelength and zero group velocity. The Reynolds number is actually conserved after re lexion (the wavelength decreasing in inverse proportion to the wave-induced velocity) so that, if the incident Reynolds number is large, turbulent efects are expected to occur in the relected ray tube. This situation has been investigated theoretically using a weakly nonlinear approach (Dauxois and Yo ung, 1999;, experimentally (Ivey et al., 2000;McPhee-Shaw and Kunze, 2002;Peacock and Weidman, 2005) and with nu merical simulations (Slinn and Riley, 1998;Zikanov and Slinn, 2001). Dauxois and Young (1999) predicted the generation of a harmonic component, which was recently conirmed by the laboratory experiment of Peacock and We idman (2005), and the formation of an array of counter-rotating vortices along the wall. Slinn and Riley (1998) showed that a strongly turbulent boundary layer can form along the wall with eicient mixing property and, according to the laboratory experiment of McPhee-Shaw and Kunze (2002), this results in irreversible mass transport along the isopycnals and therefore away fr om the boundary, in so-called nepheloid layers ( Figure 28) (see Staquet and Sommeria, 2002, for further discussion).
In a closed domain in which one wall (at least) is inclined to the vertical, the solu tions of the linear equation are not longer regular solutions (seiches) but have a fr actal structure. In this situation, the incident ray tubes focus at each relexion on the slop ing wall. It has been shown by Maas and Lam (1995) that the focused tubes actually converge toward an "attractor" , namely toward a localized region in space (Figure 29a). A laboratory experiment has been performed to check the validity of these theoretical ideas (Maas et al., 1997): the occurence of the attractor is clearly visible in Figure 29b. The occurence of turbulence and mixing along the attractor is still an open question.

The internal tide
The internal tide is the ield of internal gravity waves generated in the ocean by the interaction of the barotropic tide (the common tide) with the bathymetry. Munk and  Laboratory experiment of the attractor in a tank stably-stratiied with salted water. Horizontal lines of dye have been introduced at initial time, whose subsequent motions attest of the internal wave formation. To make the attractor more visible, the image of the low at initial time has been substracted to that at a time where the attractor has formed. The localized presence of dye lines attest of localized luid motions (from Maas et al., 1997). Wunsch (1998) suggested that the breaking and mixing of the internal tide could fun damentally contribute to the global thermal equilibrium of the ocean: the cold bottom water masses in the thermohaline circulation indeed raise toward the surface through mix ing; and mixing in the abyssal ocean, away from the boundaries, mostly occurs through internal wave breaking. Several questions therefore stem from this argument : what is the rate of conversion from the barotropic tide into the baroclinic tide? How and where does the internal tide lose energy? In particular, where does breaking occur? Since the conversion of energy is largest at continental slopes, we considered the generation of the internal tide at a continental slope.
The numerical modelling of this problem is not simple because the internal tide dy namics are nonlinear, non-hydrostatic, three-dimensional and a bathymetry is involved. Moreover, a continental slope typically extends from a depth of about 200m down to the abyssal plain, at 4000m. We are currently studying this problem with the code developed by John Marshall at MIT (Marshall et al., 1997) with the investigation of breaking and mixing as the main objective. Figure 30a displays the internal tide velocity ield for a linear hydrostatic computation performed by Gerkema (2002). This simpler context allows one to focus on the generation and propagation of the internal tide. The Brunt-VaisaUi frequency is constant in this situation and equal to 2 w-3 rad/s, the Coriolis frequency is that of mid-latitudes (lo-4 rad/s) and the topography has a simple shape. Figure 30a has two remarkable features: the internal tide is generated at the shelf break and manifests itself as a single ray relecting at the bottom and at the surface. Simple qualitative arguments can be provided to account for these features. Because the tide must jump over the slope, its velocity ield acquires a vertical velocity over the slope which is largest on top of the slope. Let us consider a frame of reference attached to the barotropic tide. In this fr ame the topography has an oscillatory motion at the tidal frequency. This is analogous to the motion of an oscillating body in a stably-stratiied medium already discussed. An internal gravity wave ield is thus generated in the region of largest vertical velocity (if the oscillating frequency is comprised between the local values of N and f), whose energy propagates along a cross with angle set by the dispersion relation. In the present case, only parts of the cross can form because the body (the topography) is not spatially bounded, which propagate toward the abyssal plain and toward the shelf.
The inclusion of nonlinear non-hydrostatic efects strongly modify this regular pattern ( Figure 30b ). As discussed by Gerkema et al. (2005), an internal wave beam in an ininite medium is a solution of the fully nonlinear Boussinesq equations (see also Tabaei and Akylas, 2003). Hence, nonlinear efects are expected to occur where the beam meets a boundary. Figure 30b shows indeed that small-scale structures are observed in the generation region and where the beam meets the abyssal plain. One can explain the former observation by noticing that the beam is tangent to the slope in the generation region and is therefore aligned with the barotropic low, so that forcing by the latter low is optimum. In the relexion region, higher harmonics are generated as a result of the superposition of the incident and relected beams. As discussed in Section 6.1, the generation of a second harmonic, predicted theoretically by Dauxois and Young (1999), has been observed in a laboratory experiment of Peacock and Weidman (2005)  are observed however in the experiment, as opposed to the present simulation, because the internal wave beam is not fo rced.
It follows that higher frequencies are generated once the beam relects, and this is attested by fr equency spectra of the kinetic and potential energy ( Figure 31). The slope spectrum is close to -2 and is made of quasi-linear waves, suggesting that the spectra may be of the Garrett-Munk type (see next Section).
A practical application of this study is that two well-distinct regions appear to be turbulent, as further discussed in Gerkema et al. (2005), so that any parameterization of the internal tide may be seeked fo r independently in each region. Of course, this point should be complemented by simulations in a more realistic context. 7 Statistical properties of "wave turbulence" As we recall it in Section 3.3, internal gravity waves transport the momentum and energy of their source, without any transport of matter as long as their dynamics are linear and non-dissipative (see, f.i., Andrews et al., 1987). The degree of nonlinearity depends upon the length scales involved and is quantiied by the wave steepness (deined in Section 4.1).
In the ocean, s « 1 for vertical wavelengths comprised between a few tens and a few hundreds of meters and the wave dynamics are weakly nonlinear5 (Garrett and Munk, 1979). Ve locity and temperature spectra measured in the ocean at these scales remarkably display the same dependency on the vertical wavenumber and frequency Kolmogorof type (v k-513) fo r a vertical wavelength smaller than :o 10m. In fr ame (c) , the spectrum stems from a laboratory experiment of mixing of a density interface by an oscillating grid (Hannoun and List, 1988) . For a certain parameter regime, internal gravity waves develop at the interface and break. In frame (d) , the result fr om a two dimensional (in a vertical plane) direct numerical simulation of a large scale breaking internal gravity wave in a constant N stratiied luid at rest is presented (Bouruet Aubertot et al., 1996) . The spectra are averaged over 14 Brunt-VaisaHi periods (BVP) around the breaking time (t :o 35 BVP) . An identical k;3 spectrum is obtained when the computation is performed in three dimensions (Carnevale et al., 2001) .

Mixing by breaking internal gravity waves
The importance of mixing stems fr om the fact that small-scale processes, such as internal gravity waves, have an impact upon their environment only through irreversible processes, such as mixing. Hence, in large scale circulation models, these small-scale motions are usually represented through their mixing properties, via a turbulent difusivity. Internal gravity waves deposit momentum when they break so that this efect should also be taken into account through a (isopycnal) transport coeicient.
It is not possible to ignore mixing from these small scale motions. As we already stressed it, in the ocean fo r instance, the cold bottom water masses raise through mixing and it is very likely that the thermohaline circulation, with typical time scale of one thousand years, would stop without mixing, with typical time scale of one minute. As well, in the atmosphere, the climate dynamics strongly depend upon mixing processes, which arise in convective clouds (Bony et al., 1995) .
In this section, we only address basic concepts of mixing. A more extensive presenta tion along with measurements of mixing in a few laboratory and numerical experiments and in geophysical luids can be fo und in Staquet (2004) .

Basic concepts
Motions in a stably-stratiied luid are associated with kinetic and available potential energy (APE) . The APE is partly reversibly converted into kinetic energy through the buoyancy lux, and partly irreversibly converted into background potential energy (BPE) through mixing (Lorenz, 1955;Thorpe, 1977;Winters et al., 1995) . The BPE is the potential energy of the background density ield Pb(z), namely the density ield that would exist in the luid if it were brought to rest. The BPE is therefore the minimum potential energy of the luid. The physical mechanism to account for the irreversible conversion rom APE into BPE is simple: when luid motions occur, luid particles of diferent density are brought into contact and mix through molecular efects. This point can be discussed further. Mixing between two particles depends upon the time the particles are in contact. In other words, the Reynolds number of the luid motion that advects the particles comes into play: if luid particles are transported by large scale motions, with large Reynolds number, the difusive time will be large relative to the advective time and mixing hardly occurs (but stirring does) . By contrast, small scale motions will lead to more mixing. These ideas are illustrated in Staquet (2000) for the in the oceanic thermocline (from Gregg, 1987). b) Stratospheric temperature spectra measured from instrumented balloons, showing a k;3 range with a Kolmogorof k;5 /3 tail at high kz (provided by F. Dalaudier). c) Normalized internal wave spectra, computed from luoresceine concentration proiles in a laboratory experiment of mixing of a density interface by an oscillating grid. 6b is the buoyancy jump across the density interface, lo is an integral scale of turbulence and R j is a Richardson number deined by ,b l0ju §, with u0 being the rms of the horizontal velocity measured in a homogeneous luid at the same distance from the oscillating grid as the density interface (from Hannoun and List, 1988). d) Two-dimensional direct numerical simulations of a large scale breaking internal gravity wave (ky is the vertical wavenumber) (from Bouruet-Aubertot et al., 1996).
At the scale of the luid system, mixing is quantiied by the rate of change of the BPE. This rate of change is always positive, a consequence of the second law of thermodynam ics. This rate is also the dissipation rate of the APE and will be referred to as the mixing rate EAPE · The eiciency of mixing is characterized by the ratio of this rate of change to the amount of energy brought to the system per time unit. If a statistically steady state is reached, the amount of energy brought to the system per time unit is equal to the rate of energy dissipated per time unit in the system, namely EAPE +E KE· The eiciency of mixing is therefore EAPE and is referred to as the lux Richardson EAPE +EKE number in the oceanographic community (see Toole, 1998, for a discussion of mixing in the oceanographic context) . The crucial point is therefore to compute EAPE· 8.2 How to compute the mixing rate EAPE?
An exact expression for the APE is presented in Holliday and Mcintyre (1981), in terms of an expansion that holds for general N(z). For constant N only the irst term of the expansion remains, giving APE=-! / < p 1 2 > , where < > denotes a volume 2up dz average and p(z) is the horizontally averaged density proile. Note that for a constant N, the scale of the vertical gradient of the density proile is theoretically ininite so that this case corresponds to the linear limit. From a practical point of view, this expression of the APE is also used when the displacement of the luid particles is small with respect to the scale of the vertical density gradient p / dz. In this case, fr om the linearized Boussinesq equations, the dissipation rate of APE is simply -: p� dz < l 7 p' l 2 > . This method cannot be used when the particle displacements are not small with respect to the scale of the density gradient (but the Boussinesq approximation still holds) .
The stably-stratiied shear layer is one example. In this case, Winters et al. (1995) proposed an exact analytical expression for EAPE, based upon the background density proile. Hence, one has to compute this proile and two methods have been proposed for this purpose. The irst method was proposed by Thorpe (1977) and extended to two and three dimensions by Winters et al. (1995). The background density proile is obtained by sorting adiabatically the density ield so as to get a stable proile (note that the energy released by the sorting process is the APE) . rom a practical point of view, one sorts the values of the density ield at each grid point so that a stable proile is reached; each plane of the numerical domain (supposed to be parallelepipedic) , starting from the lowest altitude, is then illed with the heaviest particles and the density is averaged over each plane to get a monotonic density proile. Tseng and Ferziger (2001) proposed another method to access the background density proile, which is more eicient numerically than the sorting method and relies upon the probability density function of the density ield.
As said above, these theoretical approaches have been developed to compute the mixing eiciency. Mixing is also quantiied by a (diapycnal) turbulent difusivity Kp , dp b deined as EAPE/9 = -Kp z · In the linear limit, p :: : Pb for all times and one recovers < IVP'I 2 > the expression widely used in oceanography: Kp = -, ( 5 /dz)2 .

Conclusion
We tried in these Notes to broadly describe the very rich dynamics of internal gravity waves. The very study of those waves started in the mid-sixties, when Ball (1964) realized that two surface waves could resonantly interact with an internal gravity wave . Much has been learned about these waves, especially from a theoretical point of view, although the dispersion relation of inertia-gravity waves has just been veriied experimentally (Peacock and Tabaei, 2005). And much remains to be done when natural media are considered. The sources of internal gravity waves in the atmosphere, other than the wind blowing over orography, are still the subject of active research. For instance, convective clouds at low latitudes generate such waves, which may interact in turn with the clouds. As well, the analysis of atmospheric data (FASTEX campaign) revealed that energetic large-scale low-frequency waves, with velocity comprised between 8 and 10 m/s, are emitted by the tropospheric jet (Plougonven et al., 2003). The important role of internal gravity waves in mixing the abyssal water masses has been conjectured but not proved. And the role of these waves in the dynamics of the radiative zone of the Sun is totally unknown. From a modelling point of view, weakly nonlinear theories have been developed to account fo r the statistical properties of those waves (see Staquet and Sommeria, 2002). But, from a practical point of view, no reliable parameterization of the waves exists, which would relate the transport coeicients they induce to the dynamical parameters of the large scale motions which create them. Such parameterizations are crucially needed in general circulation models of the atmosphere and the ocean.
10 Acknowledments The author thanks V. Armenio and S. Sarkar for their invitation to the CISM Summer School on "Environmental stratiied lows" in Udine, for very fruitful discussions and for providing the opportunity to write these lecture notes.