## Abstract

An attempt is made to define the thermodynamics of internal gravity waves breaking in the middle atmosphere on the basis of the energy conservation law for finite fluid volumes. Consistent with established turbulence theory, this method ultimately determines the turbulent dissipation to be equivalent to the frictional heating owing to the Reynolds stress tensor. The dynamic heating due to nonconservative wave propagation, that is, the energy deposition, consists of two terms: namely, convergence of the wave pressure flux and a residual work term that is due to the wave momentum flux and the vertical shear of the mean flow. Only if both heating terms are taken into account does the energy deposition vanish, by definition, for conservative quasi-linear wave propagation. The present form of the energy deposition can also be deduced from earlier studies of Hines and Reddy and from Lindzen.

The role of the axiomatically defined heating rates in the heat budget of the middle atmosphere is estimated by numerical experiments using a simple general circulation model (SGCM). The authors employ the theory of Lindzen to parameterize saturation of internal gravity waves, including the first theorem of Eliassen and Palm, to define the wave pressure flux. It is found that in the climatological zonal mean, the residual work and the simulated turbulent dissipation give maximum cooling/heating rates of −5 K day^{−1} and +2.5 K day^{−1} in the summer mesosphere/lower thermosphere. These values are small but not negligible against the major contributions to the heat budget.

Filtering out gravity wave perturbations in the thermodynamic equation reveals that the simulated dissipation, which is due to the shear of the planetary-scale flow, does not generally represent the total dissipation. The latter turns out to be dominated by the shear associated with gravity wave perturbations. Taking advantage of Lindzen, the total frictional heating can be calculated and is quantitatively consistent with in situ measurements of Lübken.

## 1. Introduction

Since the work of Lindzen (1981), there is quantitative evidence that the breakdown of internal gravity waves (IGWs) in the middle atmosphere accounts for the strong zonal drag that balances the Coriolis force associated with the summer-to-winter pole meridional circulation. To a first approximation, the reversed meridional temperature gradient in the upper mesosphere/lower thermosphere (MLT) can be interpreted as a dynamical consequence of the IGW-induced meridional circulation that is associated with adiabatic cooling/heating in the summer/winter MLT. However, it has been noted in several studies that the direct heating rates associated with IGW breakdown such as energy deposition, turbulent dissipation, and turbulent diffusion give important contributions to the heat budget of the MLT as well (e.g., Lübken 1992; 1997a,b; Luo et al. 1995). Several concepts to parameterize these heating terms can be found in the literature (e.g., Schoeberl et al. 1983; Strobel et al. 1985; Akmaev 1994; Fritts and Luo 1995; Hines 1997). Corresponding numerical experiments are particularly important in order to understand the thermal structure of the MLT.

Parameterizations of turbulence and wave disturbances yield contributions to the momentum and thermodynamic equations of motion of the resolved mean flow. A model so-defined is tuned by varying the free parameters so as to achieve agreement with observations. However, having a huge number of free parameters is not necessarily useful, especially in cases when their physical interpretations are not clear. It is, rather, desirable to take advantage of physically motivated constraints in order to reduce the level of ambiguity of a given parameterization scheme. For example, the energy conservation law applied to the mean flow can be used as such a constraint. Strictly speaking, for any given subscale process, this law ultimately determines the thermodynamic effects in terms of the momentum flux tensor and the sensible heat flux. Thus, the benefit is not only a minimum number of free parameters but in addition an energetically consistent formulation of the equations of motion.

With regard to IGW parameterization schemes employed in general circulation models (GCMs), it appears worthwhile to consider possibilities to formulate thermodynamic effects like turbulent dissipation and energy deposition in a most straightforward manner using the energy conservation law. In the present study, a corresponding attempt is outlined. While the planetary-scale flow is explicitly resolved, IGWs and turbulence are treated as unresolved disturbances to be parameterized. Furthermore, a strict separation between the turbulent scale (several tens of meters in the vertical) and the scale of IGWs (typically a few hundred kilometers in the horizontal and several hundred meters in the vertical) is implicitly assumed. As pointed out by McIntyre (1989), it is open to question whether such a scale separation can be justified by available observational data. On the other hand, the assumption is inherent in almost every present-day IGW parameterization (e.g., Matsuno 1982; Lindzen 1981; Holton 1982; Gardner 1996; Hines 1997). In particular, while turbulence is commonly parameterized in terms of momentum and heat diffusion, the vertical fluxes associated with IGW propagation are nondiffusive. The importance of the thermodynamic heating rates as defined by the energy conservation law will be estimated by numerical experiments. In this study, a simple GCM (SGCM) ranging from the surface to the lower thermosphere with IGW saturation being parameterized following Lindzen (1981) and Holton (1982) is employed.

