Nonlinear interaction of gravity and acoustic waves

Abstract

Here we have investigated the possibility of an inertio-acoustic wave-mode to be unstable with regards to gravity mode perturbations through non-linear triad interactions in the context of a shallow non-hydrostatic model. We have considered highly truncated Galerkin expansions of the perturbations around a resting, hydrostatic and isothermal background state in terms of the eigensolutions of the linear problem. For a single interacting wave triplet, we have shown that an acoustic mode cannot amplify a pair of inertio-gravity perturbations due to the high mismatch among the eigenfrequencies of the three interacting wave-modes, which requires an unrealistically high amplitude of the acoustic mode in order for pump wave instability to occur. In contrast, it has been demonstrated by analysing the dynamics of two triads coupled by a single mode that a non-hydrostatic gravity wave-mode participating in a nearly resonant interaction with two acoustic modes can be unstable to small amplitude perturbations associated with a pair of two hydrostatically balanced inertio-gravity wave-modes. This linear instability yields significant inter-triad energy exchanges if the nonlinearity associated with the second triplet containing the two hydrostatically balanced inertio-gravity modes is restored. Therefore, this inter-triad energy exchanges lead the acoustic modes to yield significant energy modulations in hydrostatic inertio-gravity wave modes. Consequently, our theory suggests that acoustic waves might play an important role in the transient phase of the three-dimensional adjustment process of the atmosphere to both hydrostatic and geostrophic balances.

Keywords:

• wave–wave interaction
• triad interaction
• hydrostatic adjustment
• acoustic waves
• inertia–gravity waves

1. Introduction

With the advent of high-resolution atmospheric models, there has been a renewed interest in the study of the normal modes of the non-hydrostatic atmospheric dynamics. Kasahara and Qian (2000) and Qian and Kasahara (2003) studied the linear normal mode function theory of the shallow non-hydrostatic model (White et al., 2005) in the contexts of spherical and the beta-plane geometries, respectively, and the theory has been augmented with the account of the non-traditional Coriolis terms (Kasahara, 2003a, 2003b). The linear normal mode function theory has also been presented by Kasahara (2004) for the full deep non-hydrostatic case, as well as by Kasahara and Gary (2006) for the Boussinesq system in which acoustic modes are absent. However, as the governing equations of the atmospheric dynamics are non-linear, a more accurate account of the normal mode theory should also include the effect of the non-linearity on the wave dynamics.

Recently, Raupp et al. (2019) extended the work of Kasahara and Qian (2000) by analysing both linear and weakly non-linear energetics of inertia–gravity and inertia–acoustic modes. In the latter case, they studied the dynamics of a single resonant triad interaction involving an inertia–gravity wave and two inertia–acoustic modes and showed that in this kind of resonant interaction an inertia–gravity wave essentially acts as a catalyst mode for the energy exchanges between the two inertia–acoustic waves, in the sense that it enables the interaction to occur and controls both the interaction period and the impacts of the energy modulations on the perturbed dynamical field variables. Therefore, comparing this finding with the previous investigations on the non-linear atmospheric wave theory in the hydrostatic context (Duffy, 1974; Domaracki and Lossch, 1977; Loesch and Deininger, 1979; Ripa, 1983a, 1983b; Vanneste and Vial, 1994; Raupp et al., 2008), it is clear that the role of an inertia–gravity mode in a resonant interaction involving inertia–acoustic modes is similar to the role of a Rossby mode in a resonant interaction involving two inertio-gravity waves.

In this article, we extend the work of Raupp et al. (2019) by further investigating the non-linear dynamics of the shallow non-hydrostatic equations. In particular, we are interested here in analysing the possibilities of an acoustic mode to be unstable to gravity wave perturbations, including the study of off-resonant wave triplets as well as the dynamics of two connected wave triads. In the latter analysis, we study whether a resonant triad involving inertia–acoustic and non-hydrostatic inertia–gravity waves can be unstable with regard to interacting triads (not necessarily resonant) of inertia–gravity modes. A motivation for this analysis stems from recent findings in the non-linear wave literature (Janssen, 2003; Smith and Lee, 2005; Bustamante et al., 2014) pointing out that although in a single interacting wave triad the resonance relation among the mode eigenfrequencies is crucial for significant energy exchanges to occur in the limit of weak non-linearity if one relaxes such assumption of weak non-linearity to take into account off-resonant wave triads, the mismatch among the wave frequencies within interacting triads might be important for the energy flow throughout the whole modal space. Bustamante et al. (2014) showed in a reduced dynamical system of two triads coupled by two modes (four-wave system) that, for moderate values of the modal amplitudes, the energy leakage of a triad increases as the mismatch among the wave eigenfrequencies of one triad approaches the frequency of amplitude (energy) modulations of the other wave triplet. This synchronization between the non-linear frequency and the linear mismatch frequency between different interacting triads has been labelled by the authors as precession resonance. Bustamante et al. (2014) also demonstrated the important role of precession resonance mechanism for increasing the efficiency of the energy flow throughout the whole system of several connected triads. Another mechanism that has been shown to yield significant energy transfers throughout the whole modal space in a diversity of wave problems is the modulational instability (Connaughton et al., 2010).

Another motivation for this study refers to the fact that, although acoustic modes are eigensolutions of compressible non-hydrostatic models of the atmospheric flow, due to the highly restrictive computational constraints related to their numerical treatment with explicit schemes (Pielke, 2002; Thuburn, 2011), in the numerical weather prediction models these acoustic waves are treated as noise in the sense that they are either filtered out (Davies et al., 2003; Klein, 2009) or subjected to strong damping associated with semi-implicit numerical schemes (Giraldo et al., 2010; Klemp et al., 2018). For example, Daley (1988) proposed a filter for acoustic modes based on normal mode expansion, in the same spirit of the method proposed by Tribbia (1979) in the hydrostatic context to filter out inertio-gravity waves.

However, our analysis of the highly truncated spectral model of the shallow non-hydrostatic equations (five-wave system) has demonstrated that an inertia–gravity mode participating in a resonant triad interaction with two inertia–acoustic modes can be unstable to small amplitude perturbations corresponding to a pair of lower frequency inertio-gravity modes. Since the higher the time frequency of an inertio-gravity wave the more pronounced the non-hydrostatic effect of vertical acceleration on it, our results suggest that ultra-high frequency acoustic modes can potentially yield amplitude (energy) modulations in hydrostatically balanced inertio-gravity waves through inter-triad energy exchanges. Therefore, this theoretical description suggests that acoustic modes excited by localized and explosive heating associated with convective storms might play an important role in both hydrostatic and geostrophic adjustments, as it will be discussed in Section 5.

The remainder of this article is organised as follows. In Section 2, we present the model equations, the pseudo-energy conservation and the linear eigenmodes. Section 3 presents the general solution of the non-linear problem based on its expansion in terms of the linear eigenmodes. In Section 3, we also show some energy constraints for the coupling coefficients of any interacting triad as a consequence of the pseudo-energy conservation. Section 4 analyses the reduced dynamics of one-wave triad and two triads coupled by a single mode to investigate the possibility of acoustic modes to excite hydrostatic inertio-gravity waves. The main conclusions are discussed in Section 5.

2. The model

2.1. Governing equations

