An Efficient Solution Method for the Buoyancy Wave Equation at Variable Wind and Temperature

By Room, Rein Zirk, Marko

ABSTRACT To solve a horizontally spectral, vertically discrete buoyancy wave equation in conditions of arbitrary wind and temperature distribution with height, a novel method is applied, which consists of a presentation of the solution in the form of a cumulative product of complex decrease factors. For decrease factors, a nonlinear, inhomogeneous, two-member recurrence formula follows that is initiated, assuming the radiative condition at the top. Singularities of the wave equation, corresponding to a critical layer in the vicinity of evanescent wind, are eliminated by turbulent friction. The estimation of minimal vertical resolution is derived, enabling solution stability and accuracy. The areas of application of the developed numerical scheme are the high- precision modeling of orographic waves for arbitrary orography in general atmospheric stratification conditions and testing of adiabatic kernels of numerical weather prediction models.

(ProQuest: … denotes formulae omitted.)

1. Introduction

This paper describes a numerical method for solving the linear mountain wave equation. The approach takes its beginning in early investigations by Queney (1948) and Scorer (1949) and is similar to the solution method by Smith (1980), Shutts (1998), Room and Mannik (1999), and Broutman et al. (2003). In this approach, the mountain wave equation is first reduced to an ordinary differential equation (the so-called wave amplitude equation) in a vertical coordinate by applying Fourier transform in the horizontal (and in time, if a nonstationary problem is studied), and then the amplitude equation is solved either analytically (ibid) or numerically, like in the layered approach by Smith et al. (2002). In this paper, a new finite- difference approach is introduced to solve the buoyancy wave amplitude equation in vertical, which allows for arbitrary vertical variation of stratification and directional shear in the ambient flow, and thus expands the area of practical applications of the linear mountain wave theory a while.

The proposed solution method is based on the presentation of the wave amplitude in the form of a cumulative product of the decrease factors. Those complex factors have clear physical content: the modulus presents the decrease of the wave amplitude per single layer of the discrete model, whereas its argument is the phase angle increment per layer. The obtained cumulative solution is then rescaled to satisfy the lower boundary condition. The decrease factors are computed via a straightforward recurrence formula, whose initial value is specified from the radiative condition at the top. The numerical algorithm proves economical and fast, as the computation of the recurrence will take little effort and the major computational time goes to the preparation of the equation coefficients and summation of the obtained orthogonal modes over the wavenumbers to get the solution in the ordinary physical coordinates. Thus, extremely high vertical resolutions (up to 1 m) and small horizontal grid lengths (10100 m) are achievable, if such extreme resolutions are required by stratification and orographie conditions.

As the model deals with variable winds, problems arise inevitably with critical wave vectors and critical levels, corresponding to singularities of the wave equation. This problem is solved with the inclusion of turbulent friction in the forcing terms, yielding singularity removal. The use of turbulent viscosity for the wave equation regularization purpose was proposed already by Lin (1955), Jones (1967), and Hazel (1967). However, the theoretical estimates of maximum stable vertical grid step will show that the requirements to high vertical resolution remain in the vicinity of critical levels. This is the point where the numerical efficiency of the method becomes crucial, enabling the application of sufficiently high spatial resolution where appropriate.

While there exist various “ready” wave equations, they do differ rather substantially in appearance, depending on the details of the dynamical model, like the used coordinate system, etc. To avoid potential ambiguities in initial definitions, we start with a short introduction of the wave equation (which is used in this investigation) from the Miller-Pearce-White (MPW) model [a certain nonhydrostatic, semielastic pressurecoordinate model developed by Miller (1974), Miller and Pearce (1974), Miller and White (1984), and White (1989)].

2. Continuous spectral wave equation

When using the nondimensional log pressure coordinate zeta = ln(p^sub 0^/p), [implies] p = p^sub 0^e^sup -zeta^, (p^sub 0^ = 1000 hPa is the mean sea level pressure) instead of the common pressure p, the linearized version (Room 1998) of the MPW model reads

… (1a)

… (1b)

… (1c)

and

… (1d)

Dynamic fields are the omega velocity omega = Dp/Dt (with D/Dt as the material derivative), the temperature fluctuation T, the nonhydrostatic geopotential fluctuation phi, and the horizontal velocity fluctuation v = {u, v}. The gradient operator [backward difference] and Laplacian [backward difference]^sup 2^ are strictly horizontal. Reference (background) state of the atmosphere is given by a stationary, horizontally homogeneous wind vector U(zeta) = {U(zeta), V(zeta)} with U^sub zeta^ = [partial differential]U/ [partial differential]zeta, background temperature T^sup 0^(zeta), a gas constant for dry air R, and constant Coriolis parameter f = |f|, where f is a vertical vector. The scale height H and the stability parameter theta are