Any parameterization of subscale processes implies certain approximations to be fulfilled by the perturbation flow. The complete set of approximations can be identified by comparing the assumed model equations with the full equations of motion after subscale perturbations have been filtered out. For the present concept of gravity wave parameterization in association with the energy conservation law, a corresponding analysis is given in the appendix. This consideration also sheds some light on a specific question concerning the heat budget of the summer MLT, namely, what process balances the strong turbulent dissipation inferred from in situ measurements (Lübken 1997a).

The remainder of the paper is organized as follows. Section 2 recapitulates how the energy conservation law can be used to axiomatically define the dissipation in turbulent flow. In addition, the energy deposition analog associated with IGW propagation is derived and further specified for a family of gravity waves. A corresponding idealized GCM is defined in section 3 and some numerical experiments are presented in section 4. In particular, section 4c presents our model estimate of the total frictional heating that is implicitly—rather than explicitly—accounted for in the model. Section 5 gives a concluding discussion.

## 2. General remarks on parameterization of turbulence and gravity waves in global circulation models

Let us consider the horizontal momentum equation and the thermodynamic equation of motion for the planetary-scale atmospheric flow with turbulence and IGWs being represented by parameterizations:

The three-dimensional velocity field and nabla operator are denoted by **v**_{3} and ∇_{3} while **v** and ∇ consist of horizontal components only. The symbol *z* represents the height above sea level and **e**_{z} is the unit vector in vertical direction. The notation is otherwise standard except those abbreviations explained below.

The thermodyamic equation (2) is written in terms of sensible heat *h* = *c*_{p}*T.* The sum of heating by solar insolation, cooling by longwave radiation and latent heating is represented by *c*_{p}*Q.* The Reynolds stress tensor and the heat flux due to turbulence are denoted by Σ_{turb} and **J**. These terms are usually parameterized by diffusion. A common representation for **J** is

where potential temperature, Prandtl number, and horizontal and vertical diffusion coefficient are represented by Θ, *P*_{r}, *K,* and *ϑ,* respectively. In the shallow atmosphere approximation we get from (2.3)

In the design of GCMs, it is usual to parameterize turbulent friction rather than turbulent stress. As a result, conventional formulations of horizontal momentum diffusion violate the symmetry property of Σ_{turb}. This difficulty has been discussed in Becker (2001) and a consistent stress tensor formulation has been proposed that may be written as

In (2.5), the tensor product is indicated by an open circle and the transposed tensor by the superscript T. Further, *a*_{e} is the radius of the earth. In the shallow atmosphere approximation, the corresponding friction force yields

In analogy to the standard vertical momentum diffusion contained in (5) and (6), the stress tensor due to IGW propagation can be written as

where **F** is the gravity waves' averaged vertical flux of horizontal momentum. This leads to the standard momentum deposition

The remaining terms *ε* and *E* on the rhs of (2) represent turbulent dissipation (e.g., Hinze 1959) and energy deposition, respectively. Our model equations (1)–(8) are completed by the hydrostatic approximation, the continuity equation, and boundary conditions. Furthermore, both *ε* and *E* need to be specified by means of parameterization. Alternatively, we require that the energy conservation law (e.g., Lindzen 1990; Serrin 1959) applies for the planetary-scale flow, hence

In (2.9), *G* denotes an arbitrary fluid volume, *e* = *c*_{υ}*T* is internal energy, **I** is the idemfactor, *w* is vertical velocity, and *g* is the gravity acceleration. The symbol *F*_{p} denotes the vertical flux of sensible heat due to the gravity waves' pressure perturbations [see also Eq. (A.14)]. In analogy to the definition of molecular dissipation, it follows for turbulence

and for gravity waves

We emphasize that (10) is a well-known approximation valid for small-scale turbulence (e.g., van Mieghem 1973; Smagorinsky 1993). Since (Σ_{turb}**∇**_{3}) · **v**_{3} represents the transfer between kinetic energy of the mean flow and turbulent kinetic energy, Eq. (10) simply states that this kinetic energy transfer is downscale. That is, Eq. (10) implies a cascade of kinetic energy that is maintained by (Σ_{turb}**∇**_{3}) · **v**_{3} at the largest turbulent scale and dissipation by molecular friction at the Kolmogorov scale. This is quite consistent with the established picture of turbulent dissipation. Note that any parameterization of Σ_{turb} must guarantee that the dissipation *ε* is positive definite in order to ensure consistency with the second law of thermodynamics.

As long as Σ_{turb} is not defined, the dissipation cannot consistently be accounted for. Therefore, GCMs do usually not employ (10). Alternatively, we get from (5) (Becker 2001):