In this article, we adopt the shallow non-hydrostatic model on a mid-latitude f-plane as it is the simplest context bearing the existence of both gravity and acoustic waves. The shallow non-hydrostatic model consists of a relaxation of the compressible primitive equations by adding the vertical acceleration term in the vertical momentum equation (White et al., 2005). We consider small-amplitude perturbations embedded in a resting, hydrostatic and isothermal background state. In this setting, the governing equations for the perturbations in Cartesian coordinates are: (1) $$\frac{\partial u\prime }{\partial t}-fv\prime +\frac{1}{{\rho }_0}\frac{\partial p\prime }{\partial x}=-\mathbf{V}\prime ·\nabla u\prime -w\prime \frac{\partial u\prime }{\partial z}+\frac{\rho \prime }{{\rho }_0^2}\frac{\partial p\prime }{\partial x},$$(1) (2) $$\frac{\partial v\prime }{\partial t}+fu\prime +\frac{1}{{\rho }_0}\frac{\partial p\prime }{\partial y}=-\mathbf{V}\prime ·\nabla v\prime -w\prime \frac{\partial v\prime }{\partial z}+\frac{\rho \prime }{{\rho }_0^2}\frac{\partial p\prime }{\partial y},$$(2) (3) $$\frac{\partial w\prime }{\partial t}+\frac{1}{{\rho }_0}\left(\frac{\partial p\prime }{\partial z}+\frac{g}{C_s^2}p\prime -\theta \prime \right)=-\mathbf{V}\prime ·\nabla w\prime -w\prime \frac{\partial w\prime }{\partial z}+\frac{\rho \prime }{{\rho }_0^2}\frac{\partial \rho \prime }{\partial z}+g{\left(\frac{\rho \prime }{{\rho }_0}\right)}^2,$$(3) (4) $$\begin{array}{c}\frac{1}{C_s^2}\frac{\partial p\prime }{\partial t}-\frac{{\rho }_{0}g}{C_s^2}w\prime +{\rho }_0\left(\nabla ·\mathbf{V}\prime +\frac{\partial w\prime }{\partial z}\right)\\ =-\frac{1}{C_s^2}\left(\mathbf{V}\prime ·\nabla p\prime +w\prime \frac{\partial p\prime }{\partial z}\right)\\ -\rho \prime \left(\nabla ·\mathbf{V}\prime +\frac{\partial w\prime }{\partial z}\right)+\frac{1}{C_s^2}\frac{T\prime }{T_0}\left(\frac{\partial p\prime }{\partial t}-{\rho }_{0}gw\prime \right),\end{array}$$(4) (5) $$\frac{\partial \theta \prime }{\partial t}+{\rho }_0{N}^{2}w\prime =-\left(\mathbf{V}\prime ·\nabla \theta \prime +w\prime \frac{\partial \theta \prime }{\partial z}\right)+\frac{g}{C_s^2}\frac{T\prime }{T_0}\left(\frac{\partial p\prime }{\partial t}-{\rho }_{0}gw\prime \right),$$(5) (6) $$\frac{p\prime }{p_0}=\frac{\rho \prime }{{\rho }_0}+\frac{T\prime }{T_0}.$$(6) where $\theta \prime =\frac{g}{C_s^2}p\prime -g\rho \prime$ is related to entropy perturbation, and all the remaining variables, symbols and operators have their usual meanings and are defined in Table 1. The superscript prime denotes the perturbation fields, whereas the subscript $0$ denotes the background state quantities, defined by $u_0=v_0=w_0=0,T_0=\text{const}$ and (7) $${\rho }_0(z)={\rho }_0(0)e^{-\frac{z}{H}}$$(7) (8) $$p_0(z)=p_0(0)e^{-\frac{z}{H}}$$(8) with $H=\frac{R{T}_0}{g}$ representing the scale height of the isothermal atmosphere. Only the perturbation terms required to describe the non-linear triad interactions among the wave modes have been retained in the equations above, namely, the leading-order (quadratic) non-linear terms in the equations (1)–(5) and the linear terms of the equation of state.

Table 1. Definition of variables, symbols and operators.

x	Eastward distance
y	Northward distance
z	Altitude
t	Time
u	Zonal wind
v	Meridional wind
w	Vertical wind
$\mathit{V}={(u,v)}^T$	Horizontal wind vector
p	Pressure
ρ	Density
T	Temperature
R	Gas constant for dry air
cp	Specific heat at constant pressure
cv	Specific heat at constant volume
$\gamma =c_p/c_v$
$\kappa =R/c_p$
$T_0=243.878$ K	Temperature of the isothermal basic state
$L_x=4\times {10}^7\hspace{0.17em}\text{cos}({\varphi }_0)$ m	Zonal period
$L_y={10}^7$ m	Meridional period
zT = 18,000 m	Model’s top
$\Omega =2\pi /(24\times 60\times 60)$ s^−1	Earth’s rotation rate
${\varphi }_0=\pi /4$	Central latitude of $f-$ plane
$f=2\Omega \hspace{0.17em}\text{sin}({\varphi }_0)$	Coriolis parameter
g = 9.8 m s^−2	Earth’s gravitational acceleration
$C_S=\sqrt{\gamma R{T}_0}$	Linear sound wave speed
$\mathit{V}·\nabla h=u{\partial }_{x}h+v{\partial }_{y}h$
$\nabla ·\mathit{V}={\partial }_{x}u+{\partial }_{y}v$	Horizontal divergence

2.2. Pseudo-energy conservation

An useful tool to describe the dynamics of a Hamiltonian system in the spectral space refers to pseudo-energy conservation, since it is a conserved quantity that is quadratic to leading-order in terms of the perturbation variables (Ripa, 1981; Shepherd, 1990). Pseudo-energy is the energy related to the departure from a reference steady state and, for the isothermal background state considered here, the pseudo-energy of a compressible non-hydrostatic model is (Andrews, 1981; Shepherd, 1993): (9a) $$\mathcal{E}=\iiint \{\frac{1}{2}\rho (u\prime {}^2+v\prime {}^2+w\prime {}^2)+\rho c_p{T}_0{h}_1(\theta \prime /g{\rho }_0)+\rho c_p{T}_0{h}_2(p\prime /p_0)\}\mathit{dxdydz},$$(9a) (9b) $$h_1(x)=x-\hspace{0.17em}\text{log}\hspace{0.17em}\left(1+x\right),$$(9b) (9c) $$h_2(x)={(1+x)}^{\kappa }((1-\kappa )x^{\kappa }+\kappa x^{-(1-\kappa )}-1).$$(9c) where $\kappa =R/C_p.$ Following the studies of Ripa (1983a) and Vanneste and Vial (1994), to describe the non-linear wave interactions it is suitable to explicit the pseudo-energy conservation in terms of its quadratic and higher-order dependencies in terms of the field variables. In this way, by Taylor expanding the functions h₁ and h₂ and using $\rho ={\rho }_0+\rho \prime ,$ we can express the pseudo-energy invariance as: (10) $${\mathcal{E}}^{(2)}+{\mathcal{E}}^{(3)}+\mathcal{O}{(u\prime ,v\prime ,w\prime ,p\prime ,\theta \prime )}^4=\text{const}$$(10) where (11a) $${\mathcal{E}}^{(2)}=\iiint \left\{\frac{1}{2}{\rho }_0(u\prime {}^2+v\prime {}^2+w\prime {}^2)+\frac{\theta \prime {}^2}{2{\rho }_0{N}^2}+\frac{p\prime {}^2}{2{\rho }_0{C}_s^2}\right\}d\mathit{xdydz},$$(11a) (11b) $${\mathcal{E}}^{(3)}={\iiint }_V\left\{\frac{1}{2}\rho \prime (u\prime {}^2+v\prime {}^2+w\prime {}^2)+\frac{p\prime {}^3}{2{\rho }_0{p}_0{C}_s^2}+\frac{\rho \prime \theta \prime {}^2}{2{\rho }_0^2{N}^2}-\frac{\theta {\prime }^3}{3N^2{\rho }_0^{2}g}-\frac{4+\kappa }{6{\rho }_0{C}_s^2{p}_0}p\prime {}^3\right\}d\mathit{xdydz}$$(11b)

2.3. Linear eigenmodes