…

where gamma = gamma(zeta) is the height-dependent kinematical turbulent viscosity coefficient. We apply the simplified turbulent viscosity model with the viscous terms acting only in the horizontal. This is because the main purpose of viscosity here is to regularize singularities of the wave equation on critical levels (see farther on).

Using Fourier presentations

…

where i = [the square root of]-1 is the imaginary unit, the amplitudes … are functions of zeta, k = {k, l} is the wave vector, and v^sup 0^^sub k^ is the corresponding eigenfrequency, system (1) transforms to the spectral normalmode equations

… (2a)

… (2b)

… (2c)

and

… (2d)

where

… (2e)

is the intrinsic frequency or eigenfrequency. The eigenfrequency is complex in the presence of turbulent viscosity, which will produce a weakening of the free orthogonal modes in time, or, in the case of a stationary solution, a downstream and upward weakening of the wave amplitude in comparison with the inviscid case.

It is straightforward to derive from (2) a single scalar equation for the spectral amplitude of omega velocity (short notation … is used in the following):

… (3)

with zeta-dependent coefficients

…

where N = [the square root of]Rtheta/H is the Brunt-Vaisala frequency. The spectral wave equation [(3)] was derived for the general nonstationary case with nonzero v^sup 0^^sub k^ [with the aim of solving the nonstationary problem with a future application of (3)]. From this point on, the treatment will be confined strictly to the stationary solution with v^sup 0^^sub k^ = 0. In the stationary case, the bottom boundary condition for (3) is the no- penetration condition in the spectral presentation

omega(0) = ik . U(0)p^sub k^(0), (4)

where p^sub k^(0) is the Fourier amplitude of the mean surface pressure p^sub s^(x, y) = p^sub 0^ exp[-h(x, y)/H(0)], corresponding to the orographic height distribution h(x, y). The upper boundary condition can be formulated as the radiation condition (Baines 1995), assuming that the atmosphere is topped at zeta > zeta^sup 0^ by an inviscid (gamma = 0) homogeneous layer with constant temperature T = T^sup 0^ and wind U = U^sup 0^. In such a layer, U^sub zeta^, V^sub zeta^ = 0, rho, tau, beta = 0, v = constant, and lambda = lambda^sup 0^ = constant, and Eq. (3) simplifies to

…

The solution of this equation, satisfying the group energy upward propagation (i.e., the radiation) condition, is

omega ~ e^sup (-1/2-mu)zeta^, where (5)

…

The upper mu corresponds to an evanescent wave and the lower mu represents a free wave with kappa providing the proper phase sign and enabling upward propagation of the group energy. The top boundary condition is then a requirement that the solution in the top layer is presented by (5).

3. Discrete spectral buoyancy wave equation

Let us introduce a staggered vertical grid with full levels zeta^sub i^, and half levels zeta^sub i^sub +1/2^^,

0 = zeta^sub 1/2^ < zeta^sub 1^ < . . . < zeta^sub i-1^ < zeta^sub i-1/2^ < zeta^sub i^ < zeta^sub i+1/2^ < . . . < zeta^sub M^ < zeta^sub M+1/2^ < [infinity],

where M is the number of discrete layers. Layers are centered at zeta^sub i^ and layer boundaries are at levels zeta^sub i^sub +-1/ 2^^. The discrete omega field is located at half levels: omega^sub i^sub +1/2^^, whereas differences are located at full levels:

Deltaomega^sub i^ = omega^sub i+1/2^ – omega^sub i-1/2^.

Auxiliary functions rho, tau, alpha, beta, and lambda are considered analytical functions of zeta and thus are defined for each zeta^sub i^, zeta^sub i^sub +1/2^^ analytically.

The discrete approximation of (3) is L^sup +^^sub i+1/ 2^Deltaomega^sub i+1^ – L^sup -^^sub i+1/2^Deltaomega^sub i^ + Deltazeta^sup 2^^sub i+1/2^lambda^sub i+1/2^omega^sub i+1/2^ = 0, (6)

with coefficients

…

and

…

where Deltazeta^sub i^ = zeta^sub i+1/2^ – zeta^sub i-1/2^, Deltazeta^sub i+1/2^ = zeta^sub i+1^ – zeta^sub i^.