The frictional heating *ε*_{z} due to vertical momentum diffusion is well known in the context of boundary layer mixing (e.g., van Mieghem 1973). It represents turbulent dissipation generated by gravity wave breakdown in the MLT if we identify *ϑ* with the corresponding turbulent diffusion coefficient. From observational estimates (Lübken 1997a; Hocking 1999) we can expect that this heating rate plays an important role in the climatology of the MLT, though not commonly accounted for in present-day IGW parameterizations. Horizontal diffusion terms are recapitulated here for the sake of completeness only. Indeed, in the region of gravity wave breaking, the vertical terms represent the major contributions from turbulent diffusion.

The energy deposition (11) owing to IGWs consists of two terms, namely, convergence of the wave pressure flux and a residual work *W*_{res} that is directly due to the wave momentum flux **F** and the shear of the mean flow. This form of the energy deposition may also be infered from Hines and Reddy (1967), Lindzen (1973, 1990, chapter 8), or from Hines (1997, section 9). Our assumptions made so far to define wave–mean flow interaction and dissipation can be identified with corresponding approximations applied to the equations of motion with gravity wave perturbations being filtered out. A corresponding analysis is given the appendix.

For convenience, we assume a family of individual gravity waves. Hence, the momentum flux may be written as

Unit vectors in zonal and meridional direction are abbreviated by **e**_{x} and **e**_{y}. For each wave, **e**_{αj} gives the direction of horizontal phase propagation and *c*_{j} > 0 is the corresponding phase speed. Thus, *j*^{−1}_{0}*F*_{j}**e**_{αj} is the vertical flux of horizontal momentum generated by an IGW with horizontal phase velocity *c*_{j}**e**_{αj}. Owing to the first theorem of Eliassen and Palm valid for linear wave propagation [e.g., Lindzen 1990, Eq. (8.14)], the pressure flux can be expressed in terms of the momentum fluxes, horizontal phase speeds, and the background horizontal wind:

According to the second theorem of Eliassen and Palm [e.g., Lindzen 1990, Eq. (8.24)], the vertical momentum flux is constant with height as long as diabatic processes are absent. It is evident that also the energy deposition must vanish for conservative wave propagation. Obviously, this holds neither for *W*_{res}, nor for the pressure flux convergence −*ρ*^{−1}∂_{z}*F*_{p}. Both terms can be different from zero regardless of whether the momentum flux **F** changes with height or not. This paradox has already been noted in Hines and Reddy (1967) and Lindzen (1973). Owing to (11), the complete energy deposition yields

and obviously vanishes for conservative wave propagation. It is also interesting to look at the signs of *W*_{res} and *E.* For instance, let us consider a single saturated gravity wave with westward phase velocity in the winter MLT. Because of ∂_{z}*u* < 0 for the background zonal wind we have *W*_{res} = −*ρ*^{−1}**F** · ∂_{z}**v** < 0. Nevertheless, (15) is essentially positive, as it should be: the wave deposits its zonal momentum (∂_{z}*F*_{j} < 0) while **v** · **e**_{αj} − *c*_{j} < 0 below the critical level, hence *E* > 0. An analogous consideration is valid for the summer MLT if we consider a single saturated gravity wave with eastward phase velocity and ∂_{z}*u* > 0 for the background zonal wind.

To further complete our simple gravity wave representation, it is assumed that each wave gives an independent contribution to the vertical diffusion coefficient *ϑ.* Thus, we write analogously to (13)

The first term on the rhs of (16) represents the IGW-induced turbulent diffusion coefficient (typically a few hundred square meters per second). With respect to the definition of the model experiments presented below, (16) also includes a weak background diffusion *ν*_{bg} chosen as in Becker and Schmitz (1999, their Fig. 3). In the MLT, the additional term is negligible against the IGW-induced diffusion coefficient.

In summary, given the momentum flux (13) and the turbulent diffusion coefficient (16), the complete set of direct heating rates accompanying IGW saturation is defined via (11), (12c), (14) and the second term on the rhs of (4) with the Prandtl number *P*_{r} as additional free parameters only.