The linearized version of equations (1)–(5) for the mid-latitude f-plane approximation ($f=\text{const}$) has linear wave-mode solutions of the form (Qian and Kasahara, 2003): (12a) $$u_a(x,y,z,t)={\rho }_0^{-1/2}\left(\frac{-{\omega }_{a}f+k_{a}g{H}_{e}i{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\xi }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12a) (12b) $$v_a(x,y,z,t)={\rho }_0^{-1/2}i{\xi }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12b) (12c) $$w_a(x,y,z,t)={\rho }_0^{-1/2}ig{H}_e\left(\frac{i{\omega }_a{l}_a-k_a{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\eta }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12c) (12d) $$p_a(x,y,z,t)={\rho }_0^{1/2}g{H}_e\left(\frac{i{\omega }_a{l}_a-k_a{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\xi }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12d) (12e) $${\theta }_a(x,y,z,t)={\rho }_0^{1/2}g{H}_e\left(\frac{i{\omega }_a{l}_a-k_a{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\Theta }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12e) where the vertical structure functions ${\xi }_a(z),{\eta }_a(z)$ and ${\Theta }_a(z)$ are (13) $$\eta (z)=A\hspace{0.17em}\text{sin}(\lambda z),$$(13) (14) $$\xi (z)={\left(\frac{1}{C_s^2}-\frac{1}{g{H}_e}\right)}^{-1}A\left[\lambda \hspace{0.17em}\text{cos}(\lambda z)-\Gamma \hspace{0.17em}\text{sin}(\lambda z)\right],$$(14) (15) $$\Theta (z)=\frac{N^2}{\omega }A\hspace{0.17em}\text{sin}(\lambda z),$$(15)

In the equations above, A is an arbitrary constant, $\Gamma =\frac{g}{C_S^2}+\frac{1}{2{\rho }_0}\frac{d{\rho }_0}{dz}=\frac{1-2\kappa }{2H}$ is the parameter of adiabatic expansion (Eckart, 1960), with the second equality being valid for the isothermal background state considered here; H_e is the separation constant, also known as equivalent height (Taylor, 1936), and is related to the vertical eigenvalue λ through the relation (16) $${\lambda }^2=\left(\frac{1}{g{H}_e}-\frac{1}{C_s^2}\right)(N^2-{\omega }^2)-{\Gamma }^2,$$(16) with the eigenfrequencies ω satisfying the following dispersion relation (17) $${\omega }^5-C_S^2\left\{k^2+l^2+{\lambda }^2+\frac{1}{4H^2}+f^2{C}_S^{-2}\right\}{\omega }^3+C_S^2\left\{f^2({\lambda }^2+\frac{1}{4H^2})+N^2(k^2+l^2)\right\}\omega =0,$$(17) whose solutions represent a vortical mode ω = 0 and two pairs of eastward and westward propagating inertio-gravity and inertio-acoustic modes.

Considering the commonly adopted rigid-lid boundary conditions (18) $$\eta =0\text{at}z=0\text{and}\text{at}z=z_T,$$(18) with $z_T\approx 18$ Km yields the quantization of the vertical eigenvalue spectrum according to (19) $${\lambda }_m=\frac{m\pi }{z_T},m=1,2,3,\dots$$(19)

Similarly, for periodic solutions in the (x, y) directions, the horizontal wavenumbers are quantized according to the relations (20a) $$k_j=\frac{2\pi j}{L_x},j=0,\pm 1,\pm 2,\dots ,$$(20a) (20b) $$l_n=\frac{2\pi n}{L_y},n=0,\pm 1,\pm 2,\dots ,$$(20b) with L_x and L_y representing, respectively, the length of the zonal circle along the latitude 45° and the distance from the poles to the equator, that is, $L_x=(4\times {10}^7\times \hspace{0.17em}\text{cos}\frac{\pi }{4}$) m and $L_y\approx {10}^7$ m. Therefore, each particular linear eigenmode of system (1)–(5), which is labeled by the subscript a in equation (12(12a) $$u_a(x,y,z,t)={\rho }_0^{-1/2}\left(\frac{-{\omega }_{a}f+k_{a}g{H}_{e}i{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\xi }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12a) ), must be distinguished by its zonal, meridional and vertical quantum indexes j, n and m, respectively, along with its oscillation type that can be labelled by an index r.¹ A special solution is characterised by $\eta \equiv 0$ ($\Theta \equiv 0$). This mode is labelled as external mode and is characterised by (21a) $$\xi (z)=A{e}^{-\Gamma z}$$(21a) (21b) $$H_e=\frac{C_S^2}{g}$$(21b)

For this mode, the dispersion relation (17(17) $${\omega }^5-C_S^2\left\{k^2+l^2+{\lambda }^2+\frac{1}{4H^2}+f^2{C}_S^{-2}\right\}{\omega }^3+C_S^2\left\{f^2({\lambda }^2+\frac{1}{4H^2})+N^2(k^2+l^2)\right\}\omega =0,$$(17) ) gives only a vortical mode and a pair of eastward and westward propagating inertio-gravity modes. An important feature of the non-hydrostatic wave dynamics is that the equivalent height H_e is no longer constant for all the eigenmodes having the same vertical index m, as it is the case for the hydrostatic primitive equations in which the equivalent height depends only on the vertical wavenumber. Rather, in this model the equivalent height differs from each eigenmode $a=(j,n,m,r)$ as it relies on the eigenfrequency ω, apart from its dependence on the vertical wavenumber λ. Indeed, equation (16)(16) $${\lambda }^2=\left(\frac{1}{g{H}_e}-\frac{1}{C_s^2}\right)(N^2-{\omega }^2)-{\Gamma }^2,$$(16) shows that for the acoustic modes, whose eigenfrequencies are such that ${\omega }^2>N^2,$ it follows that $H_e>H_{\mathit{ext}},$ with $H_{\mathit{ext}}=\frac{C_S^2}{g}$ indicating the equivalent height of the external mode. In contrast, the gravity wave oscillation regime (${\omega }^2<N^2$) is characterized by $H_e<C_s{}^2/g.$

Figure 1 shows the dispersion curves of inertio-acoustic and inertio-gravity waves associated with the meridional wavenumber n = 1, for the first three baroclinic modes m = 1, 2, 3. Figure 1 shows the eastward branch of both wave types. To highlight the non-hydrostatic effect on the inertio-gravity waves, Fig. 2 displays only the inertio-gravity wave dispersion curves associated with the meridional wavenumber n = 1 and the first three baroclinic modes, together with the corresponding dispersion curves obtained with hydrostatic approximation. The eigenfrequency of hydrostatic inertio-gravity waves is obtained by determining the equivalent height H_e from the simplified version of equation (16)(16) $${\lambda }^2=\left(\frac{1}{g{H}_e}-\frac{1}{C_s^2}\right)(N^2-{\omega }^2)-{\Gamma }^2,$$(16) for ${\omega }^2\ll N^2,$ and then computing the eigenfrequency by the well-known shallow-water equation dispersion relation ${\omega }_{h_{k,l}}^{(\pm 1)}=\pm \sqrt{f^2+g{H}_e(k^2+l^2)}.$ Figure 2 shows that the non-hydrostatic effect on the inertio-gravity waves only becomes noticeable for zonal wavenumbers j > 400 ($k>9\times {10}^{-5}{m}^{-1}$). An important point to be observed in Fig. 1 is that, for a discrete spectrum of vertical eigenmodes resulting from a finite top z_T in equation (18)(18) $$\eta =0\text{at}z=0\text{and}\text{at}z=z_T,$$(18) , there is a large time-scale separation between inertio-acoustic and inertio-gravity modes. This high time–frequency separation increases for higher vertical wavenumbers and has important consequences for the nature of the non-linear interactions between these wave types. For example, this time–frequency separation prevents an acoustic mode to be unstable to a pair of gravity wave modes in a single triad interaction, as it will be shown in Section 4.

Fig. 1. Dispersion curves of inertia-acoustic (AI) and inertia-gravity (IG) waves corresponding to the first three baroclinic modes m = 1, 2, 3. All the curves are referred to the meridional index n = 1.

Fig. 2. Similar to Fig. 1, but only for the inertia-gravity waves (IG) and their corresponding dispersion curves obtained by hydrostatic approximation (H).

The exact conservation of ${\mathcal{E}}^{(2)}$ given by equation (11a)(11a) $${\mathcal{E}}^{(2)}=\iiint \left\{\frac{1}{2}{\rho }_0(u\prime {}^2+v\prime {}^2+w\prime {}^2)+\frac{\theta \prime {}^2}{2{\rho }_0{N}^2}+\frac{p\prime {}^2}{2{\rho }_0{C}_s^2}\right\}d\mathit{xdydz},$$(11a) in the linear case implies that the linear eigenmodes satisfy the following orthogonality relation (Kasahara and Qian, 2000; Qian and Kasahara, 2003): (22) $$i({\omega }_a-{\omega }_b){〈{\mathbf{u}}_{\mathbf{a}},{\mathbf{u}}_{\mathbf{b}}〉}_{{\mathcal{E}}^{(2)}}=0,$$(22) where ${\mathbf{u}}_{\mathbf{a}}$ and ${\mathbf{u}}_{\mathbf{b}}$ represent two arbitrary eigenvectors whose components are defined by equation (12(12a) $$u_a(x,y,z,t)={\rho }_0^{-1/2}\left(\frac{-{\omega }_{a}f+k_{a}g{H}_{e}i{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\xi }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12a) ) and ${〈,〉}_{{\mathcal{E}}^{(2)}}$ refers to the inner product in terms of pseudo-energy ${\mathcal{E}}^{(2)},$ given by $${〈{\mathbf{u}}_{\mathbf{a}},{\mathbf{u}}_{\mathbf{b}}〉}_{{\mathcal{E}}^{(2)}}={\int }_0^{z_T}{\int }_0^{L_x}{\int }_0^{L_y}\left[{\rho }_0\left(u_a^{*}{u}_b+v_a^{*}{v}_b+w_a^{*}{w}_b\right)+\frac{p_a^{*}{p}_b}{{\rho }_0{C}_s^2}+\frac{{\theta }_a^{*}{\theta }_b}{{\rho }_0{N}^2}\right]d\mathit{ydxdz}.$$ with the superscript ‘*’ indicating the complex conjugate.

3. General solution

3.1. Modal expansion

Now we use the orthogonality and completeness of the linear eigenmode functions described in the previous section to expand the solution of our nonlinear system (1)–(5) in a series (23) $$\mathbf{u}=\left[\begin{array}{c}u\prime (x,y,z,t)\\ v\prime (x,y,z,t)\\ w\prime (x,y,z,t)\\ p\prime (x,y,z,t)\\ \theta \prime (x,y,z,t)\end{array}\right]=\sum _a{A}_a(t){\mathbf{u}}_{\mathbf{a}}+C.C.,$$(23) where ‘C.C.’ indicates the complex conjugate of what is preceding, $A_a(t)$ refers to the complex-valued spectral amplitudes and the vector ${\mathbf{u}}_{\mathbf{a}}={[u_a,v_a,w_a,p_a,{\theta }_a]}^T$ represents the eigenvector function of a particular mode given by equation (12(12a) $$u_a(x,y,z,t)={\rho }_0^{-1/2}\left(\frac{-{\omega }_{a}f+k_{a}g{H}_{e}i{l}_a}{{\omega }_a^2-k_a^{2}g{H}_e}\right){\xi }_a(z)e^{i{k}_{a}x+i{l}_{a}y-i{\omega }_{a}t}$$(12a) ). Expansion series above is an exact solution of system (1)–(5) provided the mode amplitudes A_a satisfy (24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24)

In equation (24)(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) , $$E_a={〈{\mathbf{u}}_{\mathbf{a}},{\mathbf{u}}_{\mathbf{a}}〉}_{{\mathcal{E}}^{(2)}}$$ is the intrinsic energy (i.e. the squared norm) of the a-th mode and ${\delta }_{\mathit{abc}}={\omega }_a-{\omega }_b-{\omega }_c$ represents the mismatch among the linear eigenfrequencies of each mode triplet abc. If ${\delta }_{\mathit{abc}}=0,$ the triad is said to be resonant; the coupling constants ${\sigma }_a^{bc}$ are the projection of the non-linear terms due to the action of two modes b and c onto another mode a, viz., (25) $${\sigma }_a^{bc}={〈{\mathbf{u}}_a,\mathbb{B}({\mathbf{u}}_{\mathbf{b}},{\mathbf{u}}_{\mathbf{c}})+\mathbb{B}({\mathbf{u}}_{\mathbf{c}},{\mathbf{u}}_{\mathbf{b}})〉}_{{\mathcal{E}}^{(2)}}$$(25) where $\mathbb{B}$ is the bilinear operator containing the non-linear terms of system (1)–(5), and hence all the information on the non-linearity of our model equations is contained in these coefficients. The interacting triads are those whose coupling constants are non-zero. The orthogonality of the (x, y) basis functions $e^{\mathit{ikx}+\mathit{ily}}$ requires the wave modes of an interacting triad to satisfy (26a) $$k_a=k_b+k_c$$(26a) (26b) $$l_a=l_b+l_c$$(26b)

In contrast, the vertical coupling integrals involved in equation (25)(25) $${\sigma }_a^{bc}={〈{\mathbf{u}}_a,\mathbb{B}({\mathbf{u}}_{\mathbf{b}},{\mathbf{u}}_{\mathbf{c}})+\mathbb{B}({\mathbf{u}}_{\mathbf{c}},{\mathbf{u}}_{\mathbf{b}})〉}_{{\mathcal{E}}^{(2)}}$$(25) appear in the form $${\int }_0^{z_T}{\rho }_0^{-\frac{1}{2}}\hspace{0.17em}\text{cos}[({\lambda }_a\pm {\lambda }_b\pm {\lambda }_c)z]dz.$$

Thus, due to the presence of the ‘weight function’ ${\rho }_0^{-\frac{1}{2}}$ in the vertical coupling constants, unlike the horizontal wavenumbers, there is no an excluding selection rule imposed by the vertical structures of the triad components. However, as ρ₀ is a monotonically decreasing function of z, if ${\lambda }_a\pm {\lambda }_b\pm {\lambda }_c\ne 0$ the vertical coupling integrals will be small, so that the triads whose wave modes do satisfy the condition (27) $${\lambda }_a\pm {\lambda }_b\pm {\lambda }_c=0,$$(27) together with conditions (26)(26a) $$k_a=k_b+k_c$$(26a) are believed to undergo the strongest interactions, although condition (27(27) $${\lambda }_a\pm {\lambda }_b\pm {\lambda }_c=0,$$(27) ) is no longer excluding. However, this non-excluding nature of the vertical coupling integrals allows a triad of modes whose vertical wavenumbers do not satisfy condition (27(27) $${\lambda }_a\pm {\lambda }_b\pm {\lambda }_c=0,$$(27) ) to be excited by a primary triplet that does so, as it will be shown in the next section. Therefore, hereafter we will assume only the condition (26(26a) $$k_a=k_b+k_c$$(26a) ) to be met when referring to interacting triads.

3.2. Energy constraints for the interacting triads

As explained by Ripa (1981) in the context of barotropic Rossby waves and internal gravity waves in a vertical plane and Ripa (1983a) and Vanneste and Vial (1994) for the equatorial beta-plane and spherical geometry shallow-water equations, respectively, the conserved quantities which are quadratic to lowest order in terms of wave disturbances lead to relations among the coupling constants of an interacting triad. Thus, in what follows we shall apply their approach in our non-hydrostatic context. In this way, substituting the mode expansion (23(23) $$\mathbf{u}=\left[\begin{array}{c}u\prime (x,y,z,t)\\ v\prime (x,y,z,t)\\ w\prime (x,y,z,t)\\ p\prime (x,y,z,t)\\ \theta \prime (x,y,z,t)\end{array}\right]=\sum _a{A}_a(t){\mathbf{u}}_{\mathbf{a}}+C.C.,$$(23) ) into equations (11a)(11a) $${\mathcal{E}}^{(2)}=\iiint \left\{\frac{1}{2}{\rho }_0(u\prime {}^2+v\prime {}^2+w\prime {}^2)+\frac{\theta \prime {}^2}{2{\rho }_0{N}^2}+\frac{p\prime {}^2}{2{\rho }_0{C}_s^2}\right\}d\mathit{xdydz},$$(11a) and (11 b) yields: (28a) $${\mathcal{E}}^{(2)}=\sum _a{E}_a|A_a{|}^2$$(28a) (28b) $${\mathcal{E}}^{(3)}=\sum _a\sum _b\sum _{c}2S_{\mathit{abc}}\mathrm{Re}(A_a{A}_b^{*}{A}_c^{*}{e}^{i{\delta }_{\mathit{abc}}t})$$(28b)

Equations above show that, as a consequence of the orthogonality relation (22(22) $$i({\omega }_a-{\omega }_b){〈{\mathbf{u}}_{\mathbf{a}},{\mathbf{u}}_{\mathbf{b}}〉}_{{\mathcal{E}}^{(2)}}=0,$$(22) ), the leading-order (quadratic) pseudo-energy has a diagonalised representation in terms of the linear eigenmodes, whereas the cubic energy ${\mathcal{E}}^{(3)}$ is expanded in terms of all interacting triads, with coefficients S_abc being given by (29) $$\begin{array}{c}{S}_{\mathit{abc}}={\int }_0^{z_T}{\int }_0^{L_x}{\int }_0^{L_y}\{\frac{1}{2}{\rho }_a(u_b^{*}{u}_c^{*}+v_b^{*}{v}_c^{*}+w_b^{*}{w}_c^{*})+\frac{p_a{\rho }_b^{*}{\rho }_c^{*}}{2{\rho }_0{p}_0{C}_s^2}+\frac{{\rho }_a{\theta }_b^{*}{\theta }_c^{*}}{2{\rho }_0^2{N}^2}\\ -\frac{{\theta }_a*{\theta }_b^{*}{\theta }_c^{*}}{3N^2{\rho }_0^{2}g}-\frac{4+\kappa }{6{\rho }_0{C}_s^2{p}_0}{p}_a{p}_b^{*}{p}_c^{*}+CP\}d\mathit{xdydz}\end{array}$$(29) where CP means the same term as what is preceding but with cyclic permutations among the subscripts abc. On the other hand, from equation ((24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) 24(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) )(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) it follows (30a) $$E_a\frac{d|A_a{|}^2}{dt}=\sum _b\sum _c{\sigma }_a^{bc}\mathrm{Im}(A_a^{*}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t})$$(30a) (30b) $$\frac{d(A_a{A}_b^{*}{A}_c^{*})}{dt}=\mathcal{O}(A^4)$$(30b)

Taking the time derivative of equations (28a,b)(28a) $${\mathcal{E}}^{(2)}=\sum _a{E}_a|A_a{|}^2$$(28a) and using equation (30(30a) $$E_a\frac{d|A_a{|}^2}{dt}=\sum _b\sum _c{\sigma }_a^{bc}\mathrm{Im}(A_a^{*}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t})$$(30a) ) we get (31) $$\sum _a\sum _b\sum _c[{\sigma }_a^{bc}+S_{\mathit{abc}}{\delta }_{\mathit{abc}}]2i\mathrm{Im}(A_a^{*}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t})+\mathcal{O}(A^4)=0$$(31)

Equation above shows that total pseudo-energy is no longer conserved for an arbitrary truncation of (23(23) $$\mathbf{u}=\left[\begin{array}{c}u\prime (x,y,z,t)\\ v\prime (x,y,z,t)\\ w\prime (x,y,z,t)\\ p\prime (x,y,z,t)\\ \theta \prime (x,y,z,t)\end{array}\right]=\sum _a{A}_a(t){\mathbf{u}}_{\mathbf{a}}+C.C.,$$(23) ).² However, the rate of change of total pseudo-energy for an arbitrary truncation of our model equations is of $O(A^4).$ Thus, the smaller the disturbance amplitude the smaller the variation of total pseudo-energy for a truncated version of the model. Conversely, irrespective of the modal truncation, a necessary condition for equation (31)(31) $$\sum _a\sum _b\sum _c[{\sigma }_a^{bc}+S_{\mathit{abc}}{\delta }_{\mathit{abc}}]2i\mathrm{Im}(A_a^{*}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t})+\mathcal{O}(A^4)=0$$(31) to hold is that the sum inside the square brackets must vanish identically for all interacting triads, namely (32) $${\sigma }_a^{bc}-{\sigma }_b^{ac}-{\sigma }_c^{ab}=-{\delta }_{\mathit{abc}}{S}_{\mathit{abc}}$$(32)

For resonant triads (${\delta }_{\mathit{abc}}=0$) or triads containing only vortical modes, the constraint above reduces to $${\sigma }_a^{bc}-{\sigma }_b^{ac}-{\sigma }_c^{ab}=0.$$

For these kinds of wave triplets, the quadratic component ${\mathcal{E}}^{(2)}$ of total pseudo-energy is exactly conserved. Consequently, for these triads the coupling constants ${\sigma }_a^{bc},{\sigma }_b^{ac}$ and ${\sigma }_c^{ab}$ have always the same sign, and the wave mode with the largest absolute value coupling constant (say, mode a) will always receive energy from or lose energy to the remaining triad components. Moreover, in these resonant interactions, condition ${\sigma }_a^{bc}-{\sigma }_b^{ac}-{\sigma }_c^{ab}=0,$ together with the resonance relation ${\omega }_a-{\omega }_b-{\omega }_c=0,$ implies that this mode having the largest absolute value coupling coefficient will always be the one with the largest absolute eigenfrequency. This mode is usually labeled as pump wave from the context of plasma physics (Weiland and Wilhelmsson, 1977). However, for non-resonant triads the pump mode, or the potentially unstable mode of the triad does not necessarily have the largest absolute coupling coefficient, according to relation (32(32) $${\sigma }_a^{bc}-{\sigma }_b^{ac}-{\sigma }_c^{ab}=-{\delta }_{\mathit{abc}}{S}_{\mathit{abc}}$$(32) ).

Indeed, Table 2 summarises a sample of interacting triads involving acoustic and gravity wave types. It can be noted that the acoustic mode is always the pump wave in near resonant triads involving acoustic and gravity modes (i.e. triads containing two acoustic modes and one gravity wave). In contrast, for interacting triads involving two gravity modes and one acoustic wave, which are characterised by a large frequency mismatch δ_abc, either an acoustic mode (e.g. Triads 3, 5 and 9) or an gravity wave (e.g. Triads 6, 7, 8, 14 and 15) can be the unstable mode (pump mode) of the triplet. However, the large time–frequency mismatch associated with this triad interaction type inhibits the energy exchanges among the mode components, since it requires an unrealistically high amplitude of the pump mode in order for instability to occur, as it will be shown in the next section.

Table 2. Representative examples of interacting triads involving inertia-acoustic (IA) and inertia-gravity modes (IG).

Triad	Mode a	Mode b	Mode c	ωa	ωb	ωc	δabc	$i{\sigma }_a$	$i{\sigma }_b$	$i{\sigma }_c$
1	(1, 1, 1, IA)	(339, 1, 2, IG)	(340, 1, 1, IA)	5.89E–02	4.10E–03	6.30E–02	3.31E–05	–3.11E–01	–2.16E–02	–3.33E–01
2	(170, 1, 2, IG)	(169, 1, 2, IG)	(339, 1, 2, IG)	2.09E–03	2.08E–03	4.10E–03	–6.94E–05	–3.37E + 02	–3.35E + 02	–6.73E + 02
3	(169, 1, 1, IG)	(171, 1, 1, IG)	(340, 1, 1, IA)	3.89E–03	3.93E–03	6.30E–02	5.52E–02	4.47E–01	4.45E–01	2.61E + 00
4	(1, 1, 2, IA)	(6, 1, 1, IA)	(7, 1, 3, IA)	1.11E–01	5.89E–02	1.65E–01	–4.95E–03	–2.53E–03	–1.33E–03	–3.80E–03
5	(2, 1, 1, IG)	(5, 1, 1, IG)	(7, 1, 3, IA)	1.18E–04	1.59E–04	1.65E–01	1.65E–01	–1.55E + 00	–1.55E + 00	–2.38E–04
6	(10, 1, 1, IG)	(329, 1, 2, IG)	(339, 1, 1, IA)	2.58E–04	3.99E–03	6.30E–02	5.87E–02	4.78E + 01	–5.50E + 01	–5.27E + 00
7	(2, 1, 1, IA)	(338, 1, 2, IG)	(340, 1, 1, IG)	5.89E–02	4.09E–03	7.43E–03	–5.55E–02	–1.90E + 00	2.73E + 01	2.58E + 01
8	(5, 1, 1, IG)	(350, 1, 2, IG)	(355, 1, 1, IA)	1.59E–04	4.23E–03	6.34E–02	5.90E–02	4.53E + 01	–5.06E + 01	–5.22E + 00
9	(62, 1, 1, IA)	(107, 1, 3, IG)	(169, 1, 1, IA)	5.90E–02	8.96E–04	5.99E–02	5.09E–09	–4.35E–01	–6.61E–03	–4.42E–01
10	(93, 1, 1, IA)	(421, 1, 1, IG)	(514, 1, 1, IA)	5.92E–02	8.90E–03	6.81E–02	–5.08E–08	1.06E + 00	1.60E–01	1.22E + 00
11	(180, 1, 3, IA)	(478, 1, 2, IG)	(658, 1, 3, IA)	1.66E–01	5.67E–03	1.72E–01	3.36E–08	–1.74E + 00	–5.94E–02	–1.80E + 00
12	(87, 1, 1, IG)	(252, 1, 3, IG)	(339, 1, 2, IG)	2.03E–03	2.09E–03	4.10E–03	–1.61E–05	–1.53E + 03	–1.60E + 03	–3.13E + 03
13	(170, 1, 1, IG)	(169, 1, 1, IG)	(339, 1, 2, IG)	3.91E–03	3.89E–03	4.10E–03	–3.69E–03	–4.82E + 01	–4.64E + 01	–9.50E + 01
14	(1, 1, 1, IA)	(339, 1, 1, IG)	(340, 1, 1, IG)	5.89E–02	7.41E–03	7.43E–03	–5.89E–02	–4.32E–01	–8.10E–01	–1.83E + 00
15	(169, 1, 2, IG)	(–168, 1, 1, IG)	(1, 1, 1, IA)	2.08E–03	3.86E–03	5.89E–02	5.29E–02	–2.23E + 00	–2.45E + 00	2.90E–01

From left to right, the triad members, their respective linear eigenfrequencies, the mismatch ${\delta }_{\mathit{abc}}={\omega }_a-{\omega }_b-{\omega }_c$ and the corresponding coupling constants are given. Each mode is characterised, from left to right, by its zonal, meridional and vertical quantum indexes (j, n, m) and its wave type. The eigenfrequencies ω and mismatchs δ_abc are measured in Hertz. Triads 1,2,4,9,10 and 11 are nearly resonant.

4. Analysis of highly truncated spectral solutions

Given the general theoretical framework on the non-linear interaction among the wave modes of our non-hydrostatic model employed in the previous section, in this section we will further investigate highly truncated versions of the interaction equation (24)(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) to analyse the possibility of acoustic modes to excite inertio-gravity waves. First we will consider the most elementary form of the interaction equations: a single interacting wave triplet. Then we augment our analysis for considering two coupled interacting wave triads.

4.1. Single triad interaction

If one truncates the modal expansion (23(23) $$\mathbf{u}=\left[\begin{array}{c}u\prime (x,y,z,t)\\ v\prime (x,y,z,t)\\ w\prime (x,y,z,t)\\ p\prime (x,y,z,t)\\ \theta \prime (x,y,z,t)\end{array}\right]=\sum _a{A}_a(t){\mathbf{u}}_{\mathbf{a}}+C.C.,$$(23) ) to consider a single interacting triad of modes (a, b, c), equation (24)(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) now reads: (33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) (33b) $$E_{b}A{\prime }_b={\sigma }_b^{ac}{A}_a{A}_c^{*}{e}^{i{\delta }_{\mathit{abc}}t},$$(33b) (33c) $$E_{c}A{\prime }_c={\sigma }_c^{ab}{A}_a{A}_b^{*}{e}^{i{\delta }_{\mathit{abc}}t}.$$(33c)

In the case of a resonant interaction (${\delta }_{\mathit{abc}}=0$), it is well known that if one of the wave modes holds most part of the initial energy of the triplet, this mode will only be unstable if it is the pump wave of the triad (Craik, 1988, Chapter 8). Here, we extend such linear stability analysis for an arbitrary value of the mismatch δ_abc to encompass all the possibilities of interacting triads involving acoustic and gravity modes in this model.

To study the stability of one-wave mode in the single-triad interaction equations above, it is suitable to make the transformation of variable $A_b(t)={\tilde{A}}_b(t)e^{i{\delta }_{\mathit{abc}}t}.$ Inserting this transformation, equations (33a,b,c)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) become (34a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c,$$(34a) (34b) $$E_b(A{\prime }_b+i{\delta }_{\mathit{abc}}{A}_b)={\sigma }_b^{ac}{A}_a{A}_c^{*},$$(34b) (34c) $$E_{c}A{\prime }_c={\sigma }_c^{ab}{A}_a{A}_b^{*},$$(34c) where the ‘$\sim$’ has been omitted to avoid cumbersome notation. Let us assume that Mode a holds almost the total energy of the triad initially, that is, $|A_a(0)|\gg |A_{b,c}(0)|.$ With this assumption, equation (34)(34a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c,$$(34a) can be approximated by their linearized version around the amplitude of Mode a as follows (35a) $$A_a(t)=A_a(0)=\text{const},$$(35a) (35b) $$A\prime {\prime }_b+i{\delta }_{\mathit{abc}}A{\prime }_b=-\frac{{\sigma }_c^{ab}{\sigma }_b^{ac}}{E_c{E}_b}|A_a{|}^2{A}_b.$$(35b)

Thus, as the coupling constants are purely imaginary numbers, in order for instability to occur two conditions must be satisfied: (36a) $$\mathrm{Im}({\sigma }_b^{ac})\mathrm{Im}({\sigma }_c^{ab})>0,\text{and}$$(36a) (36b) $${\delta }_{\mathit{abc}}^2<\frac{4\mathrm{Im}({\sigma }_b^{ac})\mathrm{Im}({\sigma }_c^{ab})}{E_b{E}_c}|A_a{|}^2.$$(36b)

Otherwise, the solution is stable, and no amplification of Modes b and c occurs. Condition (36a) says that Mode a must be the pump mode of the triad. For an exact resonant interaction, this condition is the only requirement for instability to occur. Conversely, according to condition (36b), the minimal value of the amplitude of Mode a for instability to occur increases linearly as the absolute value of the mismatch parameter increases.

Consequently, for an interacting wave triplet composed of two inertio-gravity waves and one inertio-acoustic mode, the amplitude regime of the pump mode required to yield instability might be so high to be observable in the real atmosphere due to the high time–frequency mismatch among the triad components in this case. For example, for the case of Triad 3 of Table 2 composed of mesoscale acoustic and gravity modes, the frequency mismatch is $\frac{{\delta }_{\mathit{abc}}}{2\Omega }\approx 378,$ giving a threshold amplitude associated with vertical wind perturbations of the order of 200 m/s. Thus, for an acoustic mode amplitude value yielding dynamical field perturbations with realistic values, there is no pump wave instability for this triplet. Figure 3 shows the result of a time integration of the full nonlinear version of the three-wave interaction equation (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) for an amplitude value of the pump mode of Triad 3 of Table 2 that yields realistic values of the perturbation field variables. Figure 3 confirms that, for the mode amplitudes associated with realistic values of atmospheric flow disturbances, the high-frequency mismatch associated with the interacting triads involving an acoustic mode and two gravity waves strongly inhibits the energy exchanges among the modes of such triads.

Fig. 3. Time evolution of the mode quadratic energies associated with the solution of the three-wave interaction equations (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) for the modes of Triad 3 of Table 2. The present triad is non-resonant.

In contrast, an acoustic-inertia wave mode can undergo pump wave instability in near resonant triads involving another acoustic-inertia wave and an inertio-gravity mode. In fact, for the parameters of Triad 1 of Table 2 and the value of $\left|A_a\right|$ chosen to yield a vertical wind magnitude of $\approx 2.4$ m/s, parametric instability does occur. Numerical integration of the full equation (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) in this referred amplitude regime illustrates the expressive energy modulations undergone by the three-wave modes (Fig. 4) of this triplet.

Fig. 4. Time evolution of the mode quadratic energies associated with the solution of the three-wave interaction equations (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) for the modes of Triad 1 of Table 2. The present triad is nearly resonant.

Raupp et al. (2019) analysed the dynamics of a single resonant triad interaction with a wave triplet similar to the one illustrated in Fig. 4. They analysed the analytical solution of the interaction equation (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) for the exact resonant case ${\delta }_{\mathit{abc}}=0$ with $|A_a(t=0)|=0$ and discussed the consequences of the mode energy modulations for the physical space solution $[u\prime ,v\prime ,w\prime ,p\prime ,\theta \prime ]$ in view of the energy partition of each mode type on the kinetic, available elastic and available potential forms. Here we will not repeat this analysis and, rather, will stick to the modal space dynamics to further investigate the possibility of an acoustic mode to excite inertio-gravity waves. As in the full system (24(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) ) a single-wave mode may participate of several connected wave triplets, to investigate in a simplified fashion the possibility of each of the wave modes excited by the pump acoustic wave instability shown in Fig. 4 to excite other gravity wave modes, we shall augment our analysis of the phase-space dynamics to consider two triads coupled by a single mode.

4.2. Two triads coupled by one mode

Let us now consider a truncated version of modal expansion (23(23) $$\mathbf{u}=\left[\begin{array}{c}u\prime (x,y,z,t)\\ v\prime (x,y,z,t)\\ w\prime (x,y,z,t)\\ p\prime (x,y,z,t)\\ \theta \prime (x,y,z,t)\end{array}\right]=\sum _a{A}_a(t){\mathbf{u}}_{\mathbf{a}}+C.C.,$$(23) ) that considers five modes (a,b,c,d,e) whose wavenumbers and eigenfrequencies satisfy the relations (37a) $$k_a=k_b+k_c,$$(37a) (37b) $$l_a=l_b+l_c,$$(37b) (37c) $${\lambda }_a=\pm {\lambda }_b\pm {\lambda }_c,$$(37c) (37d) $${\omega }_a={\omega }_b+{\omega }_c+{\delta }_{\mathit{abc}},$$(37d) (37e) $$k_c=k_d+k_e,$$(37e) (37f) $$l_c=l_d+l_e,$$(37f) (37g) $${\omega }_c={\omega }_d+{\omega }_e+{\delta }_{\mathit{cde}}.$$(37g)

Notice that we have imposed the interaction condition for the vertical wavenumbers only to the primary triad (a,b,c). This is because the condition for the vertical wavenumbers is no longer excluding, as previously discussed, so that other modes whose coupling constants are non-zero, although small, can be excited by the mode coupling the two triads. In this situation, equation ((24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) 24(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) )(24) $$E_{a}A{\prime }_a(t)=\sum _b\sum _c{\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t}$$(24) now reads (38a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(38a) (38b) $$E_{b}A{\prime }_b={\sigma }_b^{ac}{A}_a{A}_c^{*}{e}^{i{\delta }_{\mathit{abc}}t},$$(38b) (38c) $$E_{c}A{\prime }_c={\sigma }_c^{ab}{A}_a{A}_b^{*}{e}^{i{\delta }_{\mathit{abc}}t}+{\sigma }_c^{de}{A}_d{A}_e{e}^{-i{\delta }_{\mathit{cde}}t},$$(38c) (38d) $$E_{d}A{\prime }_d={\sigma }_d^{ce}{A}_c{A}_e^{*}{e}^{i{\delta }_{\mathit{cde}}t},$$(38d) (38e) $$E_{e}A{\prime }_e={\sigma }_e^{cd}{A}_c{A}_d^{*}{e}^{i{\delta }_{\mathit{cde}}t}.$$(38e)

Thus, to analyse the stability of the triad interaction (a,b,c) to small amplitude perturbations associated with the interacting modes (c,d,e), let us first assume that $|A_c(t=0)|\gg |A_{d,e}(t=0)|.$ In this case, Mode c evolves independently of Modes d and e, and its amplitude obeys the three-wave equations (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) . Furthermore, by explicitly expressing A_c, A_d and A_e in terms of their real and imaginary parts, the linearized version of equation (38d(38d) $$E_{d}A{\prime }_d={\sigma }_d^{ce}{A}_c{A}_e^{*}{e}^{i{\delta }_{\mathit{cde}}t},$$(38d) ,e) around the amplitude of Mode c can be written as (39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) where $\mathbf{X}={[A_{d_r},A_{d_i},A_{e_r},A_{e_i}]}^T$ and the matrix M(t) is defined by (40) $$M(t)=\left[\begin{array}{cc}{0}_{2\times 2}& M_d(t)\\ M_e(t)& 0_{2\times 2}\end{array}\right]$$(40) with $M_d(t)$ being defined as (41) $$M_d(t)=\left[\begin{array}{cc}-\frac{\mathrm{Im}({\sigma }_d^{ce})}{E_d}(A_{c_r}\hspace{0.17em}\text{sin}({\delta }_{\mathit{cde}}t)+A_{c_i}\hspace{0.17em}\text{cos}({\delta }_{\mathit{cde}}t))& \frac{\mathrm{Im}({\sigma }_d^{ce})}{E_d}(-A_{c_i}\hspace{0.17em}\text{sin}({\delta }_{\mathit{cde}}t)+A_{c_r}\hspace{0.17em}\text{cos}({\delta }_{\mathit{cde}}t))\\ \frac{\mathrm{Im}({\sigma }_d^{ce})}{E_d}(A_{c_r}\hspace{0.17em}\text{cos}({\delta }_{\mathit{cde}}t)-A_{c_i}\hspace{0.17em}\text{sin}({\delta }_{\mathit{cde}}t))& \frac{\mathrm{Im}({\sigma }_d^{ce})}{E_d}(A_{c_i}\hspace{0.17em}\text{cos}({\delta }_{\mathit{cde}}t)+A_{c_r}\hspace{0.17em}\text{sin}({\delta }_{\mathit{cde}}t))\end{array}\right]$$(41)

If the matrix coefficient M(t) is periodic, Floquet theorem (Arnol’d, 1989; Vanneste and Vial, 1994; Majda, 2003) can be used to study the stability of system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ). Conversely, as we are considering arbitrary time–frequency mismatchs δ_abc and δ_cde, in order for M(t) to be exactly periodic (and, therefore, Floquet theory be applicable), two conditions must be met: (i) the mismatch δ_abc and the non-linear oscillation frequency of the spectral amplitudes A_a, A_b and A_c must be co-measurable; and (ii) the resulting oscillation frequency of $A_c(t)$ be also co-measurable with the second triad mismatch δ_cde. As these two conditions are very restrictive, a more general way to analyse the stability of the aforementioned linear system is to estimate its maximal Lyapunov exponent (MLE) (Zounes and Rand, 1998), for which the MLE being positive means instability.³ The MLE has been evaluated using the method described in Benettin et al. (1976), whose implementation is available in (Datseris, 2018).

For the five-wave system composed of the modes of Triads 1 and 2 of Table 2, the MLE of system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) associated with the linearised dynamics of the gravity modes (170, 1, 2) and (169, 1, 2) is ${\Lambda }_L=1.35\times {10}^{-5}{s}^{-1},$ corresponding to a growth rate of $1/20.5h^{-1},$ which seems compatible with the typical time scale of internal-gravity waves. In this case, the time modulation of $A_c(t)$ refers to the solution shown in Fig. 4. The instability of the gravity wave-mode (339, 1, 2) to the other two gravity modes of Triad 2 of Table 2 is illustrated in the numerical integration of system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) shown in Fig. 5, which shows the growth of the gravity wave harmonics (170, 1, 2) and (169, 1, 2). Since the gravity mode (339, 1, 2) representing Mode c in this example is the pump wave of the triplet (c,d,e), this instability appears more likely to be a pump wave instability than a modulational instability explored by Connaughton et al. (2010) in the Rossby wave context.