4. Solution factorization technique

Considering tentatively that the ith layer is homogeneous and inviscid, the solution in this layer is given by the exponential form (5), allowing one to define the decrease factor of omega in the layer as

…

where chi^sub i^ is the (complex) phase shift in the layer of unit depth. The solution can be presented at discrete half levels as a cumulative product of decrease factors

… (7)

with omega^sub 1/2^ as the surface value. In the general nonhomogeneous case, we will seek solution in the same form, although the layers are not homogeneous anymore [moreover, the presentation of (7) does not need such a restrictive precondition but supports instead an assumption of continuous and differentiable altering of reference atmosphere inside layers]. Loading (7) into wave equation [(6)] yields a two-point, nonlinear, nonhomogeneous recurrence for c^sub i^:

… (8)

The recurrence direction proves to be an essential property here. The proper solution omega decreases exponentially with height, which yields stable recurrence at the stepping from the top to the bottom, in the direction of decreasing i and increasing omega, and unstable for the opposite stepping direction. The start value C^sub M^ has to be specified from the top boundary condition, which can be formulated, first solving the homogeneous problem. Then, the remaining factor omega^sub 1/2^ in (7) can be specified from the bottom boundary condition [(4)], which becomes in the discrete case

… (9)

a. Special case of the homogeneous inviscid atmosphere

Homogeneous stratification in a discrete model assumes both a homogeneous background state and homogeneous layering:

… (10)

In this case, beta^sub 1+1/2^ = 0, while L+- and lambda in (6) become height-independent:

…

Considering the inviscid atmosphere, gamma = O, we will have v = k . U (i.e., lambda becomes a real constant). Assuming that in the case of height-independent coefficients, the solution of (8) also becomes height-independent,

… (11)

(8) reduces to a quadratic equation for c, identical at each ith level:

…

Solutions of this equation depend on the parameter

… (12)

Excluding an option Q < O, which is theoretically possible yet can appear at very low vertical resolutions Deltazeta ~1 only, the two options do remain, Q >/= 1 and O < Q < 1, yielding evanescent and free-wave radiative solutions of the decrease factor:

… (13)

… (14)

The solution [(7)] becomes in the case of decrease factor [(13)]

… (15)

presenting the discrete approximation to the exact solution in a homogeneous, inviscid atmosphere.

At sufficiently high vertical resolution, Deltazeta [much less than] 1 (which can always be considered a valid assumption), Q reduces approximately to

…

If, in addition,

… (16)

then (14) simplifies to

…

whereas the dicrete solution [(15)] converts to the exact solution of the homogeneous layer [(5)].

Using the solution of the homogeneous layer above, the “radiative” top initialization of the recurrence (8) in the nonhomogeneous case is provided by

…

with mu from (14) and Q from (12) as

…

where H^sub M^, N^sub M^, and Deltazeta^sub M^, and V^sub M^ correspond to the topmost layer M of the discrete model.

b. Requirements to vertical resolution

At high vertical resolution, if Deltazeta [much less than] 1, the decrease factors c^sub i^ are supposed to be, and L^sup +-^ are, close to unity in absolute terms, yielding a constriction of the free term in (8):

… (17)

[rather close to the condition (16)]. In the case of an inviscid atmosphere, lambda^sub j^ + 1/2 can become infinite for critical wave vectors k* (specific for each level i), for which v^sup 2^^sub i^ – f^sup 2^ = (k* . U^sub i^)^sup 2^ – f^sup 2^ = 0, yielding Deltazeta^sub i+1/2^ [arrow right] 0 in (17). Due to their discrete nature, the wave vectors are not (except in the case of the special choice of U^sub i^) strictly critical, but a lot of them can be nearly critical, (k* . U^sub i^)^sup 2^ – f^sup 2^ [asymptotically =] 0, which would cause very large lambda and unfeasibly small Deltazeta. Especially important for numerical solution accuracy are levels with evanescent wind U [arrow right] 0, as the corresponding critical wave vectors are located in the maximum area of spectral amplitudes. The levels near the evanescent wind U^sub i^ [arrow right] 0 are therefore called the critical layer. In the case of inviscid atmosphere, it is impossible to satisfy (17) for any computationally considerable size of Deltazeta in the critical layer. Fortunately, introduction even of a rather moderate turbulent friction would regularize lambda, as v^sup 2^^sub i^ – f^sup 2^ = (k . U^sub i^ – igammak^sup 2^)^sup 2^ -f^sup 2^ cannot become zero anywhere anymore.