We note that the dissipation (12c) and the energy deposition (15) bear some resemblence to other existing formulations briefly reviewed by McIntyre (1989, section 1). For example, turbulent dissipation has often been written as *ε*_{z} = *ϑN*^{2}/*β*_{1} [see Eq. (1) in the paper of McIntyre]. With *N*^{2} = *g*Θ^{−1}∂_{z}Θ, equivalence to (12c) is achieved if the free parameter *β*_{1} is equal to the Richardson number *N*^{2}(∂_{z}**v**)^{−2} [e.g., Lübken 1992, Eq. (21)]. A previous formulation of the energy deposition yields *E* = *ϑN*^{2}/*β*_{2} [see Eq. (3) in the paper of McIntyre]. This expression is consistent with (15) under the following conditions: *β*_{2} = 1, only a single gravity wave is considered, and the turbulent diffusion coefficient is calculated according to Lindzen (1981). An expression equivalent to the energy deposition (15) has also been employed by Akmaev (1994) in a global circulation model.

## 3. Model description

We employ the Kühlungsborn mechanistic general circulation model (KMCM) previously used in Becker and Schmitz (1999, 2001). The dynamical core of the model is a standard GCM—T29 spectral resolution, 60 hybrid levels extending from 990 up to 0.0003 mb. Turbulence and gravity waves are accounted for according to the preceding section. Radiative heating is represented by temperature relaxation toward a zonally symmetric equilibrium temperature *T*_{E} (as in Fig. 3a) that in the middle atmosphere resembles the radiatively determined state (Shine 1987). In the troposphere, *T*_{E} resembles a typical boreal winter state with increased meridional gradients (see Fig. 1a of Becker and Schmitz 2001). The relaxation time and the horizontal diffusion coefficient *K* are prescribed as functions of hybrid coordinate *η.* They are plotted in Fig. 1 as functions of the equivalent pressure defined as *η* × 1013 mb.

Horizontal momentum diffusion is applied on the model's hybrid surfaces using its definition derived in *z* coordinates. This may cause transformation errors in regions of steep orography. Therefore *K* is set equal zero below *η* = 0.8, and it is continuously raised to its midtropospheric value of about 5.5 × 10^{4} m^{2} s^{−1} between *η* = 0.8 and *η* = 0.6. In the upper troposphere, wavenumber 29 of horizontal vorticity is damped within about 9 days. The minimum damping time constant slopes down to about 0.6 h in the upper mesosphere. A Prandtl number of 2 is chosen for both horizontal diffusion and vertical diffusion owing to IGW saturation.

Cumulus heating in the deep Tropics and latent heating in the regions of the warm ocean currents are mimicked via prescribed heat sources and self-induced heating, respectively. Boundary layer mixing is based on Prandtl's mixing length concept in association with the Richardson criterion. The model furthermore includes a standard envelope orography. A detailed description of the model formulation with respect to the troposphere is given in Becker and Schmitz (2001).

Momentum flux and turbulent diffusion coefficient owing to gravity waves are parameterized analogously to Holton (1982). However, in the present 3D model experiments, IGW propagation is considered in eight equally spaced azimuths (*j*_{0} = 8), that is, *α*_{j} := (*j* − 1)*π*/4 for *j* = 1, … , 8. Assuming that horizontal variations of the resolved flow are negligible on the scale of IGWs, the IGW parameterization scheme is applied at each horizontal model grid point. Furthermore, the specification of IGW parameters is horizontally uniform. Individual breaking levels *z*_{bj} and critical levels *z*_{cj} are implicitly defined via

where *H* is a scale height and *ũ*_{j} represents the initial amplitude of gravity wave *j.* Between *z*_{bj} and *z*_{cj}, the individual wave contributions *F*_{j} and *D*_{j} follow directly from Lindzen (1981). Wave obliteration is complete above *z*_{cj} while a smooth attack of *D*_{j} is assumed between *z*_{0} and *z*_{bj}. The resulting formulas may be written as

The initial amplitudes *ũ*_{j}, the phase speeds *c*_{j}, and the horizontal wavenumbers *k*_{j} are assumed to have elliptical distributions in the horizontal plane:

The level *z*_{0} of gravity wave excitation is located in the upper troposphere (at approximately 180 mb), and the present choice of the scale height is *H* = 7.5 km. The remaining free parameters are *k*_{0} = 2*π*/340 km, *s*_{k} = 0.1, *c*_{0} = 10 m s^{−1}, *s*_{c} = 0.65, *ũ*_{0} = 0.91 m s^{−1}, and *s*_{u} = −0.16. Figure 2 illustrates the resulting azimuthal distributions of horizontal wavenumbers, phase speeds, and amplitudes. The present choice of parameters is consistent with Lindzen's estimate that gravity waves with eastward phase speeds of a few tens of meters per second and those with phase speeds close to zero are of primary importance. In particular, due to the climatological background wind, the former propagate into the summer MLT where they break at higher altitudes and more efficiently in comparison to the latter, which dominate gravity wave breaking in the winter MLT. This asymmetry between summer and winter will also become evident from computational results presented below.

## 4. Numerical experiments