Fig. 5. Numerical solution of the linearized system (39) composed of modes (170,1,2,IG) and (169,1,2,IG) of Triad 2. These modes are parametrically forced by the mode (339, 1, 2, IG) of Triad 1. This solution presents a maximal Lyapunov exponent ${\lambda }_L=1.35\times {10}^{-5}{s}^{-1}.$

In fact, if we consider the coupling between Triads 1 and 13 of Table 2 through the same gravity mode (339, 1, 2), that is, if we replace Triad 2 to Triad 13 (which has a larger time–frequency mismatch) in our five-wave system, the solution of (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) is stable (MLE ${\Lambda }_L=-4.85\times {10}^{-10}{s}^{-1}$) and no amplification of a pair of other gravity waves occurs (Fig. 6). Figure 7 shows a similar situation when Modes d and e are gravity modes and Mode c is an inertio-acoustic mode. In this case, the planetary acoustic mode of Triad 1 of Table 2 couples this triad with the two inertio-gravity modes of Triad 15, and the MLE of system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) is positive but very close to zero. In fact, the MLE in this case is ${\Lambda }_L=6.59\times {10}^{-11}{s}^{-1},$ corresponding to a growth rate of $\approx 1/180da{y}^{-1},$ which is no longer compatible with the observed time scale of gravity waves. Therefore, in this situation, as in the previous case shown in Fig. 6, system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) is stable and there is no amplification of the two gravity modes. The numerical integration of the non-linear five-wave system (38(38a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(38a) ) for the corresponding stable cases of the linearized system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) shown in Figs. 6 and 7 confirms that there is no energy leakage from Triad 1 towards the other gravity modes, so that the time evolution of Modes (a,b,c) is identical to that predicted by the three-wave problem shown in Fig. 4 (figures not shown).

