Nonlinear interaction of gravity waves in a nonisothermal and dissipative atmosphere

Starting from a set of fully nonlinear equations, this paper studies that two initial gravity wave packets interact to produce a third substantial packet in a nonisothermal and dissipative atmosphere. The effects of the inhomogeneous temperature and dissipation on interaction are revealed. Numerical experiments indicate that significant energy exchange occurs through the nonlinear interaction in a nonisothermal and dissipative atmosphere. Because of the variability of wavelengths and frequencies of interacting waves, the interaction in an inhomogeneous temperature field is characterised by the nonresonance. The nonresonant three waves mismatch mainly in the vertical wavelengths, but match in the horizontal wavelengths, and their frequencies also tend to match throughout the interaction. Below 80 km, the influence of atmospheric dissipation on the interaction is rather weak due to small diffusivities. With the further propagation of wave above 80 km, the exponentially increasing atmospheric dissipation leads to the remarkable decay and slowly upward propagation of wave energy. Even so, the dissipation below 110 km is not enough to decrease the vertical wavelength of wave. The dissipation seems neither to prevent the interaction occurrence nor to prolong the period of wave energy exchange, which is different from the theoretical prediction based on the linearised equations. The match relationship and wave energy evolution in numerical experiments are helpful in further investigating interaction of gravity waves in the middle atmosphere based on experimental observations.