In the following we inspect some perpetual January simulations using the SGCM defined in the preceding sections. In the troposphere, full forcing of planetary waves is applied in the way described in Becker and Schmitz (2001). Temporal averages are computed on the basis of 900-day time series using a sampling rate of 1 day. In sections 4a–c we will consider the model climatology and the heat budget in the standard run where all thermodynamic effects owing to tubulence and IGWs are accounted for. In section 4d, some sensitivity experiments are presented to estimate the importance of the residual work *W*_{res} [Eq. (11)] and the dissipation rates *ε*_{z} and *ε*_{h} [Eq. (12)]. Furthermore, the model sensitivity to the choice of the Prandtl number is assessed.

### a. Climatology

Figures 3 and 4 illustrate the zonal mean model climatology in the standard run. The simulated thermal structure of the middle atmosphere (see Fig. 3a) captures the prominent substantial differences between the climatological temperature (black contours) and the assumed equilibrium temperature *T*_{E} (white contours). The model yields a reasonable summer mesopause with a minimum temperature of less than 140 K above the pole as well as a realistic boreal winter stratopause. Note that both features are not in any way captured by the structure of *T*_{E} but result from the nonlinear dynamics of gravity waves, planetary waves, and the accompanying diffusion processes. The simulated zonal-mean zonal wind (Fig. 3b) has reasonable magnitude. The westerly wind maximum reproduces the prominent vertical tilt toward the North Pole in the stratosphere at low altitudes. Also consistent with observations, the tilt is reversed at higher altitudes. The residual mass streamfunction (Fig. 3d) is calculated in quasigeostrophic approximation. In the stratosphere, it is consistent with observational analyses (Rosenlof 1995). Comparison with the Eulerian mass streamfunction (Fig. 3c) proves that in the MLT, the residual circulation is driven by gravity waves rather than planetary waves. Furthermore, as a result of gravity waves and the vertical profile of the horizontal diffusion coefficient (Fig. 1b), the meridional fluxes of momentum and heat owing to transient and stationary planetary waves are damped out in the upper mesosphere (Fig. 4). The stationary wave structure and the longitudinal distribution of transients (both not shown) are quite reasonable, as in our prior model version (Becker and Schmitz 2001).

### b. Heat budget

The climatological zonal-mean heat budget is displayed in Fig. 5. First of all, we note the almost perfect balance between adiabatic heating plus advection by the resolved flow (Fig. 5f) and the total heating (Fig. 5e), which we take as the sum of temperature relaxation, dissipation, energy deposition, and heat diffusion. The overall thermal driving by temperature relaxation (Fig. 5a) is consistent with other estimates for the net radiative forcing (e.g., Gille and Lyjak 1986; Berger and von Zahn 1999). Figures 5a–d reveal that in the MLT, the total heating is a residual of several contributions to heat budget. Particularly in the region of gravity wave breaking, diffusion of heat, and energy deposition (Figs. 5c,d) are of the same order of magnitude as temperature relaxation.

The energy deposition (Fig. 5c) resembles the observed asymmetry between the summer and the winter mesosphere as deduced from rocketborne measurements (Lübken 1997a,b), hence reflecting corresponding hemispheric differences in gravity wave saturation as anticipated by Lindzen (1981). However, the energy deposition does not represent the turbulent dissipation, even though such an interpretation is sometimes suggested. The simulated dissipation is displayed in Fig. 5b (using a contour interval of 0.5 K day^{−1}, whereas the contour interval is 3 K day^{−1} in all other panels of Fig. 5). This simulated dissipation also reproduces the expected asymmetry between summer and winter. The maximum dissipation is about 2.5 K day^{−1}. This value is consistent for instance with Fritts and Luo (1995, see also references therein), who noted that a typical value of 2 K day^{−1} for the summer MLT lies within the wide range of experimental and theoretical estimates. Nevertheless, the simulated dissipation is too weak by almost one order of magnitude in comparison with the more recent observational results of Lübken (1997a). This apparent discrepancy is analyzed in the appendix and further addressed in the the next section.

As discussed in section 2, the energy deposition consists of convergence of the wave pressure flux plus a residual work that is negative in the MLT. Further, the dissipation or the heat diffusion consists of two terms, each of which is due to vertical and horizontal diffusion. Figure 6 depicts these individual parts of the heat budget. It is confirmed that frictional heating by IGW-induced vertical momentum diffusion dominates the dissipation in the summer MLT (Figs. 6a,b). The overall energy deposition is mainly due to convergence of the pressure flux (Fig. 6d). Nevertheless, a maximum cooling rate of about −5 K day^{−1} due to *W*_{res} (Fig. 6c) is not negligible. Comparing Fig. 5b to 5d or Figs. 6a,b to 6e,f reveals that vertical diffusion of potential temperature dominates the overall thermal effects from turbulence that are explicitly accounted for in the model. In the present simulation we obtain maximum cooling rates of about −30 K day^{−1} in summer and −15 K day^{−1} in winter (Fig. 6e). These values may crucially depend upon the assumed Prandtl number (see section 4d). Horizontal diffusion of temperature gives contributions of −6 to +12 K day^{−1} in the summer MLT (Fig. 6f).