Fig. 6. Numerical solution of the linearized system (39) composed of modes (170,1,1,IG) and (169,1,1,IG) of Triad 13. These modes are parametrically forced by the mode (339, 1, 2, IG) of Triad 1. This solution presents a maximal Lyapunov exponent ${\lambda }_L=-4.85\times {10}^{-10}{s}^{-1}.$

Fig. 7. Numerical solution of the linearized system (39) composed of modes (169,1,2,IG) and (–168,1,1,IG) of Triad 15. These modes are parametrically forced by the mode (1, 1, 1, IA) of Triad 1. This solution presents a maximal Lyapunov exponent ${\lambda }_L=6.59\times {10}^{-11}{s}^{-1}.$

A result of numerical integration of the full five-wave system for the unstable case of system (39(39) $$\mathbf{X}\prime =M(t)\mathbf{X}$$(39) ) is shown in Figs. 8 and 9. In this numerical integration, the initial condition is similar to those considered in Fig. 4, in which the pump acoustic wave (340, 1, 1) holds almost total initial energy of the system, with only a small perturbation distributed among the remaining modes b, c, d and e. From Figs. 8 and 9 one notices that, after Mode a excites Modes b and c through pump wave instability, part of Mode c energy leaks to the gravity modes (170, 1, 2) and (169, 1, 2) (Modes d and e), yielding the excitation of these modes. In this case, Mode c acts as the pump mode in the triad interaction with Modes d and e. Consequently, the time evolution of the mode energies of the five-wave system in this case shows that, apart from the considerable energy modulations of the gravity modes d and e at the expense of Mode c, the coupling of Mode c with the two gravity modes d and e also yields a multiplication of the periods associated with the time evolution of the energies of Modes a and b, disturbing the energy modulations of these modes from their exactly periodic nature predicted by the three-wave equations (33)(33a) $$E_{a}A{\prime }_a={\sigma }_a^{bc}{A}_b{A}_c{e}^{-i{\delta }_{\mathit{abc}}t},$$(33a) .