Introduction
Atmospheric wave perturbations include tidal, planetary and gravity waves.These atmospheric waves propagate from the lower atmosphere to the mesosphere and lower thermosphere (MLT) and deposit their energy and momentum in the MLT due to their instability and dissipation, which has substantial dynamical effects on the large-scale circulation and thermal structure of the MLT (Lindzen, 1981;Sably, 1984).Due to the exponential decrease of the atmospheric density, atmospheric waves propagating into the MLT usually have large amplitudes, thus strong wave-wave and wave-flow interactions may frequently occur in the MLT (Fritts and Vincent, 1987;Liu and Hagan, 1998;Pancheva, 2000;Zhou et al., 2000;Jacobi et al., 2006;Chang et al., 2011;Yue et al., 2012;Huang et al., 2013a).Nonlinear interaction between planetary waves and tides is proposed as a mechanism responsible for the observed variability of tides because the beat between tides and generated secondary waves leads to modulation of tidal amplitudes at periods of planetary waves (Teitelbaum and Vial, 1991;Kamalabadi et al., 1997;Jacobi et al., 2001;Beard et al., 1999;Pancheva, 2000;Liu et al., 2007;Huang et al., 2013b).Especially, nonlinear interaction between tides and quasi-stationary planetary waves is regarded as a significant excitation source of nonmigrating tides (Hagan and Published by Copernicus Publications on behalf of the European Geosciences Union.Roble, 2001;Mayr et al., 2003;Liu et al., 2007;Chang et al., 2009).Interaction of tides with gravity waves can also result in the variability of tidal amplitude and phase in the MLT (Liu and Hagan, 1998;Mayr et al., 2003;Ortland and Alexander, 2006;A. Z. Liu et al., 2013).Nonlinear interaction of gravity waves redistributes energy and momentum among different spectral components, which may be a cause of the universal spectrum of gravity waves observed in the MLT (Tsuda et al., 1989;Yeh and Liu, 1985;Smith et al., 1987;Hines, 1991).Hence, wave-wave interactions are an important dynamics factor of the MLT region.
Gravity waves are small scale perturbations, which is different from planetary waves and tides with global scale.Dominant sources of gravity waves are convective storm, airflow over topography, wind shear and geostrophic adjustment in the lower atmosphere (Fritts and Alexander, 2003).In the MLT, body forcing accompanying wave dissipation and breaking, wave-wave interaction and auroral heating can locally generate upward and downward gravity waves (Snively and Pasko, 2003;Vadas et al., 2003;Li et al., 2007;Huang et al., 2012;Vadas and Liu, 2009).Satellite (Preusse et al., 2009;Alexander and Teitelbaum, 2011) and ground-based observations (Tsuda et al., 1990;Hickey et al., 1998;Vincent and Alexander, 2000;Dou et al., 2010;Yue et al., 2013;Zhang et al., 2013) show that small scale gravity waves have also broad temporal and spatial scales.When these gravity waves with different scales interact with each other, a sum or difference new wave is usually excited (Klostermeyer, 1984(Klostermeyer, , 1990;;Widdel et al., 1994;Rüster, 1997;Wüst and Bittner, 2006;Huang et al., 2009Huang et al., , 2012)), while numerical experiments with GCM presented by Chang et al. (2011) demonstrate that in interactions among large scale planetary waves and tides, both sum and difference new waves are generated, and evidence of sum and difference secondary waves is also observed in radar measurements (Jacobi et al., 1998;Beard et al., 1997Beard et al., , 1999;;Pancheva, 2006;Huang et al., 2013b), moreover, tidal modulation due to the beat is generally clear (Teitelbaum and Vial, 1991;Beard et al., 1997Beard et al., , 1999;;Huang et al., 2013b).Hence, relative to interactions among planetary waves and tides, nonlinear interaction of gravity waves shows some special characteristics.
In fact, the phenomenon of interaction for gravity waves in deep water was theoretically demonstrated by Phillips (1960).He pointed out that two waves with phases k 1 • r − ω 1 t and k 2 • r − ω 2 t could act as forcing terms in the linearised equations to force a new mode with their sum or difference phase.If the new mode is a free wave, resonant interaction takes place, and significant energy exchange among the wave triad occurs.The wave vector k and frequency ω matching conditions in resonance are expressed as where the subscripts j = 1, 2 and 3 denote three interacting waves, respectively.As free modes, the frequencies and wavenumbers of three interacting waves require to obey the dispersion relation ω j = (k j ), which can be written in a two-dimensional atmosphere as following expression (Hines, 1960) where k x and k z are the horizontal and vertical components of the wave vector k, respectively; N is the buoyancy frequency; v a is the acoustic speed; and ω a is the acousticcutoff frequency.As functions of atmospheric temperature, the buoyancy and acoustic-cutoff frequencies and acoustic speed vary in height.Hereby, whether or not resonant interaction can take place depends on the dispersion relation of gravity waves.For instance, resonant interaction cannot arise for deep water gravity waves because it is impossible to find three waves satisfying the resonant conditions (1) and ( 2) under the dispersion relation of ρ = (gk) 1 2 , where g is the acceleration due to gravity (Phillips, 1960).
Based on the linearised equations under the condition of weak nonlinear phenomena of waves in fluids, the properties of resonant and nonresonant interactions of oceanic and atmospheric internal gravity waves were extensively studied, and the effects of interaction on spectral energy transfer and universal spectrum formation were discussed (Phillips, 1960;Yeh and Liu, 1981;Müller et al., 1986;Inhester, 1987;Dong and Yeh, 1988;Fritts et al., 1992;Yi and Xiao, 1997;Stenflo and Shukla, 2009).Fritts et al. (1992) presented general interaction equations and coefficients for resonant and nonresonant interactions of gravity waves in a compressible atmosphere with rotation and small wind shear.The influence of viscous dissipation on nonlinear interaction was also examined (Yeh and Dong, 1989;Vanneste, 1995;Yi and Xiao, 1997).Model studies further revealed the characteristics of interaction of gravity waves in a nonlinear regime (Dunkerton, 1987;Vanneste, 1995;Yi, 1999;Zhang and Yi, 2004;Huang et al., 2009Huang et al., , 2011)).Fully nonlinear experiments indicate that energy mainly transfers from the high frequency primary wave to the excited wave in resonant and nonresonant interactions (Huang et al., 2009(Huang et al., , 2011)), which is completed in several periods of the primary waves (Vanneste, 1995).Because of spatially localised interaction for gravity wave packets, energy exchange is permanent rather than periodic from theoretical study in a linear regime (Huang et al., 2009(Huang et al., , 2011)).Huang et al. (2009) introduced a detuning degree of interaction which was proposed to measure whether there is an effective energy exchange in nonlinear interaction.Resonant excitation is reversible due to the satisfaction of the resonant conditions while nonresonant excitation is not so (Huang et al., 2012).Significant third-order resonance of gravity waves can arise in the atmosphere, which occurs through direct interaction of waves that satisfy the corresponding third-order resonant conditions instead of through an intermediate forced mode involved in interaction suggested by a previous study (Huang et al., 2013c).Because of extreme complexity, theoretical study of interaction based on the linearised equations used the assumptions that the background atmosphere was isothermal, and the molecular kinematic diffusivity was set to be constant if the effect of viscosity was examined.The theoretical study predicts that there is an amplitude threshold on interaction in the presence of molecular viscous dissipation, which means that the atmospheric dissipation may prevent interaction for small amplitude waves.And if interaction takes place, the dissipation makes energy exchange time larger, at this case, interacting waves do not obey the dispersion relation any more (Dong and Yeh, 1988;Yeh and Dong, 1989;Yi and Xiao, 1997).Hereby, the dissipative effect on interaction in the MLT is believed to be important since the molecular kinematic diffusivity exponentially increases with height (Yi and Xiao, 1997;Yi, 1999).Similarly, for the sake of examining the match relationship among interacting waves and exploring the dependence of energy exchange upon the detuning degree of interaction, previous numerical studies focused on resonant and nonresonant interactions in an isothermal atmosphere (Huang et al., 2009(Huang et al., , 2011(Huang et al., , 2013c)), and only the effect of constant molecular kinematic diffusivity was investigated (Yi, 1999;Huang et al., 2007).Even so, numerical results based on fully nonlinear model indicate that the restriction of amplitude threshold on interaction in the presence of viscous dissipation seems to be rather loose, and all interacting waves still satisfy the dispersion relation (Huang et al., 2007), which is different from the prediction of theoretical study from the linearised equations (Dong and Yeh, 1988;Yeh and Dong, 1989;Yi and Xiao, 1997).This indicates that interaction of gravity waves is in need of nonlinear investigation.
In the realistic atmosphere, the background temperature shows an increase or decrease with height in different atmospheric layers, respectively, which leads to a complex profile of buoyancy frequency.The molecular kinematic and thermal diffusivities exponentially increase with the decreasing atmospheric density, and strong eddy diffusion exists in the mesopause region (Balsley et al., 1983;Hocking, 1987;Lübken, 1997;Xu et al., 2009a).Intensely dissipative gradient can decrease the vertical wavelength of gravity waves propagating in the thermosphere (Zhang and Yi, 2002;Vadas and Fritts, 2005).However, we lack an adequate understanding at present of the influences of atmospheric inhomogeneous temperature and dissipation on nonlinear interaction of gravity waves.It can be inferred that interaction should mainly be nonresonant in a nonisothermal atmosphere because it is difficult for interacting waves to satisfy the resonant conditions for a long period.In addition, there are shear winds in the realistic atmosphere, which has important influences on interaction of gravity waves.Not only does the Doppler-shift due to wind shear change the wavelengths and frequencies of interacting waves, but also wave propagating along (against) the wind shear absorbs (transfers) energy from (to) the background atmosphere because of the Reynolds stress work (Huang et al., 2010), thus the effect of background wind on interaction is more complex.In this paper, we will numerically study interaction of gravity waves in a nonisothermal and dissipative atmosphere by neglecting the background wind, and explored the effect of inhomogeneous temperature and dissipation on interaction.The influence of wind fields will be focused on in the further work.
In Sect.2, the fully nonlinear numerical model is briefly described, and a nonisothermal atmosphere with molecular and eddy diffusion is introduced.In Sect.3, wave excitation in a nonisothermal atmosphere is exhibited through numerical experiments, and match relationship of interacting waves is investigated.The effect of dissipation that resulted from molecular and eddy diffusion is examined in Sect. 4. A summary is provided in Sect. 5.

Governing equations and numerical scheme
The numerical study utilises a set of primitive hydrodynamic equations describing conservations of mass, momentum and energy, and the equation of state for ideal gas, which can be written as Here, we use a two-dimensional numerical model in the x-z-plane, where x and z are the horizontal and vertical (positive upward) coordinates, respectively.The operator D Dt = ∂ ∂t + v • ∇ is the substantial derivative.t is the time; v is the velocity vector with the horizontal and vertical components of u and w; and ρ, T and p are the atmospheric density, temperature and pressure, respectively.R = 287 J kg −1 K −1 is the special gas constant for air; γ = c p /c v (c p = 1005 J kg −1 K −1 and c v = 718 J kg −1 K −1 are the specific heats at constant pressure and volume, respectively).µ is the molecular dynamical diffusivity, which can be expressed as µ = ρ 0 υ, where υ is the molecular kinematic diffusivity, and ρ 0 is the background density.λ m is the molecular thermal conductivity, and the molecular thermal diffusivity may be expressed as k m = λ m ρ 0 c v .η e and k e are the eddy momentum and thermal diffusivities, respectively.The molecular and eddy diffusion, as the dissipative terms, are written at the right sides of Eq. ( 5).
In order to precisely simulate the long-lasting propagation and interaction of gravity waves in the atmosphere, a numerical scheme should be of high accuracy and fine stability.Here, a composite difference scheme with secondorder temporal and spatial precision is applied to discretizing the Eq. ( 5), which was described in detail by Huang et al. (2007), and a corresponding three-dimensional model was extended to investigate the propagation characteristics of gravity waves (Huang et al., 2008).In the numerical scheme, a uniform Eulerian mesh with staggered grids is used to eliminate the checkerboard error.To avoid boundary reflection, the lateral boundaries are set to be periodic, and projected characteristic line boundaries are employed at the top and bottom boundaries (Hu and Wu, 1984;Zhang and Yi, 1999).
Considering an explicit scheme used in the projected characteristic line boundaries, the time step should be restricted by the Courant condition: where the acoustic speed v a = √ γ RT , and v = √ u 2 + w 2 .We take t = 0.5 t c in numerical experiments.

Initial background
In the two-dimensional model, we assume that the initial background atmosphere is homogeneous in the horizontal.In order to investigate the effects of the vertically inhomogeneous temperature and dissipation on wave-wave interaction, we avoid the wave-flow energy exchange and Doppler-shift due to wind shear by assuming that the background atmosphere is windless.Figure 1 shows the background temperature and density, which are prescribed using the MSISE-90 model (Hedin, 1991) for a location of latitude 20 • N and longitude 0 • , a date of 15 July, and a local time of 12:00.Based on the background temperature shown in Fig. 1, the buoyancy frequency is evaluated, and is also shown in Fig. 1.
In a dissipative atmosphere, the molecular kinematic and thermal diffusivities in units of m 2 s −1 may be expressed as (Banks and Kockarts, 1973) Then, the molecular dynamical diffusivity µ is in units of kg m −1 s −1 , and the molecular thermal conductivity λ m has the units in J m −1 K −1 s −1 .Here, the profiles of the molecular kinematic and thermal diffusivities are almost the same as those of previous studies (Hickey and Yu, 2005;Ortland and Alexander, 2006;Xu et al., 2009b;Yu et al., 2009;Liu et al., 2010).Observational analyses show that the eddy diffusion is strong in a relatively narrow height range in the mesopause region (Balsley et al., 1983;Hocking, 1987;Lübken, 1997).We simply chose a Gaussian profile to represent the eddy momentum diffusivity, which is written as follows (Xu et al., 2009b;X. Liu et al., 2013) The maximum eddy momentum diffusivity is taken to be η max = 100 m 2 s −1 at the height of 100 km, which is comparable to the values derived from radar inferences of turbulent energy dissipation rates (Hocking, 1987).The profile of total momentum diffusivity υ + η e is shown in Fig. 1.On the whole, the total momentum diffusivity exhibits an approximately exponential increase with height.The eddy thermal diffusivity is calculated based on the Prandtl number P l = η e k e , and the Prandtl number is chosen to be P l = 0.7.The profile of total thermal diffusivity k m + k e is similar to that of the total momentum diffusivity υ + η e , thus it is not presented in Fig. 1.

Initial wave packets
In the initial background atmosphere, we introduce two discrete gravity wave packets as the initial wave perturbations, for which the horizontal velocity disturbances have the following Gaussion forms where u cj (j = 1 and 2) is the maximum horizontal wind amplitude; x cj and z cj are the initial geometric centre positions of wave packets in the horizontal and vertical directions, respectively; δ xj and δ zj represent the horizontal and vertical half-widths of wave packets, which are chosen to be δ xj = λ xj and δ zj = λ zj , where λ xj and λ zj are the horizontal and vertical wavelengths at the centre positions of wave packets, respectively.We choose a downward propagating wave 1 as an initial wave with frequency of ω = 1.80 × 10 −3 Rads −1 and horizontal wavelength of λ x = 50 km.According to the dispersion relation of gravity waves, the vertical wavelength of wave 1 is calculated to be λ z = 4.50 km at its centre height of 80 km.An upward propagating wave 2 that has frequency of ω = 0.32 × 10 −3 Rads −1 and horizontal wavelength of λ x = 130 km is chosen as the other initial wave, and then its vertical wavelength is λ z = −2.43km at the height of 60 km, where the negative sign of the vertical wavelength represents downward phase propagation.The initial parameters of waves 1 and 2 are listed in Table 1.The waves propagating to the MLT generally have large amplitudes, thus the dynamical or convective instability frequently occurs (Lindzen, 1981;Fritts and Alexander, 2003).The studies showed that the body forcing due to local instability and breaking of waves could create high frequency gravity waves with large vertical scales (Fritts et al., 2002;Vadas et al., 2003).Besides this, the upgoing high frequency waves are likely to experience the reflection in shear wind and temperature gradient in the MLT, and then turn to propagate downwards (Huang et al., 2010).Hence, we select the downward propagating wave 1 as a high frequency perturbation.These two initial waves listed in Table 1 have typical temporal and spatial scales observed in the middle atmosphere (Tsuda et al., 1990;Alexander and Teitelbaum, 2007;Snively et al., 2010;Vincent et al., 2013;Zhang et al., 2013).
In the numerical computation, considering that the spatial scales of the interacting waves will vary with wave propagation, the grid sizes are set to be x = 3.5 km and z = 0.2 km in the horizontal and vertical directions, and the computational domains are chosen to be 0 km ≤ x ≤ 2520 km and 0 km ≤ z ≤160 km, respectively.

Effect of inhomogeneous temperature
Firstly, we investigate the effect of the inhomogeneous temperature on wave-wave interaction.Hence, in case 1, we neglect the atmospheric dissipation, including the molecular and eddy diffusion.The inhomogeneous temperature is shown in Fig. 1.

Wave excitation in interaction
Figure 2 exhibits the propagation and interaction of gravity waves in a nonisothermal atmosphere by using square root density-weighted horizontal velocity disturbances, which are calculated from the expression of u (z) = ρ 0 (z) ρ 0r 1 2 u (z), where ρ 0r is the background density at a reference level of z = 60 km.When waves 1 and 2 meet with each other at t = 3 h, the nonlinear interaction occurs.After 4 h of propagation, a new wave (wave 3) excited through the interaction clearly appears and staggers with waves 1 and 2. At t = 10 h, the upward propagating wave 3 departs from waves 1 and 2. In this process, the maximum value of the square root density-weighted horizontal velocity disturbances of wave 3 reaches 0.34 m s −1 at t = 10 h due to obtaining energy through the nonlinear interaction, while the velocity disturbance peak of wave 1 drops from 1.19 m s −1 at the beginning time to 0.47 m s −1 at t = 10 h as a result of wave energy transfer and coverage expansion.The maximum velocity disturbance of wave 2 keeps almost invariable in the propagation and interaction.
The wavenumber spectra of the interacting waves can be obtained by making a Discrete Fourier Transform (DFT) on the square root density-weighted horizontal velocity disturbances over the whole computational domain.The calculated wavenumber spectra are normalised by its maximum spectral magnitude.Figure 3 shows the normalised wavenumber spectra of the interacting waves at several selected times.It can be clearly seen from Fig. 3 that the spectral peaks of waves 1-3 maintain unchanged in the horizontal direction, however, obviously vary with time in the vertical direction.
According the corresponding wavenumbers at the spectral peaks, we calculate the dominant wavelengths of waves 1-3.The dominant horizontal wavelengths of waves 1-3 are 50.4km, 132.63 km and 36.52 km throughout the propagation and interaction process, respectively.The tiny differences between the numerically calculated and initially specified wavelengths of waves 1 and 2 are due to limited spectral resolution in the numerical computation.Figure 4 shows the evolution of the absolute vertical wavelengths and frequencies of the interacting waves with height.In Fig. 4, the wavelengths and frequencies denoted by the asterisks are marked at the height of their horizontal velocity disturbance peaks, and these parameters for wave 3 are presented from 3 h when this excited wave appears as a very weak perturbation.The wave frequency is calculated from the dispersion relation by using the background temperature and buoyancy frequency at the height of the disturbance peak.The vertically inhomogeneous background temperature leads to a remarkable variation of the vertical wavelengths and frequencies of waves 1-3 with height.In particular, wave 1 propagates through the largest height range due to its fastest vertical group speed, the variability of its wavelengths and frequencies is the most complex, while for wave 2 the opposite is the case.Because the buoyancy and acoustic cutoff frequencies and sound speed involved in the dispersion relation (3) have close relation to the background temperature, the response of wave frequency to the inhomogeneous temperature is distinguished from that of the vertical wavelength.These are different from the results in an isothermal atmosphere in which the wavelengths and frequencies of the interacting waves are constant in both the resonant and nonresonant interactions (Huang et al., 2009(Huang et al., , 2011)).In a nonisothermal atmosphere, it is usually difficult for the interacting waves to satisfy the resonant conditions for a long time due to the variability of their wavelengths and frequencies, thus the interaction should occur mainly through the nonresonance.In case 1, it is easy to demonstrate that waves 1-3 do not satisfy the resonant condition during the interaction period form 3 h to 6 h.When the new wave appears at t = 3 h, the sum and difference detuning degrees are calculated to be about 0.15 and 0.86, respectively.Hence, the sum nonresonant interaction arises due to a small sum detuning degree.The vertical wavelength (λ z3 = −4.21km) wave 3 is smaller than the matching value of −4.85 km derived from the wavelengths λ z1 = 4.71 km and λ z2 = −2.39km of waves 1 and 2, and its frequency ω 3 = 2.05 × 10 −3 Rad s −1 is slightly larger than the sum frequency 2.01×10 −3 Rad s −1 of waves 1 and 2 (ω 1 = 1.70 × 10 −3 Rad s −1 and ω 2 = Fig. 5. Evolution of wave energies in propagation and interaction for cases 1 (left) and 2 (right).The cyan, green and red lines denote the energies of waves 1-3, respectively; and the blue line denotes their total energy.0.31 × 10 −3 Rad s −1 ).At t = 6 h, the wavelength and frequency of λ z3 = −3.90km and ω 3 = 2.08 × 10 −3 Rad s −1 for wave 3 do not equal the matching values of −4.84 km and 1.98 × 10 −3 Rad s −1 , either.As we expected, wave 3 is excited through the sum nonresonant interaction.Even so, we can note an interesting phenomenon that the horizontal wavelength of wave 3 equals the matching wavelength of 36.52 km, which indicates that the interacting waves meet the matching condition in the horizontal scale.In an isothermal atmosphere, the nonresonant wave triad mismatches in the vertical scale but matches in in the horizontal scale, and their frequencies tend to match, moreover, such a relationship is explained based on the energy exchange in the nonresonance (Huang et al., 2012).Therefore, the match relationship of the interaction in a nonisothermal atmosphere is consistent with those in an isothermal atmosphere.Considering the atmospheric temperature is inhomogeneous mainly in the vertical direction, the approximate match in the horizontal wavelengths is an important characteristic of the interaction of gravity waves in the realistic atmosphere.

Wave energy exchange
Since the wavenumber spectra of waves 1-3 are apart from each other, as shown in Fig. 3, waves 1-3 can be extracted by applying a band-pass filter, respectively.Using the extracted perturbation quantities, we calculate the energies of waves 1-3 at each integral hour by integrating the wave energy den- over the whole computational domain, where p is the pressure perturbation, and w is the vertical velocity disturbance.Figure 5 exhibits the evolution of wave energies in the propagation and interaction.The wave energy experiences a rapid exchange between 3 h and 6 h, and mainly transfers from wave 1 to wave 3, which is in good agreement with previous numerical results in an isothermal atmosphere (Huang et al., 2009(Huang et al., , 2011)).At the beginning time, the initial energies of waves 1 and 2 are 1.36 × 10 5 J and 1.28 × 10 5 J, respectively.After 10 h, Fig. 6.Propagation paths of interacting waves.The cyan, green and red lines denote the paths of waves 1-3, respectively.The asterisk denotes the position of wave energy centre at each integer hour, and the position of wave 3 is presented from t = 3 h.the energy of wave 2 approximately equals its initial value.The energy of wave 1 decreases to 9.72×10 4 J, while wave 3 obtains the energy of 2.91×10 4 J through the nonlinear interaction, which is about 21.4 % of the initial energy of wave 1.Compared with the interactions in an isothermal atmosphere under the condition of the same initial amplitudes of waves 1 and 2 (Huang et al., 2009(Huang et al., , 2011)), this energy exchange is strong.The total energy of waves 1-3 shows a slight decrease with loss of about 3.67 % at t = 10 h due to the small wave energy conversion to the large-scale background movement in long-lasting nonlinear propagation and the energy loss in the numerical computation and filtering.Hence, the total energy of these interacting waves is approximately conservative during the propagation and interaction.
Figure 6 shows the energy propagation paths of the interacting waves.The asterisk is marked at the position of wave energy center.Wave 2 propagates a rather short vertical distance in which the background temperature has a slight decrease, thus its path is almost straight.The paths of waves 1 and 3 are obviously affected by the nonisothermal background and wave energy exchange.When wave 3 appears at t = 3 h, its energy center locates at x 3 = 658 km between the energy centers (x 1 = 651 km and x 2 = 672 km) of waves 1 and 2 in the horizontal direction, and z 3 = 65 km slightly higher than the energy centres (z 1 = 63.2 km and z 2 = 61 km) of waves 1 and 2 in the vertical direction.After 10 h of propagation, the energy centre of wave 3 reaches x 3 = 948.5 km and z 3 = 91 km, and is close to the mesopause, as shown in Fig. 2.

Effect of atmospheric dissipation
In order to investigate the influence of the atmospheric dissipation on the nonlinear interaction of gravity waves, we simulate case 2. In case 2, besides using the same background temperature and density (shown in Fig. 1) and the same initial wave parameters (listed in Table 1) as those in case 1, the molecular and eddy momentum and thermal diffusion is added to the initial background atmosphere, and the total momentum diffusivity is shown in Fig. 1.
Figure 7 exhibits the new wave excitation in a nonisothermal and dissipative atmosphere.Compared with case 1 shown in Fig. 2, the amplitudes of waves 1 and 2 in case 2 have a slight decrease owing to the small diffusivities below 80 km, as shown in Fig. 1, while the square root densityweighted horizontal velocity disturbance peak (0.28 m s −1 ) of wave 3 at t = 10 h is evidently weaker than 0.34 m s −1 in case 1.The spectral analysis (not presented) indicates that the evolutions of the dominant horizontal and vertical wavelengths of waves 1-3 with time are consistent with those in case 1.According to the ray trace solutions presented by Vadas and Fritts (2005), when gravity waves propagate upward in the thermosphere, their vertical wavelengths may diminish because of the effect of diffusion gradient.Previous numerical studies with similar diffusivity profiles indicate that the dissipation below 100 km is not strong enough to cause a decrease of the vertical wavelengths of upward propagating gravity waves, and the vertical wavelength begins to reduce at a certain height above 110 km (Zhang and Yi, 2002).In case 2, the excited wave does not attain the height of 100 km at t = 10 h, as shown in Fig. 7, thus its wavelength is mainly affected by the background temperature, and shows the same evolution as that in case 1.When wave 3 continues to propagate upward, it should experience a vertical wavelength decrease and a stronger decay.Nevertheless, this is the effect of the dissipation on the wave propagation rather than the wave excitation.
We calculate the energies of waves 1-3 at each integral hour in case 2, which is also shown in Fig. 5 for the convenience of comparison.It can be noted from Fig. 5 that the wave energy exchange time does not become prolonged in the presence of the dissipation, which differs from theoretical study based on the linearised equations (Yi and Xiao, 1997).At t = 3 h, the energy (6.0 × 10 3 J) of the new wave almost equals the values of 6.05 × 10 3 J in case 1.At t = 6 h, the energy of wave 3 attains the maximum of 2.79 × 10 4 J slightly less than the corresponding energy of 2.87 × 10 4 J in case 1, and then slightly declines to 2.71×10 4 J at t = 7 h.This indicates that the dissipation is rather weak due to wave propagating in a relative low height range.After this, as the approximately exponential increase of the diffusivities, the energy of wave 3 rapidly drops to 1.46 × 10 4 J at t = 10 h, which is much smaller than the value of 2.91 × 10 4 J in case 1. Relative to case 1, one can see that the energy of wave 1 is slightly dissipated in the first hour.Because propagating in a relative low region with small diffusivities, wave 2 with small vertical scale shows rather slight energy decay though the scaledependent dissipation is inversely proportional to the square of wavelength.At t = 10 h, the total energy of three waves decreases about 12.7 %, which is larger than the energy loss of 3.67 % in case 1.On the other hand, the exponentially increasing dissipation causes slowly upward propagation of the energy centre of wave 3 because the foreside of wave 3 suffers a more intense dissipation.At t = 10 h, the energy centre of wave 3 arrives at the height of 84.2 km lower than 91 km in case 1, as shown in Figs. 2 and 7. Therefore, on the whole, the influence of the atmospheric diffusion on wave-wave interaction is rather weak below 80 km due to the small diffusivities.With the further propagation above 80 km, the wave energy dissipation becomes stronger and stronger, and this inhomogeneous dissipation leads to a slowly upward propagation of wave energy.The mesospheric dissipation can neither prolong the interaction time nor change the wavelengths of these waves, which is in agreement with the results of a constant molecular kinematic diffusivity (Huang et al., 2007).
Finally, cases 3 and 4 are performed to examine the general characteristics of the interaction in a nonisothermal and dissipative atmosphere.The initial parameters of waves 1 and 2 in cases 3 and 4 are listed in Table 2.The initial background atmospheres in cases 3 and 4 are chosen to be as those in cases 1 and 2, with and without the molecular and eddy diffusion, respectively.In the numerical computation, the horizontal grid size is adjusted to be x = 5.0 km.The new wave excitations in cases 3 and 4 are presented in Fig. 8. Similar to Figs. 2 and 7, Fig. 8 shows the effective occurrence of the interactions of gravity waves in the nonisothermal and dissipative atmosphere.Because of the approximately exponential increase of atmospheric diffusivity, at t = 10 h, the new wave in the presence of dissipation in case 4 is weaker than in case 3. Figure 9 shows the normalised wavenumber spectra of interacting waves at t = 0 and 10 h.According the wavenumber spectra, we can verify that when the new wave is excited at t = 3 h, its dominant horizontal wavelength (λ x3 = 50 km) in both cases 3 and 4 meets the matching condition with the wavelengths (λ x1 = 60 km and λ x2 = 300 km) of waves 1 and 2. The frequency of the excited wave 3 (ω 3 = 1.81 × 10 −3 Rad s −1 ) is slightly larger than the value (1.75×10 −3 Rad s −1 ) of ω 1 +ω 2 (ω 1 = 1.59 × 10 −3 Rad s −1 and ω 2 = 0.16 × 10 −3 Rad s −1 ), and the dominant vertical wavelength (λ z3 = −5.0km) of wave 3 mismatches the matching value of 5.50 km derived from the wavelengths (λ z1 = 5.71 km and λ z2 = −2.80km) of waves 1 and 2. At t = 10 h, as shown in Fig. 9, wave 3 Fig.10.Evolution of wave energies in propagation and interaction for cases 3 (left) and 4 (right).The cyan, green and red lines denote the energies of waves 1-3, respectively.
in cases 3 and 4 still has the same dominant wavelengths, which indicates that the influence of dissipation on the vertical wavelength is not obvious below 110 km.At this time, the horizontal wavelengths of waves 1-3 remain to be 60, 300 and 50 km, and their vertical wavelengths are changed to be 4.21, −2.76 and −4.10 km due to the inhomogeneous background temperature.Figure 10 shows the variation of wave energies with time in cases 3 and 4. In case 3, the energy (2.50 × 10 4 J) of wave 3 at t = 10 h is about 12.25 % of the initial energy (2.04 × 10 5 J) of wave 1.In case 4, the energy of wave 3 reaches the maximum value of 2.30 × 10 4 J at t = 6 h, and then decreases to 1.05 × 10 4 J at t = 10 h.Wave 3 undergoes a stronger energy decay in case 4 than in case 2 because it propagates to a higher height in case 4. It is clearly seen from Fig. 10 that the energy exchange time in case 4 is in good agreement with that in case 3. The primary influence of atmospheric dissipation is to dissipate the energy of the interacting waves.Therefore, the general features of interaction in cases 3 and 4, such as the match relationship, the wave energy evolution and the influence of dissipation, are consistent with those in cases 1 and 2.

Summary
In this paper, the nonlinear interaction of gravity waves in a nonisothermal and dissipative atmosphere is numerically investigated, and the whole interaction process is clearly exhibited.The numerical experiments confirm that the significant nonlinear interaction can take place in a nonisothermal and dissipative atmosphere, and considerable energy is transferred mainly from the high frequency primary wave to the excited wave in the interaction, which is consistent with previous results in an isothermal atmosphere.In a nonisothermal atmosphere, the wavelengths and frequencies of the interacting waves vary with height, thus the interaction shows a nonresonant feature.The numerical results indicate that although their wavelengths and frequencies show temporal variability, the nonresonant three waves mismatch mainly in the vertical wavelengths but match in the horizontal wavelengths, and their frequencies also tend to match throughout the interaction.The match relationship is in agreement with the nonresonance in an isothermal atmosphere.Below 80 km, the effect of atmospheric dissipation on the interaction is weak due to the small diffusivities.As the approximately exponential increase of atmospheric diffusivities, the wave energy decay becomes sharp with the upward propagation of wave in the mesopause region, which also leads to slowly upward propagation of wave energy center.But then, the dissipation in the mesopause region is not strong enough to cause a decrease in the vertical scale of wave.The primary effect of dissipation is to dissipate the energy of the interacting waves.The dissipation seems neither to prevent the interaction occurrence nor to prolong the period of wave energy exchange, which is different from the theoretical prediction based on the linearised equations.
This numerical study indicates that the nonlinear interaction of gravity waves can frequently arise in the middle atmosphere.It is well-known that the interactions among global-scale tidal and planetary waves were extensively reported based on experimental observations.However, although several resonant gravity wave triads were presented based on radar, rocket, and satellite measurements (Klostermeyer, 1984(Klostermeyer, , 1990;;Rüster, 1997;Wüst and Bittner, 2006), observational study of the interaction events of gravity waves is much less.There are two possible reasons.One is that as shown in Fig. 2, a whole interaction process spans a large height range and lasts a long period, and the amplitudes of the interacting waves exhibit drastic variation in the process, thus the observational study of the localised interaction of gravity waves requires a high level of experimental measurement ability.The other is that the characteristics of the interactions of gravity waves in the realistic atmosphere are not very clear yet.As our numerical results revealed, despite that the interaction is characterised by nonresonance in the atmosphere, three interacting waves mismatch mainly in their vertical scales.Therefore, the match relationship and evolutions of wave packages and wave energies in our numerical experiments are helpful in investigating the interaction of gravity waves in the middle atmosphere based on the present and further experimental measurements.

Fig. 1 .
Fig.1.Vertical profiles of background temperature, buoyancy frequency, density and total momentum diffusivity.The solid lines with scale on the lower horizontal axis denote the background temperature and density in left and right panels, respectively.The dashed lines with scale on the upper horizontal axis denote the buoyancy frequency and total momentum diffusivity in left and right panels, respectively.

Fig. 2 .
Fig. 2.Wave-wave nonlinear interaction in a nonisothermal atmosphere by using square root density-weighted horizontal velocity disturbances.

Fig. 4 .
Fig. 4. Variations of (left) absolute vertical wavelengths and (right)wave frequencies with time.The cyan, green and red lines represent waves 1-3, respectively.The asterisks marked at the height of disturbance peak denote the values at each integer hour, and the parameters of wave 3 are presented from t = 3 h.

Fig. 7 .
Fig.7.Wave-wave nonlinear interaction in a nonisothermal and dissipative atmosphere by using square root density-weighted horizontal velocity disturbances.