### c. Model estimates of the turbulent dissipation

Obviously, there is a quantitative discrepancy between the simulated dissipative heating in the extratropical MLT (Figs. 5b, 6a, 6b) and corresponding state-of-the-art in situ measurements. One is tempted to attribute this result to the various uncertancies of the model and to the parameterization of gravity waves in general. Or one may suspect that the observational result is not representative of the climatological zonal mean. However, the analysis presented in the appendix offers another explanation.

The simulated dissipation owing to vertical momentum diffusion (12c) includes the shear of the mean flow only. On the other hand, the total frictional heating, which corresponds to the measured quantity, includes an additional dissipation that is directly due to the vertical wind shear of the gravity waves. Using Lindzen's saturation hypothesis, this additional dissipation can be estimated, yielding *ϑN*^{2}/2. Figure 7a shows the corresponding model estimate of the total dissipation *ε*_{h} + *ε*_{z} + *ϑN*^{2}/2 [Eq. (A.11)]. Comparison with Figs. 6a and 6b yields that the gravity wave wind shear gives the dominant contribution, hence reflecting the IGW-induced vertical diffusion coefficient *ϑ* shown in Fig. 7b. The model estimate of the total dissipation is quantitatively consistent with results of Lübken (1997a,b) for both the summer and the winter hemispheres. It seems worthwhile to mention that this agreement of computational and observational result is obtained without adjusting any additional free parameter. Note also that the additional dissipation *ϑN*^{2}/2 does not enter the thermodynamic equation of the model because it is overcompensated by the adiabatic term (A.16).

### d. Sensitivity experiments

A long-term perturbation run has been performed with the residual work *W*_{res} [see Eq. (11)] being neglected. The model setup has otherwise been as in the standard run described above. Figure 8a shows the difference of the climatological zonal-mean temperature in the perturbation run from the standard run. Due to the moderate cooling rate owing to *W*_{res} (Fig. 6c), there is a significant increase of mesospheric temperatures of 3–10 K. A second perturbation run has been performed with the simulated dissipation due to vertical momentum diffusion *ε*_{z} [Eq. (12c)] discarded in the thermodynamic equation of the model. Figure 8b shows the corresponding climatological temperature difference from the standard run. Compared to the previous sensitivity experiment, the signal has reversed sign and is less intense. This finding is consistent with the fact that *ε*_{z} is positive definite and weaker than *W*_{res}. Either temperature signal in Fig. 8 is located at higher altitudes in comparison to the maximum cooling/heating rates owing to *W*_{res} or *ε*_{z}.

In a third perturbation run, the simulated dissipation owing to both vertical and horizontal momentum diffusion *ε*_{h} + *ε*_{z} [Eq. (12)] has been neglected. Now the model response is much stronger compared to the previous cases, with maximum temperature and wind signals of −18 K and ±5, … , 10 m s^{−1} (Fig. 9). It is beyond the scope of the present study to interpret these model sensitivities in detail. Rather, we want to demonstrate that neglecting any seemingly minor contribution to the heat budget can have significant or even strong consequences in the model climatology of the MLT.

One of the most crucial parameters in computational dynamics of the MLT is the Prandtl number for vertical diffusion of potential temperature [second term on the rhs of Eq. (4)]. Owing to existing estimates (Hocking 1999), the Prandtl number should be greater than 1. For our standard run we chose a value of 2. In order to assess the model sensitivity to the choice of this free parameter, another long-term integration has been performed, employing a Prandtl number of 5 for vertical diffusion in the MLT. All other model parameters have been maintained (in particular, a Prandtl number of 2 for horizontal diffusion throughout the model domain and a Prandtl number of 1 in the planetary boundary layer). Figure 10 shows the MLT temperatures and residual meridional winds in both the standard run (lhs panels) and the perturbation run (rhs panels). Figures 10a,b indicate that due to a larger Prandtl number, the summer mesopause is heated up by about 20° and shifted to lower altitudes. Furthermore, its vertical slope with latitude (e.g., Berger and von Zahn 1999) flattens, and the absolute temperature minimum sharpens. Warmer temperatures with an averaged signal of about 10° are also observed throughout the winter mesosphere. These temperature changes do in the main reflect lower diabatic cooling by vertical diffusion of potential temperature. However, Figs. 10c,d reveal that the gravity wave–driven residual circulation in the MLT is also significantly reduced, particularly in summer, indicating that gravity effects become generally weaker if we simply increase the Prandtl number.