Using in the vicinity of critical wave vectors an approximation … (17) gives

…

Minimization of the square root here with respect to … provides an estimation of the maximum vertical grid step

… (18)

which is independent of the horizontal wave vector k and thus presents a global estimate for all horizontal spectral modes. It is convenient to present the diffusion coefficient as

… (19)

The nondimensional parameter … has a simple meaning: … is the e-folding decrease period in units 1/N for the highest horizontally resolved gravity wave mode with scale ~ Deltax. In the discrete case, Deltax is the horizontal grid length; in the horizontally continuous formalism, Deltax can be estimated as the internal spatial scale of orography. From (18) and (19) we finally get

… (20)

As an example, Deltaz = HDeltazeta < 10 m for Deltax = 1 km, gamma^sup 0^ = 2 x 10^sup -2^, and f/N = 10^sup -2^. This vertical resolution limit decreases to l m for a 10-times higher horizontal resolution Deltax = 100 m.

Condition (20) should not be interpreted as an exact upper limit, but rather as a rough estimation of vertical resolution that is expected to provide the required solution accuracy. The actual vertical resolution, though based on estimation (20), must be established experimentally in every particular case. As an example, in the horizontally one-dimensional flow experiments with critical layers, the vertical grid step has to be taken up to 10 times smaller than estimation (20) inside of a critical layer and can be chosen several times larger then (20) far away from the layer. In many cases when the absolute value of the wind is large and there is moderate horizontal resolution (Deltax >/= 1 km), no friction is required at all (though the friction inclusion is actually not prohibited but rather wanted in bringing the model closer to reality).

5. Modeling examples

In the following examples, a horizontally discrete Fourier transform is applied to orography, which presents a “witch of Agnesi” profile,

… (21)

where the maximum height h^sub 0^, central coordinates x^sub 0^, y^sub 0^, and half-widths in the directions of coordinate axes a^sub x^, a^sub y^ are constant parameters. The spectral wave equation is solved using the above described approach, and the result is then inverted back to physical space. The Coriolis parameter is f = 10^sup -4^ s^sup -1^. Vertical velocity w = Homega/p is shown in all cases. The horizontal resolution in the following examples is chosen as the lowest possible while ensuring precision in all instances, which is controlled with the resolution-doubling method.

Figure 1 presents the wave pattern of w for the one-dimensional Agnesi ridge with h^sub 0^ = 100 m, a^sub x^ = 2 km, and a^sub y^ = [infinity]. Reference temperature T^sup 0^(p) presents a climatological profile. It is 280 K at the surface, has a lapse rate of 6.5 K km^sup -1^ in the troposphere, and is constantly 202 K in the stratosphere. The tropopause height is 12 km. The horizontal resolution is Deltax = 500 m and the grid domain in x-direction is 2048 points. Vertical resolution is chosen to be Deltaz = 100 m, (Deltazeta [asymptotically =] 0.01), and M = 300 levels, while the atmosphere is inviscid with gamma^sup 0^ = O at all levels. The control experiment shows that the vertical resolution increase and introduction of weak (gamma^sup 0^ = 0.01) viscosity does not alter modeling results. However, more strong friction with gamma^sup 0^ = 0.05 would damp the wave field moderately. Two wind profiles are applied. In Fig. 1a, the reference wind is unidirectional and uniform with U= 12 m s^sup -1^. In Fig. 1b, the wind U = 12 m s^sup – 1^ on the surface, has shear 0.25 m s^sup -1^ km^sup -1^ in the troposphere, and becomes constant U= 15 m s^sup -1^ in the stratosphere.

Though the buoyancy wave reflection on the tropopause and the tropospheric waveguide formation has theoretically been proven some time ago (Eliassen 1968), there has been little numerical experimentation, showing the details of the process. As seen in Fig. 1a, a reflected wave train forms already at shearless wind conditions. However, Fig. 1b demonstrates that a rather moderate wind shear will cause substantial wave reflection strengthening and wave train elongation. The wave train will increase in length along with the wind shear and can reach several thousands kilometers (depending on the turbulent friction intensity). The current examples have a special interest due to the wave train wiggling, which is observable at shearless casea and at weak shear but would disappear with further shear strengthening.