Fig. 8. Numerical solution of the five-wave system (38) composed of the modes of Triads 1 and 2 of Table 2. This figure illustrates the time evolution of the quadratic energies corresponding to Modes of Triad 1 only.

Fig. 9. This is same as Fig. 8, but illustrating the quadratic energies of the secondary gravity modes of Triad 2.

Another interesting feature regarding the unstable case is that the eigenfrequencies of the gravity modes (170, 1, 2) and (169, 1, 2) representing Modes d and e in this example are such that the linear wave dynamics of these modes is well described by hydrostatic approximation, as well as the non-linear interaction involving these inertio-gravity waves and vortical modes. As the non-linear interaction involving inertio-gravity waves and vortical modes has been demonstrated to play an important role in the non-linear geostrophic adjustment (Majda and Embid, 1998; Vanneste, 2004; Vanneste and Yavneh, 2004), the theoretical results described here suggest that acoustic modes might be important for both hydrostatic and geostrophic adjustment processes in the atmosphere.

5. Summary and conclusions

Non-linear triad interactions involving inertio-acoustic and inertio-gravity waves are studied here in the context of the mid-latitude f-plane shallow non-hydrostatic equations for a background state at rest and characterised by a hydrostatic balance and an isothermal temperature profile. In this context, we have adopted highly truncated Galerkin expansions in terms of the eigensolutions of the linear problem. For a single-triad interaction, we have shown that the interacting triplets involving two inertio-gravity waves and one acoustic mode require a likely unrealistic modal amplitude regime in order for pump wave instability to occur. Consequently, for direct triad interactions, we have shown that an inertio-acoustic wave mode can only be unstable to perturbations associated with a pair of acoustic/gravity modes.