## 5. Concluding discussion

This study has recapitulated that the dynamic heating of the planetary-scale flow due to energy deposition by internal gravity waves consists of two terms (Hines and Reddy 1967; Lindzen 1973). In addition to convergence of the wave pressure flux, a second term follows from the energy conservation law and represents a residual work *W*_{res} that is directly due to the wave momentum flux and the vertical shear of the mean flow. For a discrete spectrum of linearly propagating IGWs, the wave pressure flux can be written according to Eliassen and Palm's first theorem, and the complete energy deposition is ultimately defined in terms of the momentum fluxes and phase speeds of the gravity waves. It can readily be shown that there is heating in the region of wave obliteration above the mesospheric jet. No energy deposition occurs in the case of conservative wave propagation.

In the literature, a clear distinction between energy deposition and turbulent dissipation is often lacking. However, it has been emphasized by Lindzen (1981, 1984) that turbulent diffusion and convergence of the vertical wave fluxes are distinct processes associated with IGWs breaking in the middle atmosphere. Even though momentum and energy deposition occur as a consequence of a diabatic process—for example, damping of the wave amplitudes by turbulent diffusion—they represent wave–mean flow interactions and must be regarded as dynamic processes. In a global circulation model with energetically self-consistent formulation of the equations of motion, the turbulent dissipation owing to gravity wave breakdown corresponds to the frictional heating associated with wave-induced turbulent momentum diffusion.

Summarizing these conceptional considerations, once momentum fluxes, phase speeds, and turbulent viscosity for a family of gravity waves are given, the direct heating rates from IGW saturation can well be defined without additional free parameters (except the Prandtl number) on the basis of elementary contraints. This method has been tested in numerical experiments using a considerably realistic SGCM for the lower and middle atmosphere. It has been found that in the climatological zonal mean, the residual work or the turbulent dissipation amount to maximum cooling/heating rates of −5 or +2.5 K day^{−1}. These values are small against the maximum radiative or adiabatic heating rates, but they are generally not negligible. Sensitivity experiments reveal that a neglect of the residual work or the simulated frictional heating owing to vertical momentum diffusion can significantly affect the model climatology of the MLT (Fig. 8). On the other hand, the choice of the Prandtl number for wave-induced turbulent diffusion appears to be particularly crucial to the simulated general circulation of the MLT (Fig. 10).

In the appendix we have outlined all approximations implicitly assumed by the formalism presented in section 2, and particularly by application of the energy conservation law for the mean flow. There is one issue to be mentioned that is important in any parameterization of momentum flux, pressure flux, and turbulent diffusion due to IGWs: The frictional heating owing to turbulent momentum diffusion is calculated from the resolved planetary-scale flow only rather than taking the shear due to IGWs into account as well. This leads to an underestimation of the simulated dissipation rate. The analytical consideration proves that the model formulation is nevertheless appropriate since frictional heating owing to gravity wave shear is overcompensated by adiabatic cooling generated by the waves themselves. Furthermore, employing Lindzen's gravity wave saturation theory, the total dissipation can be estimated [Eq. (A.11)], and the corresponding model result (Fig. 7a) describes in situ measurements of Lübken (1997a,b) in a satisfactory way.

With respect to the topic of the present study, the chief limitation of the idealized GCM lies in the quasi-linear parameterization of gravity wave saturation. No direct interaction among the gravity waves has been considered, like in the study of Holton and Zhu (1984), for instance. It is furthermore questionable to assume that turbulent diffusion is independent of the phase of the breaking gravity wave (Hodges 1967; Chao and Schoeberl 1984). As noted by Weinstock (1988), supersaturation may have important effects as well, but it has been neglected here. Other oversimplifications of the present model may be due to using a poor gravity wave spectrum and application of the first theorem of Eliassen and Palm to define the wave pressure flux. Despite these idealizations we believe that Lindzen (1981) presented the simplest physical model that captures the essential dynamics of gravity wave saturation. In this study we have tried to stay as close as possible to this concept in order to assess the direct thermodynamic effects in the framework of a global circulation model. This approach may not at least be justified by the fact that Lindzen's theory is capable of predicting quantitatively realistic frictional heating rates, once the wave parameters have been adjusted to yield the appropriate momentum deposition. The applied definitions of energy deposition and turbulent dissipation may also be useful in combination with other gravity wave parameterization schemes.

## Acknowledgments

For valuable discussions we would like to thank R. S. Lindzen, G. E. Thomas, E. M. Volodin, U. Achatz, and M. Rapp. Also, the comments of three anonymous reviewers are gratefully acknowledged.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