Figure 2 presents a flow over an Agnesi ridge with h0^sub ^ = 100 m, a^sub x^ = 3 km, and a^sub y^ = [infinity] and with the same temperature profile as in the previous case. However, the wind is backing with height in this model, having value 10 m s^sup -1^ at the surface and decreasing linearly with height. It becomes zero at 5- (Fig. 2a) and 2.5-km (Fig. 2b) levels, which represents the central heights of respective critical layers, and decreases with height further to a constant value -2.0 m s^sup -1^ at the 5.5- and 3-km heights, respectively. The horizontal resolution is 500 m and the number of horizontal grid points is 256. The vertical grid step decreases linearly with height from Deltaz = 100 m at the surface to Deltaz = 5 m at the wind reversal level, in the case of Fig. 2a, and from Deltaz = 50 m at the surface to Deltaz = 10 m, in the case of Fig. 2b. Above these levels, the vertical grid step is kept constant (i.e., 5 and 10 m, respectively). The number of vertical levels is 240 and turbulent viscosity gamma[degrees] = 0.05. The main aim of these examples is to demonstrate the need of enhanced resolution in the vicinity of a critical layer. The coincidence with the earlier reported results (Miranda and Valente 1997; Grubisic and Smolarkiewicz 1997; Shen and Lin 1999) is excellent, demonstrating complete absorption of orographic waves at the critical layer.

Figure 3 shows horizontal cross sections of w at different heights for a wind, sheared both in amplitude and direction, at the case of a circular hill with h^sub 0^ = 300 m, a^sub x^ = a^sub y^ = 3 km. The horizontal resolution is Delta^sub 0^ = 300 = 1.1 km. The horizontal grid consists of 256 x 256 points. The model has 300 levels with constant resolution Deltaz = 50 m in the vertical. The turbulent viscosity gamma[degrees] = 0.05. The reference atmosphere is isothermal with T[degrees] = 280 K. The wind amplitude is a hyperbolic function of height with a maximum of 40 m s^sup -1^ at z = 15 km, while the wind direction rotates uniformly with height:

u^sub x^ = U(z)cos(piz/z^sub rev^), u^sub y^ = U(z)sin(piz/z^sub rev^), (22)

with

…

where z^sub rev^ = 12 km, U^sub s^ = 10 m s^sup -1^ U^sub 0^ = 120 m s^sup -1^, and z^sub max^ = 30 km.

The solution appears to be insensitive to vertical resolution doubling, which means that a constant resolution Deltaz = 50 m is sufficient here. However, the decrease of gamma[degrees] to 0.01 implies a need for a vertical resolution increase to Deltaz = 20 m.

The coincidence of the present example wave pattern with the former numerical experiment by Shutts and Gadian (1999) in similar sophisticated wind profile conditions (variable shear with height plus a uniform directional shear) is excellent. There is no analytical solution available for the presented wind profile (22) yet, but the closest available analytical model (Shutts 2003) with constant shear both in height and direction exhibits quite close behavior. A typical feature of this kind of flow regime with uniform directional wind shear is the complete buoyancy wave absorption at the wind reversal height Z^sub rev^.

6. Conclusions

The described solution factorization method presents an adequate, simple, and fast orographie wave modeling tool in the case of rather sophisticated flow regimes both in two- and three-dimensional cases. Experiments with realistic height-dependent temperature, including the tropopause, and optional sheared winds are attainable. Rather large modeling domains in association with high horizontal and vertical resolutions can be used, which makes high-resolution modeling of extended wave fields possible. As an example, in the horizontally one-dimensional case, the horizontal domain of 5000-km lengths with 500-m horizontal resolution and 1000 levels in vertical with 10-m vertical resolution would be approximately a “1-min” task on a personal computer.

The requirement for vertical resolution (20) holds good, if the reference wind does not become evanescent at some height. As numerical experimentation shows, for large winds (U = 10 m s^sup – 1^ can be considered “large”) and for moderate horizontal resolutions (Deltaz >/= 1 km), the inviscid atmospheric model can even be used without the loss of stability and accuracy. However, wind evanescence at some level, associated with the formation of a critical layer around that height, will require quite high vertical resolutions (up to Delta ~ 1 m at horizontal resolutions Deltax ~ 100 m) in the vicinity and inside of the critical layer, outmatching resolution condition [(20)] about 5 times. The extremely high requirement for resolution turns the modeling of critical layer events into a most expensive and resource-demanding computational task.