In contrast, the analysis of the dynamics of two triads coupled by a single-wave mode shows that a non-hydrostatic inertio-gravity wave mode (i.e. having an eigenfrequency such that the non-hydrostatic effect of vertical acceleration is not negligible) participating of a nearly resonant interaction with two acoustic modes is unstable in nearly resonant triad interactions with a pair of lower frequency inertio-gravity waves. In fact, for the representative example illustrated here in which the eigenfrequencies of the two excited inertio-gravity waves are nearly a half of the time–frequency of the primary wave (which couples the two triads), the numerical results of the non-linear dynamics of the five-wave system confirm that this instability yields energy modulations on the two secondary gravity modes. Since the higher the time frequency the more important the non-hydrostatic effect of vertical acceleration on the inertio-gravity waves, the results suggest that inertio-acoustic waves may induce hydrostatically balanced inertio-gravity waves to undergo episodic amplitude (energy) modulations due to inter-triad energy exchanges. On the other hand, if one of the acoustic modes of a resonant interaction involving acoustic/gravity waves couple the two triads in our reduced five-wave dynamics, our results show that the maximum Lyapunov exponent of the corresponding linearized system for the two gravity waves gives a growth time-scale of $\approx 180$ days, which is no longer compatible with the observed time-scale of internal gravity waves and, consequently, the resulting instability might likely be irrelevant in the atmosphere.