### APPENDIX

#### Filtering Out Gravity Wave Perturbations

Adapting the notation introduced in section 2 and neglecting the Coriolis force, we start with the three-dimensional momentum equation, mass conservation, and the thermodynamic equation. Turbulence is assumed to be filtered out and parameterized by diffusion, while IGW perturbations are considered as part of the resolved flow:

The turbulent dissipation is defined, as in section 2, as the frictional heating owing Σ_{turb}. Now, gravity wave perturbations are filtered out using ordinary and density weighted averages of any quantity *X*:

The average in (A.4) is assumed to extend over the typical horizontal and temporal scales of the planetary flow (i.e., a few thousands of kilometers and a few days). Owing to (A.4), (A.1)–(A.3) are transformed to

In order to achieve consistency of the two momentum equations (A.5) and (1), we must assume that the horizontal momentum flux due to IGWs averages out and that the corresponding vertical flux of vertical momentum is negligible, hence

Before identifying the conditions that lead to equivalence of (A.7) and the model's thermodynamic equation, let us consider the averaged dissipative heating

The first term on the rhs of (A.9) corresponds to the dissipation (12) explicitly accounted for in the model, that is, the simulated dissipation due to the shear of the planetary scale flow. The second term on the rhs of (A.9) results in additional dissipation owing to the vertical shear associated with gravity waves. Hence, an estimate of the total dissipation (A.9), rather than (12), should be compared with measurements of Lübken (1997). Taking advantage of Lindzen's theory, the last term on the rhs of (A.9) yields for a single saturated gravity wave

where *D*_{j} is determined by (21) and *u*^{″}_{j} represents the horizontal velocity perturbation in the direction of wave propagation. Consistent with our assumption that the gravity waves are propagating independently of each other, the dissipation owing to gravity wave shear is simply the sum of the individual contributions (A.10), hence

This expression may be regarded as the total turbulent dissipation associated with gravity wave breakdown and horizontal momentum diffusion. Note, however, that in the upper mesosphere, *ε*_{h} is neglible against each of the vertical dissipation terms. In section 4c (Fig. 7), the model estimate of (A.11) is discussed in the context of observational results.

We now proceed in evaluating (A.7) so as to accomplish agreement with (2), (10), and (11). From the first and second law of thermodynamics we have

with *s* denoting entropy per unit mass, which can be written as a function of potential temperature Θ = *T*(*p*_{00}/*p*)^{R/cp} in the case of an ideal gas. By means of (A.12), the wave sensible heat flux yields

For turbulent fluctuations, the second term on the rhs of (A.13) is usually negligible while the first is parameterized by diffusion. The situation is quite opposite in the case of gravity wave perturbations: Entropy displacements by wave–mean flow interaction are negligible and the pressure flux must be taken into account. Moreover, the latter is assumed to have a vertical component only, hence

As concerns the averaged turbulent heat flux **J**, we must apply Eq. (3) with averaged density *ρ* and averaged temperatures *T* and Θ.

Finally, we eliminate the adiabtic heating term occuring on the rhs of (A.7). For this purpose we multiply the momentum equation (A.1) by **v**^{″}_{3} and take the average. Owing to our basic assumption of quasi-linear wave propagation, we can discard terms of the order |**v**^{″}_{3}|^{3}. Then the wave kinetic energy equation yields

In line with quasi-linear propagation of IGWs we assume that and the advection of wave kinetic energy by the planetary-scale flow are negligible in (A.15), hence

Obviously, (A.16) is essentially negative in the region of gravity wave breakdown. The first term on the rhs is the residual work *ρ**W*_{res} discussed in section 2. The second term describes the rate of change of wave kinetic energy owing to vertical momentum diffusion, which is negative as well in the region of IGW saturation. Our final assumption is that

in the region of wave breaking. In other words, the dissipation owing to gravity wave shear and the loss of gravity wave kinetic energy owing to turbulent momentum diffusion average out on the scale of the planetary flow. Then, combination of (A.7), (A.9), (A.14), (A.16), and (A.17) leads to the thermodynamic equation (2) with the simulated dissipation and the energy deposition given by (10) and (11). In the case of the present idealized application of the Lindzen scheme, the rhs of (A.17) can be calculated analogously to (A.10), yielding *ρ**ϑN*^{2}/2 and thus proving that (A.17) is quantitatively valid.

## Footnotes

*Corresponding author address:* Erich Becker, Leibniz-Institut für Atmosphärenphysik an der Universität Rostock e. V., Schlossstraße 6, 18225 Kühlungsborn, Germany. Email: becker@iap-kborn.de