Like all simplified models, the developed numerical scheme has its restrictions. The first restriction is common with all linear models-it cannot assess nonlinear effects like nonlinear waves, wave breaking, and blocking. Also, due to specifics of the numerical scheme, the model is not suited for linear wave study in conditions of horizontally inhomogeneous stratification. The main area of application of the developed numerical solution is high-resolution, high-precision modeling of linear orographie waves for arbitrary low orography in vertically complex atmospheric stratification conditions. Also, the model can be applied as a test tool for the numerical accuracy study of adiabatic kernels of the nonlinear nonhydrostatic NWP models.

Acknowledgments. This investigation has been supported by the Estonian Science Foundation under Research Grant 5711.

REFERENCES

Baines, P. G., 1995: Topographie Effects in Stratified Flows. Cambridge University Press, 482 pp.

Broutman, D., J. W. Rottman, and S. D. Eckermann, 2003: A simplified Fourier method for nonhydrostatic mountain waves. J. Atmos. Sci., 60, 2686-2696.

Eliassen, A., 1968: On the meso-scale mountain waves on the rotating Earth. Geofys. Publ, 27 (6), 1-15.

Grubisie, V., and P. K. Smolarkiewicz, 1997: The effect of critical levels on 3D orographie flows: Linear regime. J. Atmos. Sci., 54, 1943-1960.

Hazel, P., 1967: The effect of viscosity and heat conduction on internal gravity waves at a critical level. J. Fluid Mech., 30, 775- 783.

Jones, W. L., 1967: Propagation of internal gravity waves in fluids with shear flow and rotation. J. Fluid Mech., 30, 439-488.

Lin, C. C, 1955: The Theory of Hydrodynamic Stability. Cambridge University Press, 155 pp.

Miller, M. J., 1974: On the use of pressure as vertical co- ordinate in modelling convection. Quart. J. Roy. Meteor. Soc., 100, 155-162.

_____, and R. P. Pearce, 1974: A three-dimensional primitive equation model of cumulonimbus convection. Quart. J. Roy. Meteor. Soc., 100, 133-154.

_____, and A. A. White, 1984: On the nonhydrostatic equations in pressure and sigma coordinates. Quart. J. Roy. Meteor. Soc., 110, 515-533.

Miranda, P. M. A., and M. A. Valente, 1997: Critical level resonance in three-dimensional flow past isolated mountains. J. Atmos. Sci., 54, 1574-1588.

Queney, P., 1948: The problem of airflow over mountains. A summary of theoretical studies. Bull. Amer. Meteor. Soc., 29, 16- 26.

Room, R., 1998: Acoustic filtering in nonhydrostatic pressurecoordinate dynamics: A variational approach. J. Atmos. Sci., SS, 654-668.

_____, and A. Mannik, 1999: Response of different nonhydrostatic, pressure-coordinate models to orographie forcing. J. Atmos. Sci., 56, 2553-2578.

Scorer, R. S., 1949: Theory of waves in the lee of mountains. Quart. J. Roy. Meteor. Soc., 75, 41-56.

Shen, B.-W., and Y.-L. Lin, 1999: Effects of critical levels on two-dimensional back-sheared flow over an isolated mountain ridge. J. Atmos. Sci., 56, 3286-3302.

Shutts, G. J., 1998: Stationary gravity wave structure in flows with directional wind shear. Quart. J. Roy. Meteor. Soc., 124, 1421- 1442.

_____, 2003: Inertia-gravity wave and neutral Eady wave trains forced by directionally sheared flow over isolated hills. J. Atmos. Sci., 60, 593-606.

_____, and A. Gadian, 1999: Numerical simulations of orographic gravity waves in flows which back with height. Quart. J. Roy. Meteor. Soc., 125, 2743-2765.

Smith, R. B., 1980: Linear theory of stratified hydrostatic flows past an isolated mountain. Tellus, 32, 348-364.

_____, S. Skubis, J. D. Doyle, A. Broad, C. Kiemle, and H. Volkert, 2002: Mountain waves over Mont Blanc: Influence of a stagnant boundary layer. J. Atmos. Sci., 59, 2073-2092.

White, A. A., 1989: An extended version of a nonhydrostatic, pressure coordinate model. Quart. J. Roy. Meteor. Soc., 115, 1243- 1251.

REIN ROOM AND MARKO ZIRK

Tartu University, Tartu, Estonia

(Manuscript received 15 September 2006, in final form 29 November 2006)

Corresponding author address: Rein Room, Institute of Environmental Physics, Tartu University, Ulikooli 18, Tartu 50090, Estonia.

E-mail: rein.room@ut.ee

Copyright American Meteorological Society Oct 2007

(c) 2007 Monthly Weather Review. Provided by ProQuest Information and Learning. All rights Reserved.