As discussed in Section 1, due to the ultra-high frequency of acoustic waves, their numerical treatment using explicit schemes implies highly restrictive computational constraints. Consequently, in non-hydrostatic numerical models adopted for meso-scale simulations the acoustic modes are either filtered out or subjected to strong damping associated with implicit schemes having time-steps much higher than the acoustic cut-off period. However, acoustic waves play an important role in the hydrostatic adjustment process as they are responsible for the vertical displacements of fluid parcels associated with the expansion of an instantaneously heated atmospheric layer (Bannon, 1995, 1996; Duffy, 2003). In fact, Chagnon and Bannon (2001) demonstrated that the steady-state solutions of anelastic and other sound filtering models exhibit significant differences from those of fully compressible models allowing acoustic modes.

Acoustic waves may be excited by thermal forcings associated with convective storms, especially the localized ones that have a duration shorter than the acoustic cut-off period given by $\frac{4\pi H}{C_S}$ (Chagnon and Bannon, 2005a, 2005b). Our simplified theoretical model suggests that these acoustic modes generated by explosive and localised storms might play an important role in the transient phase of the three-dimensional adjustment process of the atmosphere to both hydrostatic and geostrophic balances. Specifically, this role of acoustic modes in the adjustment process of the atmosphere might be due to not only their linear energy propagation as studied by Chagnon and Bannon (2005a, 2005b) but also their non-linear effect of exciting hydrostatic inertio-gravity waves as pointed out by our theoretical analysis. Fanelli and Bannon (2005) investigated the hydrostatic and geostrophic adjustments to a prescribed thermal forcing utilising a non-linear compressible model, but they considered a heating function with a duration longer than the acoustic cut-off period of $\approx 5$ min so that no acoustic waves were excited. However, similar numerical studies with a shorter time-scale forcing should be done to apply this theory to the hydrostatic and geostrophic adjustments in a more realistic fashion by considering the full expansion (23(23) $$\mathbf{u}=\left[\begin{array}{c}u\prime (x,y,z,t)\\ v\prime (x,y,z,t)\\ w\prime (x,y,z,t)\\ p\prime (x,y,z,t)\\ \theta \prime (x,y,z,t)\end{array}\right]=\sum _a{A}_a(t){\mathbf{u}}_{\mathbf{a}}+C.C.,$$(23) ). This might be the next step in the generalisation of this theory to further understand both the non-linear dynamics of the non-hydrostatic wave modes itself and its role in the non-linear hydrostatic/geostrophic adjustment, along with testing the robustness of this theory.

Acknowledgements

The authors thank Dr. Breno Raphaldini and Dr. David Ciro Taborda for their helpful suggestions. The authors also acknowledge the insightful and constructive comments of the anonymous reviewer.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

The work of André Teruya was supported by CAPES (Coordenação de Aperfeiçoamento de Pessoal de Ensino Superior) and CNPQ (Conselho Nacional de Desenvolvimento Científico e Tecnológico. The work of Carlos F. M. Raupp was supported by FAPESP (Fundação de Amparo à Pesquisa do Estado de São Paulo) through grants 2009/11643-4 and 2015/50686-1. The authors also acknowledge the support from CAPES IAG/USP PROEX (grant number 0531/2017).

Notes

1 r = 0 for a vortical mode and $r=\pm 1$ (±2) for westward/eastward inertio-gravity (inertio-acoustic) modes, for example.

2 In reality, total pseudo-energy of model equations (1)–(5) is not exactly conserved even for the full expansion case due to the $\mathcal{O}{(u\prime ,v\prime ,w\prime ,p\prime ,\theta \prime )}^4$ terms in equation (10)(10) $${\mathcal{E}}^{(2)}+{\mathcal{E}}^{(3)}+\mathcal{O}{(u\prime ,v\prime ,w\prime ,p\prime ,\theta \prime )}^4=\text{const}$$(10) .

3 To estimate the MLE (Λ_L) of $\mathbf{X}\prime =M(t)\mathbf{X},$ we obtain two solutions whose initial conditions differ by a small separation, namely ${\mathbf{X}}_1(t)$ and ${\mathbf{X}}_2(t)\equiv {\mathbf{X}}_1(t)+\delta (t),$ with $|\delta (0)|=|{\mathbf{X}}_1(0)-{\mathbf{X}}_2(0)|.$ Then, the MLE is defined by the following equation $|\delta (t)|=|\delta (0)|\hspace{0.17em}\text{exp}\hspace{0.17em}({\Lambda }_{L}t),$ representing a measure of the exponential rate of separation between two neighbouring trajectories.