1. Introduction
As in Scott & Cambon (Reference Scott and Cambon2024) and Scott (Reference Scott2025) (henceforth referred to as [I] and [II]), this paper concerns weak (small Rossby or Froude number) unbounded, homogeneous turbulence subject to rotation and stable stratification with constant Brunt–Väisälä frequency,
$N$
, and a constant and vertical rotation vector,
$\boldsymbol{\varOmega }$
, whose norm gives the rotation rate,
${\varOmega} =| \boldsymbol{\varOmega }|$
. The approach used is a spectral one, i.e. based on Fourier transforms, and a general overview of homogeneous turbulence with this viewpoint in mind can be found in the book by Sagaut & Cambon (Reference Sagaut and Cambon2018). The difference compared with [I] and [II] is that we focus on the case
$\beta =2{\varOmega} /N$
close to
$1$
, for which the weak-turbulence analysis used in the previous articles is no longer applicable.
Different approaches have been used to address the problem of homogeneous turbulence with rotation and/or stratification. These include the eddy-damped, quasi-normal, Markovian (EDQNM) closure (see Cambon & Jacquin (Reference Cambon and Jacquin1989), Godeferd & Cambon (Reference Godeferd and Cambon1994), Cambon, Mansour & Godeferd (Reference Cambon, Mansour and Godeferd1997) and Godeferd & Staquet (Reference Godeferd and Staquet2003) for studies of turbulence with both rotation and stratification). Another is direct numerical simulation (DNS), which approximates the flow as periodic using truncated Fourier series (see e.g. Orszag & Patterson Reference Orszag and Patterson1972; Canuto et al. Reference Canuto, Hussaini, Quarteroni and Zang1988). Direct numerical simulation studies have been undertaken of purely stratified turbulence (see e.g. Waite & Bartello Reference Waite and Bartello2004, Reference Waite and Bartello2006a
) and of the purely rotating case (see e.g. Smith & Waleffe Reference Smith and Waleffe1999; Bourouiba & Bartello Reference Bourouiba and Bartello2007). Of course, pure rotation or stratification does not include the case of small
$\beta -1$
considered here. Applications of DNS to stratified, rotating turbulence can be found in Coleman, Ferziger & Spalart (Reference Coleman, Ferziger and Spalart1992), Smith & Waleffe (Reference Smith and Waleffe2002), Liechtenstein, Godeferd & Cambon (Reference Liechtenstein, Godeferd and Cambon2005), Waite & Bartello (Reference Waite and Bartello2006b
), Thomas & Daniel (Reference Thomas and Daniel2020, Reference Thomas and Daniel2021) and Thomas (Reference Thomas2023). Several of these papers discuss interactions and energy transfer between wave and non-propagating (NP) modes, these being central themes of the present work. This being said, none of those studies are specifically aimed at the case of weak turbulence with small
$\beta -1$
studied here.
As discussed in the introduction of [I], in the case of weak turbulence, inertial-gravity wave oscillations lead to problems for classical DNS. Firstly, the time step must be sufficiently small as to resolve the oscillations, while the total time for evolution needs to be long enough that nonlinearity is significant. Since nonlinearity requires many wave periods for weak turbulence, the required number of time steps is considerably more that in the absence of waves. Secondly, DNS involves discretisation in spectral space. Wave oscillations in that space become increasingly rapid as time increases and the spectral resolution should be sufficiently fine that it resolves the oscillations at the times for which nonlinearity intervenes. These constraints mean that accuracy of classical DNS requires more and more computer resources as the turbulence weakens, hence the interest of studies, such as the present one, which deliberately focus on the weak-turbulence limit.
With weak turbulence in mind, the analysis of [I] expressed the flow as a sum over modes, which are solutions of the governing equations without nonlinearity or visco-diffusion (see Bartello (Reference Bartello1995) and [I]). Modes have the form
$\exp (i(\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}-s\omega (\boldsymbol{k})t))$
, where
$\omega (\boldsymbol{k})$
is the dispersion relation of inertial-gravity waves. They are indexed by the wave vector
$\boldsymbol{k}$
and the integer
$s=0,\pm 1$
. Those with
$s=\pm 1$
are inertial-gravity waves, while
$s=0$
modes have zero frequency, hence zero group velocity when regarded as waves. Thus, the latter modes are NP modes. They are often called vortical or geostrophic modes in the geophysical literature, where vortical refers to potential vorticity (see e.g. Pedlosky Reference Pedlosky1987). In any case, there are two components of the flow, wave and NP, whose physics is rather different in the weak-turbulence limit.
Mode-amplitude equations were derived in [I]. These equations describe the temporal evolution of the mode amplitudes once nonlinearity and visco-diffusion are allowed for. A spectral matrix, whose diagonal elements give the energy density of the different modes in spectral space, was defined, the total energy being a sum of wave and NP contributions. This analytical framework forms the basis of the present study and is summarised in § 2.
Assuming weak turbulence and
$\beta =2{\varOmega} /N$
not too close to
$1$
, it was shown in [I], § 3, that the NP component decouples from the wave one at leading order and a modified version of DNS was developed for that component. On the other hand, the wave component was described by wave-turbulence analysis, which has a long and ongoing history (see e.g. Benney & Saffman Reference Benney and Saffman1966; Benney & Newell Reference Benney and Newell1969; Zhakarov, Lvov & Falkovich Reference Zakharov, Lvov and Falkovich1992; Newell & Rumpf Reference Newell and Rumpf2011; Nazarenko Reference Nazarenko2011; Galtier Reference Galtier2022). To close the wave-turbulence equations, it was assumed in [I] that the amplitude of the NP component was small compared with the wave one. This restriction was lifted in [II], which allowed for comparable amplitudes of the two components. In that case, the NP component was found to induce coupling between wave modes, but itself still evolved as if the wave component did not exist.
However, both [I] and [II] assume that
$\beta -1$
is not too small, otherwise the wave modes are only weakly dispersive and wave-turbulence analysis, which assumes dispersive waves, no longer applies. The case in which
$\beta -1$
is small or zero is the subject of this article and wave-turbulence analysis makes no appearance. Instead, we derive weak-turbulence equations for the mode amplitudes and examine their consequences, at first analytically, then using numerical integration.
As in [II], we find that, owing to vanishing of a particular NP–wave coupling coefficient for
$\beta =1$
, the NP component is decoupled from the wave one at leading order, though the latter is affected by the former. We should, however, note that, being based on leading-order analysis, the theory describes what happens up to times of
$O(\varepsilon ^{-1})$
times the wave period, where
$\varepsilon$
is a small parameter characterising the amplitude of turbulence (the Rossby or Froude number). The scaling
$\varepsilon ^{-1}$
expresses the time for nonlinear evolution, but it may be that, as suggested by Thomas (Reference Thomas2023), NP decoupling does not persist to larger times.
The approach used here is, in the end, numerical and similar to classical DNS, though weak-turbulence analysis is used to suppress the oscillations due to inertial-gravity waves. As a result, it avoids the problems of direct application of classical DNS referred to above.
This paper is organised as follows. In § 2, the analytical framework using modal projection is summarised and weak-turbulence approximations of the mode-amplitude equations for small
$\beta -1$
derived. Decoupling of the NP component from the wave component is shown in Appendix A, while energy conservation by nonlinearity is the subject of Appendix B. Finally, Appendix C describes the numerical procedure used to integrate the mode-amplitude equations and results are given in § 3.
2. Analytical treatment
Consider homogeneous turbulence in a rotating, stably stratified fluid having constant Brunt–Väisälä frequency,
$N$
, and constant rotation vector
$\boldsymbol{\varOmega }$
, which is supposed parallel to gravity and has unit vector
$\boldsymbol{e}=\boldsymbol{\varOmega }/{\varOmega}$
, where
${\varOmega} =| \boldsymbol{\varOmega }|$
. Henceforth, spatial coordinates, time and velocity are non-dimensionalised using
$L ,\ (N^{2}+4{\varOmega} ^{2})^{-1/2}$
and
$L(N^{2}+4{\varOmega} ^{2})^{1/2}$
, where
$L$
is a length scale characterising the large scales of the initial turbulence.
Using a rotating frame of reference and Cartesian coordinates, as well as the summation convention, the non-dimensional Boussinesq equations of motion are
\begin{align} \frac{\partial u_{i}}{\partial x_{i}}&=0 , \end{align}
where 
and
$\eta$
denote velocity, pressure and the buoyancy variable, while
$\varepsilon _{\textit{ijk}}$
is the alternating tensor. Furthermore,
$D_{u}=\nu (N^{2}+4{\varOmega} ^{2})^{-1/2}L^{-2}$
and
$D_{\eta }=\kappa (N^{2}+4{\varOmega} ^{2})^{-1/2}L^{-2}$
, where
$\nu$
is the kinematic viscosity and
$\kappa$
the diffusivity associated with the buoyancy variable
$\eta$
. Note that 
and 
, where 
. Thus, the non-dimensional governing equations, (2.1)–(2.3), only depend on the parameters
$\beta , D_{u}$
and
$D_{\eta }$
. For simplicity’s sake, 
are denoted
${\varOmega} , N$
in what follows. Because we will only be working with non-dimensional quantities, this should not lead to confusion. Note that the right-hand sides of (2.1) and (2.3) express nonlinearity, viscosity and diffusion. The visco-diffusive terms dissipate energy.
As described in [I], in the absence of nonlinearity and visco-diffusion, (2.1)–(2.3) have particular solutions, referred to as modes, which are indexed by an integer
$s=0,\pm 1$
and a wave vector
$\boldsymbol{k}$
. Modes have the velocity and buoyancy variable
\begin{align} u_{i}&=v_{i}^{(s)}(\boldsymbol{k})\exp \big(i\big(\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}-s\omega (\boldsymbol{k})t\big)\big) ,\end{align}
\begin{align} \eta &=\eta ^{(s)}(\boldsymbol{k})\exp \big(i\big(\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}-s\omega (\boldsymbol{k})t\big)\big) ,\end{align}
where
\begin{equation} \omega (\boldsymbol{k})=\big(N^{2}\sin ^{2}\theta _{\boldsymbol{k}}+{\unicode{x1D6FA}} ^{2}\cos ^{2}\theta _{\boldsymbol{k}}\big)^{1/2} ,\end{equation}
$\theta _{\boldsymbol{k}}$
is the angle between
$\boldsymbol{k}$
and the rotation vector and
$v_{i}^{(s)} , \eta ^{(s)}$
are given by (I.2.12)–(I.2.15) (where, here and henceforth, (I.x.y) refers to equation (x.y) of [I]).
Although modes do not allow for nonlinearity or visco-diffusion, they form a complete set for
$u_{i}(\boldsymbol{x}) , \eta (\boldsymbol{x})$
satisfying (2.2). Thus, the flow can be expressed at any given time as
\begin{equation} u_{i}=\sum _{s=0,\pm 1}u_{i,s} ,\quad \eta =\sum _{s=0,\pm 1}\eta _{s} ,\end{equation}
where
\begin{align} u_{i,s}&=\int\! a_{s}\left(\boldsymbol{k},t\right)v_{i}^{(s)}(\boldsymbol{k})\exp \big(i\big(\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}-s\omega (\boldsymbol{k})t\big)\big)\textrm{d}^{3}\boldsymbol{k} ,\end{align}
\begin{align} \eta _{s}&=\int\! a_{s}\left(\boldsymbol{k},t\right)\eta ^{(s)}(\boldsymbol{k})\exp \big(i\big(\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}-s\omega (\boldsymbol{k})t\big)\big)\textrm{d}^{3}\boldsymbol{k} \end{align}
and
$a_{s}$
are the mode amplitudes at time
$t$
. The
$a_{s}$
evolve in time due to nonlinearity and visco-diffusion.
Equations (2.7)–(2.9) express the flow as a linear combination of modes indexed by
$s=0,\pm 1$
and the wave vector
$\boldsymbol{k}$
, each mode being a solution of the governing equations without visco-diffusion and nonlinearity of the form
$\exp (i(\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}-s\omega (\boldsymbol{k})t))$
. The modes with
$s=\pm 1$
are inertial-gravity waves, whose dispersion relation is (2.6). The non-dimensionalisation of time makes the period of oscillation of such waves of order
$1$
. The mode
$s=0$
has zero frequency, hence zero group velocity. The corresponding component of the flow is non-propagating, hence its designation by NP. As in [I], since
$u_{i}$
and
$\eta$
are real,
$a_{s}(-\boldsymbol{k})=a_{-s}^{*}(\boldsymbol{k})$
, where
$*$
denotes complex conjugation. Using the identities
$v_{i}^{(s)}(-\boldsymbol{k})={v_{i}^{(-s)}}^*(\boldsymbol{k}) ,\ \eta ^{(s)}(-\boldsymbol{k})={\eta ^{(-s)}}^*(\boldsymbol{k})$
and
$\omega (-\boldsymbol{k})=\omega (\boldsymbol{k})$
, (2.8) and (2.9) imply
$u_{i,s}=u_{i,-s}^{*} ,\ \eta _{s}=\eta _{-s}^{*}$
. Thus, the NP component, defined by
$u_{i,0}$
and
$\eta _{0}$
, is real, whereas
$s=+1$
and
$s=-1$
provide conjugate contributions to (2.7), whose sum is the wave component and is also real.
Since we are dealing with turbulent flow, a statistical description is to be expected. In particular averaging plays an important role. The mean flow is supposed zero, hence
$\overline{a_{s}}(\boldsymbol{k})=0$
, where the overbar denotes an ensemble average. Given statistical homogeneity, the second-order moments have the form
\begin{equation} \overline{a^{*}_{s} (\boldsymbol{k})a_{s^{\prime}}\left(\boldsymbol{k}'\right)}=A_{{\textit{ss}}^{\prime}}(\boldsymbol{k})\delta \left(\boldsymbol{k}-\boldsymbol{k}'\right)\! ,\end{equation}
where
$A_{{\textit{ss}}^{\prime}}$
is the spectral matrix, which is Hermitian and positive semi-definite. Given
$a_{s}(-\boldsymbol{k})=a_{-s}^{*}(\boldsymbol{k}) ,\ A_{{\textit{ss}}^{\prime}}(-\boldsymbol{k})=A_{-s^{\prime},-s}(\boldsymbol{k})$
. The average, non-dimensional energy density in physical space is the sum of kinetic and potential energies:
\begin{equation} \underset{\text{Kinetic}}{\underbrace{\frac{1}{2}\overline{u_{i}u_{i}}}}+\underset{\text{Potential}}{\underbrace{\frac{1}{2}\overline{\eta ^{2}}}} .\end{equation}
As shown in [I], (2.11) can be expressed as
\begin{equation} \frac{1}{2}\big(\overline{u_{i}u_{i}}+\overline{\eta ^{2}}\big)=\int\! e(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} ,\end{equation}
where
$e(\boldsymbol{k})$
is the energy density in spectral space and is given by
\begin{equation} e=\frac{1}{2}\sum _{s=0,\pm 1}A_{\textit{ss}} .\end{equation}
The diagonal elements of the spectral matrix, which appear in (2.13), are real, non-negative and express the energy densities of the different modes in spectral space. Given
$A_{{\textit{ss}}^{\prime}}(-\boldsymbol{k})=A_{-s^{\prime},-s}(\boldsymbol{k})$
, they are such that
$A_{00}(-\boldsymbol{k})=A_{00}(\boldsymbol{k})$
and
$A_{-1,-1}(-\boldsymbol{k})=A_{11}(\boldsymbol{k})$
. On the other hand, the off-diagonal elements of
$A_{{\textit{ss}}^{\prime}}$
can be complex and represent correlations between modes. The energy densities of the wave and NP components are
\begin{equation} e_{W}=\frac{1}{2}\sum _{s=\pm 1}A_{\textit{ss}} ,\quad e_{\textit{NP}}=\frac{1}{2}A_{00} ,\end{equation}
whose sum provides
$e$
. The energy densities are such that
$e_{W}(-\boldsymbol{k})=e_{W}(\boldsymbol{k}) ,\ e_{\textit{NP}}(-\boldsymbol{k})=e_{\textit{NP}}(\boldsymbol{k})$
and
$e(-\boldsymbol{k})=e(\boldsymbol{k})$
.
The mode amplitudes evolve according to (I.2.26), i.e.
where
\begin{equation} \tilde{f}(\boldsymbol{k})=\frac{1}{8\pi ^{3}}\int\! f\left(\boldsymbol{x}\right)\exp \left(-i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}\right)\textrm{d}^{3}\boldsymbol{x} \end{equation}
denotes spatial Fourier transformation and
$D_{s\hat{s}}$
is the Hermitian, positive-definite matrix defined by (I.2.20). The first term on the right-hand side of (2.15) represents nonlinearity, whereas the second expresses visco-diffusive effects, which are dissipative.
Assuming, as we do from here on, that the turbulence is weak (small Rossby or Froude number), the nonlinear term in (2.15) is small and many wave periods (large
$t$
) are required for it to be effective. The visco-diffusive term must then also be small (i.e. small
$D_{s\hat{s}}$
), otherwise dissipation will kill the turbulence before nonlinearity intervenes. Terms in the final sum of (2.15) are then small and oscillatory when
$\hat{s}\neq s$
and thus have little effect on the long-time evolution considered here. Neglecting such terms, (2.15) becomes
Since the matrix
$D_{s\hat{s}}$
is Hermitian and positive-definite, the diagonal elements appearing in (2.17) are real and positive. They can be written as
$D_{00}=D_{\textit{NP}} , D_{11}=D_{-1,-1}=D_{W}$
, where
$D_{\textit{NP}}$
and
$D_{W}$
are the NP- and wave-mode damping factors. It remains to consider the nonlinear term, i.e. the right-hand side of (2.17).
From here on, we suppose small
$\beta -1$
and derive the leading-order amplitude equations. Given
$N=(\beta ^{2}+1)^{-1/2}$
and
${\varOmega} =\beta (\beta ^{2}+1)^{-1/2} ,\ \beta =1$
leads to
$\omega (\boldsymbol{k})=2^{-1/2}$
according to (2.6). This value, independent of
$\boldsymbol{k}$
, means that the wave modes are non-dispersive when
$\beta =1$
, hence wave-turbulence theory does not apply. For small
$\beta -1$
, (2.6) implies
\begin{equation} \omega (\boldsymbol{k})\sim 2^{-1/2}\left(1+\frac{1}{2}\left(\beta -1\right)\cos 2\theta _{\boldsymbol{k}}\right)\! ,\end{equation}
which gives the constant value
$\omega =2^{-1/2}$
, noted above, when
$\beta =1$
and shows the small departures from that value when
$| \beta -1|$
is small. Using (2.18) in (2.8) and (2.9) gives
\begin{equation} u_{i,s}=w_{i,s}\exp \big(-ist/2^{1/2}\big),\quad \eta _{s}=\zeta _{s}\exp \big(-ist/2^{1/2}\big), \end{equation}
where
\begin{align} w_{i,s}&=\int\! c_{s}(\boldsymbol{k})v_{i}^{(s)}(\boldsymbol{k})\exp \left(i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}\right)\textrm{d}^{3}\boldsymbol{k} ,\end{align}
\begin{align} \zeta _{s}&=\int\! c_{s}(\boldsymbol{k})\eta ^{(s)}(\boldsymbol{k})\exp \left(i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}\right)\textrm{d}^{3}\boldsymbol{k} ,\end{align}
\begin{align} c_{s}(\boldsymbol{k})&=a_{s}(\boldsymbol{k})\exp \big(\!-i2^{-3/2}s\big(\beta -1\big)t\cos 2\theta _{\boldsymbol{k}}\big) . \\[9pt] \nonumber \end{align}
According to (2.17), given small nonlinearity and visco-diffusion,
$a_{s}$
is a slowly varying function of time. Combining this with small
$\beta -1 ,\ c_{s}$
is also slowly varying, hence so are
$w_{i,s}$
and
$\zeta _{s}$
according to (2.20) and (2.21). Note that
$w_{i,-s}=w_{i,s}^{*} ,\ \zeta _{-s}=\zeta _{s}^{*}$
because
$u_{i,s}=u_{i,-s}^{*}$
and
$\eta _{s}=\eta _{-s}^{*}$
. Note also that small
$\beta -1$
allows us to approximate
$D_{\textit{ss}} , v_{i}^{(s)}$
and
$\eta ^{(s)}$
in (2.17), (2.20) and (2.21) by their values for
$\beta =1$
(values obtained using (I.2.12)–(I.2.15) and (I.2.20) with
$N={\varOmega} =2^{-1/2}$
). From here on,
$D_{\textit{ss}} , v_{i}^{(s)}$
and
$\eta ^{(s)}$
should be understood as having these values. The resulting
$v_{i}^{(s)}$
and
$\eta ^{(s)}$
are given in equations (A13) and (A14). This approximation, which is valid at leading order, has important consequences. In particular, it leads to decoupling of the NP component from the wave one. Note that
$c_{s}(-\boldsymbol{k})=c_{-s}^{*}(\boldsymbol{k})$
follows from (2.22) and
$a_{s}(-\boldsymbol{k})=a_{-s}^{*}(\boldsymbol{k})$
.
Using (2.6), (2.18) and (2.19):
Given that
$w_{i,s}$
and
$\zeta _{s}$
are slowly varying, the sum consists of oscillatory (
$s-s^{\prime}-s^{\prime\prime}\neq 0$
) and non-oscillatory (
$s-s^{\prime}-s^{\prime\prime}=0$
) terms. Since nonlinearity is small, it takes a long time to act and we thus focus on large-
$t$
evolution. Whereas the oscillatory terms in (2.23) only induce small oscillations of
$a_{s}$
when (2.17) is integrated, the non-oscillatory ones can lead to significant cumulative evolution. Thus, we neglect the oscillatory terms, leaving only the non-oscillatory ones,
$s-s^{\prime}-s^{\prime\prime}=0$
; hence (2.17), (2.22) and (2.23) yield
where
$\delta _{0}=1$
and
$\delta _{\sigma }=0$
for
$\sigma \neq 0$
.
When
$s=0$
, it is shown in Appendix A that (2.24) gives
Hence, as in [I] and [II], the NP component evolves independently of the wave one. This is due to vanishing of the NP–wave coupling coefficient
$M_{0s^{\prime},-s^{\prime}}$
(
$s^{\prime}=\pm 1$
) for
$\beta =1$
, where
$M_{{\textit{ss}}^{\prime}s^{\prime\prime}}$
is given by (A11). The resulting decoupling of the NP component from the wave one is one of the main results of this paper, though it is far from obvious a priori. The result indicates that the treatment of the NP component in [I] and [II] continues to apply for small
$\beta -1$
. It should, however, be recalled that
$D_{\textit{ss}} ,\ v_{i}^{(s)}$
and
$\eta ^{(s)}$
here take their values for
$\beta =1$
. Thus, the NP component in the present work is the
$\beta \rightarrow 1$
limit of that of [I] and [II]. Note that decoupling of the NP component has been found to hold numerically by Thomas & Daniel (Reference Thomas and Daniel2021) for particular cases in which
$\beta$
is not close to
$1$
and the wave energy is not too large compared with the NP energy (i.e. the regime covered by [I] and [II]). As far as we are aware, that it also holds for small
$\beta -1$
is a new result.
As regards the wave component, (2.24) yields
Equations (2.20), (2.21), (2.25) and (2.26) form the evolution equations of
$c_{s}$
, which are the basis of this study. Equations (2.20), (2.21) and (2.25) show that the NP component is decoupled from the wave one at leading order, whereas (2.20), (2.21) and (2.26) indicate that the NP component induces coupling between wave modes. These conclusions are similar to those of [II], though no use of wave-turbulence analysis has been made because, as discussed earlier, it does not apply in the present case.
Let
$\varepsilon$
be a small, positive parameter measuring the amplitude of
$u_{i}$
and
$\eta$
. It is then natural to rescale variables to be of order
$1$
by defining
$\hat{u}_{i}=u_{i}/\varepsilon$
and
$\hat{\eta }=\eta /\varepsilon$
. As a result, the rescalings
$\hat{c}_{s}=c_{s}/\varepsilon ,\ \hat{w}_{i,s}=w_{i,s}/\varepsilon ,\ \hat{\zeta }_{s}=\zeta _{s}/\varepsilon$
and
$\hat{t}=\varepsilon t$
are introduced, so (2.20), (2.21), (2.25) and (2.26) become
\begin{align} \hat{w}_{i,s}&=\int\! \hat{c}_{s}(\boldsymbol{k})v_{i}^{(s)}(\boldsymbol{k})\exp \left(i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}\right)\textrm{d}^{3}\boldsymbol{k} , \end{align}
\begin{align} \hat{\zeta }_{s}&=\int\! \hat{c}_{s}(\boldsymbol{k})\eta ^{(s)}(\boldsymbol{k})\exp \left(i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x}\right)\textrm{d}^{3}\boldsymbol{k} , \end{align}
and
where
$\hat{D}_{\textit{NP}}=D_{\textit{NP}}/\varepsilon ,\ \hat{D}_{W}=D_{W}/\varepsilon$
and
$\hat{\beta }=(\beta -1)/\varepsilon$
. Equations (2.27)–(2.30) describe the time evolution of the coefficients
$\hat{c}_{s}$
and are independent of
$\varepsilon$
. Given the definition of
$\hat{t}$
, it is apparent that the time scale for nonlinear evolution is of order
$\varepsilon ^{-1}$
, a non-dimensional scaling which corresponds to the classical large-eddy turnover time (
$L/u' ,\ u'$
being a measure of the turbulent velocities) for evolution of turbulence without rotation or stratification. This is short compared with the time scale,
$O(\varepsilon ^{-2})$
, which is typical of wave-turbulence theory and characterised the evolution of the wave component in [I] and [II]. According to (2.27)–(2.30), the NP component is independent of
$\hat{\beta }$
, whereas, for the wave component,
$\hat{\beta }=(\beta -1)/\varepsilon$
is the parameter which determines closeness or otherwise to
$\beta =1$
. We expect the wave-turbulence formulation of [II] to apply when
$| \hat{\beta }| \rightarrow \infty$
.
It is shown in Appendix B that nonlinearity conserves the total energies of both the NP and wave components, resulting in
\begin{align} \frac{\rm d}{{\rm d}\hat{t}}\int\! \hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}&=-2\int\! \hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} , \end{align}
\begin{align} \frac{\rm d}{{\rm d}\hat{t}}\int\! \hat{e}_{W}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}&=-2\int\! \hat{D}_{W}(\boldsymbol{k})\hat{e}_{W}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} , \end{align}
where
$\hat{e}_{\textit{NP}}=e_{\textit{NP}}/\varepsilon ^{2}$
and
$\hat{e}_{W}=e_{W}/\varepsilon ^{2}$
are scaled energy densities. These equations indicate that the energy of the NP and wave components each decay due to their individual visco-diffusive dissipation. The sum of (2.31) and (2.32) describes the decay of the total energy.
Of course, it should not be concluded from the above energy conservation principles that nonlinearity plays no role in dissipation. Rather, when large-scale dissipation is small (i.e. small
$\hat{D}_{\textit{ss}}(\boldsymbol{k})$
for
$k=O(1)$
, where
$k=| \boldsymbol{k}|$
is the wavenumber), as is usual in turbulent flows, nonlinearity may induce energy transfer from larger to smaller scales, energy which is finally dissipated by linear mechanisms at the smallest scales. This scenario, involving transfer from large to small scales, is often referred to as a forward energy cascade, the expression ‘inverse cascade’ being used to describe energy transfer from small to large scales. In [I], it was found (numerically) that the NP component has only small dissipation, which suggests that, despite nonlinear spectral transfer, there is no forward energy cascade for that component. Since the description of the NP component is here unchanged, the same conclusion applies in the present study. On the other hand, a forward cascade of wave energy was found in [II] and we expect a similar result here.
Note that if one of the two components, NP or wave, is initially absent, it remains so according to (2.27)–(2.30). Thus, pure NP and pure wave flows are possible according to the present theory, as indeed they were in [I] and [II]. In the pure-wave case, the nonlinear term in (2.30) is zero, hence there is no nonlinear transfer and the initial
$\hat{e}_{W}$
simply decays under linear dissipation. Nonlinear effects on wave modes require the NP component, which disallows an energy cascade in the pure-wave case.
For simplicity’s sake, it is henceforth assumed that, in addition to being homogeneous, the initial turbulence has two further statistical symmetries: axisymmetry and symmetry under reflection in any plane perpendicular to the axis of gravity/rotation. Both are preserved by evolution according to the governing equations and so apply at all times. These symmetries lead to
$A_{\textit{ss}}=A_{\textit{ss}}(k,\theta _{\boldsymbol{k}}) ,\ A_{\textit{ss}}(k,\pi -\theta _{\boldsymbol{k}})=A_{\textit{ss}}(k,\theta _{\boldsymbol{k}})$
and
$A_{-1,-1}=A_{11}=e_{W}$
.
The initial (
$\hat{t}=0$
) conditions are chosen according to
\begin{align}&\qquad\quad \hat{A}_{00}(\boldsymbol{k})=k^{2}\exp \big(\!-k^{2}\big) ,\end{align}
\begin{align}& \hat{A}_{11}(\boldsymbol{k})=\hat{A}_{-1,-1}(\boldsymbol{k})=\alpha k^{2}\exp \big(\!-k^{2}\big) , \end{align}
where
$\hat{A}_{\textit{ss}}=A_{\textit{ss}}/\varepsilon ^{2}$
and
$\alpha \geq 0$
is a parameter which determines the relative importance of the wave and NP components. This choice respects the statistical axisymmetry and reflection symmetry discussed above. For simplicity’s sake, the form
$k^{2}\exp (-k^{2})$
used here is the same for both components. This form is isotropic and places the initial energy firmly in the large scales. Time evolution is expected to produce anisotropy and to transfer energy towards smaller scales, leading to inertial and dissipative ranges, and perhaps an energy cascade. Note that the spherically averaged spectra resulting from (2.33) and (2.34) have the form
$k^{4}\exp (-k^{2})$
, a form often used in theoretical studies of turbulence. Integrating (2.33) and (2.34) over
$\boldsymbol{k}$
, (2.14) leads to
\begin{equation} \int\! \hat{e}_{W}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=\frac{3}{2}\alpha \pi ^{3/2} ,\quad \int\! \hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=\frac{3}{4}\pi ^{3/2} \end{equation}
for the total initial scaled energies of the two components. Using
$\hat{e}_{\textit{NP}}=e_{\textit{NP}}/\varepsilon ^{2}$
,
\begin{equation} \varepsilon =\left(\frac{4}{3\pi ^{1/2}}\int\! e_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}\right)^{1/2} \end{equation}
gives
$\varepsilon$
in terms of the initial NP energy. This fixes the precise value of
$\varepsilon$
, up to now an order-of-magnitude parameter.
The full problem, which is integrated numerically, consists of (2.27)–(2.30) with initial conditions for
$\hat{c}_{s}$
which respect from (2.33) and (2.34). The problem, and hence the results of the next section, is independent of
$\varepsilon$
. Given (2.27)–(2.29) and (2.33), the NP component is independent of
$\alpha$
and
$\hat{\beta }$
. Once that component has been determined, (2.27), (2.28) and (2.30) form a linear system for evolution of
$\hat{c}_{s}$
(
$s=\pm 1$
). Equation (2.34) implies that the initial
$\hat{c}_{s}$
scale like
$\alpha ^{1/2}$
and it follows from linearity that
$\hat{c}_{s}(\hat{t})/\alpha ^{1/2}$
is independent of
$\alpha$
, hence the same is true of
$\hat{A}_{\textit{ss}}(\hat{t})/\alpha$
. Thus, the NP calculation has no physical parameters, whereas, when scaled appropriately, the wave one depends only on
$\hat{\beta }$
.
3. Numerical results
As usual in studies of turbulence, dissipation, represented by
$\hat{D}$
in (2.29) and (2.30), is taken to be small for the large scales. In fact, it is only included to ‘mop up’ energy before it can be transferred to higher
$k$
than resolved by the numerics, thus avoiding numerical problems such as instabilities or a pile-up of energy near a large-
$k$
numerical cutoff. As in the DNS of [I], the numerical treatment of (2.27)–(2.30) uses a hyperviscosity in which
$\hat{D}_{\textit{ss}}$
becomes
$dk^{6}$
, with small
$d$
and, for simplicity’s sake,
$d$
has the same value for the wave and NP components. As usual (see e.g. Haugen & Brandenburg Reference Haugen and Brandenburg2004), this is done with small large-scale dissipation in mind and in the hope of extending the expected inertial range to higher wavenumbers, without significantly affecting the large and inertial-range scales.
Given
$\hat{c}_{s}$
at some instant of time, (2.27) and (2.28) yield
$\hat{w}_{i,s}$
and
$\hat{\zeta }_{s}$
in physical space, where the products appearing in (2.29) and (2.30) are taken. Fourier transformation then provides the right-hand sides of (2.29) and (2.30) and hence the time derivative of
$\hat{c}_{s}$
. It will be observed that this procedure is similar to classical DNS and indeed the numerical scheme used to integrate (2.27)–(2.30), and described in Appendix C, is based on that approach. The flow is approximated as periodic, resulting in discretisation in spectral space and truncated Fourier series in place of the Fourier integrals of (2.27) and (2.28). The
$\hat{c}_{s}$
are initialised using (2.33) and (2.34) and random phases. Time stepping of (2.29) and (2.30) gives the flow at later times, where the various measures of energy described below are calculated. Numerical parameters include the hyperviscous coefficient,
$d$
. Others appear in Appendix C, where the values of all the numerical parameters used to obtain the results are given at the end.
In an ideal world the numerical calculations would be used to obtain ensemble averages. However, this would require very many runs, which, given the computational cost, is not attempted. Instead, as usual for DNS, the results given below represent just a single run for each value of
$\hat{\beta }$
. The initial conditions, and hence the subsequent flow, vary randomly from run to run (as described in Appendix C, the initial
$\hat{c}_{s}$
involve random phases). Thus, random fluctuations about the average are to be expected. The numerical scheme also introduces discretisation in spectral space.
In this section, different measures of wave and NP energy are used to present the numerical results. Each of these quantities illustrate different aspects of the results. One consists of the spectra,
$\hat{e}_{W}$
and
$\hat{e}_{\textit{NP}}$
, which represent the energy density in spectral space. Another is the total energies
\begin{equation} \hat{E}_{W}^\textit{tot}=\int\! \hat{e}_{W}(\boldsymbol{k}){\rm d}^{3}\boldsymbol{k} ,\quad \hat{E}_{\textit{NP}}^\textit{tot}=\int\! \hat{e}_{\textit{NP}}(\boldsymbol{k}){\rm d}^{3}\boldsymbol{k} ,\end{equation}
which decay according to (2.31) and (2.32), while a third is the spherically averaged spectra
\begin{align} \hat{E}_{W}(k)&=2\pi \int _{0}^{\pi }k^{2}\hat{e}_{W}\left(k,\theta _{\boldsymbol{k}}\right)\sin \theta _{\boldsymbol{k}}{\rm d}\theta _{\boldsymbol{k}} , \end{align}
\begin{align} \hat{E}_{\textit{NP}}(k)&=2\pi \int _{0}^{\pi }k^{2}\hat{e}_{\textit{NP}}\left(k,\theta _{\boldsymbol{k}}\right)\sin \theta _{\boldsymbol{k}}{\rm d}\theta _{\boldsymbol{k}} , \end{align}
which describe the distribution of energy over
$k=| \boldsymbol{k}|$
and give the total energies of the components when integrated over
$k$
.
Spectral energy fluxes are defined as follows. The wave and NP energies inside the sphere of radius
$k$
and centred on the origin evolve according to
\begin{align} \frac{{\rm d}}{{\rm d}\hat{t}}\int _{\left| \boldsymbol{k}'\right| \lt k}\hat{e}_{W}(\boldsymbol{k}'){\rm d}^{3}\boldsymbol{k}'&=-{\varPhi} _{W}(k)-2\int _{\left| \boldsymbol{k}'\right| \lt k}\hat{D}_{W} (\boldsymbol{k}')\hat{e}_{W} (\boldsymbol{k}'){\rm d}^{3}\boldsymbol{k}' , \end{align}
\begin{align} \frac{\rm d}{{\rm d}\hat{t}}\int _{\left| \boldsymbol{k}'\right| \lt k}\hat{e}_{\textit{NP}} (\boldsymbol{k}'){\rm d}^{3}\boldsymbol{k}'&=-{\varPhi} _{\textit{NP}}(k)-2\int _{\left| \boldsymbol{k}'\right| \lt k}\hat{D}_{\textit{NP}} (\boldsymbol{k}')\hat{e}_{\textit{NP}} (\boldsymbol{k}'){\rm d}^{3}\boldsymbol{k}' , \end{align}
where
${\varPhi} _{W}$
and
${\varPhi} _{\textit{NP}}$
are the spectral fluxes due to nonlinearity and the final terms represent dissipation inside the given spectral sphere. It may be noted that, according to (3.4) and (3.5), the fluxes are zero when
$k=0$
. Taking the limit
$k\rightarrow \infty$
of (3.4) and (3.5) and using (2.31) and (2.32),
${\varPhi} _{W}(\infty )={\varPhi} _{\textit{NP}}(\infty )=0$
represents conservation of energy by nonlinearity.
Beginning with the NP component, there is no dependence on
$\alpha$
or
$\hat{\beta }$
. Figure 1 shows the time evolution of the spherically averaged spectrum. The expected initial transfer of energy from large to small scales in figure 1(a) leads to the establishment of inertial and dissipative ranges. The result is a clear inertial range in which a power law,
$k^{-4}$
, appears to apply. During the later time range of figure 1(b), the exponent of this law gradually evolves from
$-4$
to
$-3$
, the value reported in [I], while the spectral peak moves towards lower
$k$
, corresponding to growth in size of the large scales and suggesting an inverse energy cascade. Energy
$\hat{E}_{\textit{NP}}^\textit{tot}$
was found to decrease by only about
$1\,\%$
between
$\hat{t}=0$
and
$\hat{t}=20$
, a final time much larger than that required for establishment of the dissipative range. This suggests that, as noted earlier and in [I], there is no forward energy cascade for the NP component in the weak-turbulence limit. These results illustrate the distinction between nonlinear energy transfer from large to small scales, which creates such scales, as in figure 1(a), and a forward cascade which delivers significant energy to the dissipative ones. Note that random variations, resulting from use of a single realisation, are apparent in figure 1 for the lower values of
$k$
. Taking the spherical average tends to suppress such fluctuations, the more so as
$k$
increases because more Fourier components contribute to the average. Discretisation is also evident.
Log–log plots of the spherically averaged NP spectrum at times from (a)
$\hat{t}=0$
to
$\hat{t}=3$
in steps of
$0.15$
and (b)
$\hat{t}=3$
to
$\hat{t}=19$
in steps of
$2$
. The dashed line represents the power law
$k^{-3}$
and the finely dashed ones
$k^{-4}$
. The arrows indicate the direction of increasing time.
Plots of the NP energy spectral flux at times from (a)
$\hat{t}=0$
to
$\hat{t}=1$
in steps of
$0.1$
and (b)
$\hat{t}=1$
to
$\hat{t}=6$
in steps of
$1$
. The horizontal axis is logarithmic, while the vertical one is linear. The arrows indicate increasing time.
Figure 2 shows the NP spectral flux,
${\varPhi} _{\textit{NP}}(k)$
, at different times. As noted earlier,
${\varPhi} _{\textit{NP}}(0)={\varPhi} _{\textit{NP}}(\infty )=0$
follows from the definition of
${\varPhi} _{\textit{NP}}$
and conservation of total NP energy by nonlinearity. Furthermore, the initial
${\varPhi} _{\textit{NP}}(k)$
should be zero because the numerical initialisation using random phases (described in Appendix C and in keeping with classical DNS) makes the third-order spectral moments zero. As noted in Appendix B, and reflecting the usual turbulence closure problem, the third-order moments determine nonlinear energy transfers between different wave vectors. When they are zero, there is no transfer and hence no energy flux. In figure 2(a),
${\varPhi} _{\textit{NP}}(k)$
is small for
$\hat{t}=0$
, but not exactly zero. This is because a single realisation, rather than an ensemble average, has been used.
As time increases, figure 2(a) shows the appearance of regions of both positive and negative
${\varPhi} _{\textit{NP}}$
. The former represent transfer from smaller to larger scales, while the latter result in the formation of the inertial and dissipative ranges found earlier. At the later times of figure 2(b), the positive peak drops away and essentially disappears. This leaves only small- to large-scale transfer at large times and again suggests that there is an inverse, rather than a forward, NP energy cascade. That the NP energy cascade is inverse has been recognised for many years in the geophysical community. As noted in [I], (2.25) is equivalent to the three-dimensional quasi-geostrophic equation. The inverse nature of the energy cascade in quasi-geostrophic turbulence was proposed in the seminal work of Charney (Reference Charney1971) and confirmed by later studies (see e.g. Vallgren & Lindborg Reference Vallgren and Lindborg2010). Also note that Charney (Reference Charney1971) suggested the
$k^{-3}$
spectrum found above.
Turning attention to the wave component, the scaled results given below are only dependent on
$\hat{\beta }$
. For each of the
$\hat{\beta }$
considered here, a calculation was performed for
$-\hat{\beta }$
. There are differences between the two, but they are sufficiently small that there is no visible difference for the results given here. The closeness of results for
$\hat{\beta }$
and
$-\hat{\beta }$
suggests that there may be a precise theoretical connection between them, but, if so, we have been unable to find one. In any case, we only give results for
$\hat{\beta }\geq 0$
, because those for
$\hat{\beta }\lt 0$
are essentially identical.
Figure 3 shows the time evolution of the scaled overall wave energy,
$\hat{E}_{W}^\textit{tot}/\alpha$
, for different values of
$\hat{\beta }$
. Following a phase in which the energy is nearly constant, there is a transition to decay. This behaviour is quite different from the near constancy of the NP energy and suggests a forward energy cascade for the wave component.
Log–log plots of the scaled wave energy as a function of time for
$\hat{\beta }=0 ,\ 5 ,\ 10 ,\ 20$
and
$50$
. The arrow indicates increasing
$\hat{\beta }$
.
Figure 4 shows the scaled spherically averaged wave spectrum,
$\hat{E}_{W}/\alpha$
, for different
$\hat{\beta }$
at the same times as figure 1. As for the NP component, there is an initial transfer of energy from large to small scales, followed by the establishment of a dissipative range and slower subsequent evolution. There is very little variation with
$\hat{\beta }$
during the first phase, whereas the dependence on
$\hat{\beta }$
is more pronounced for the second one, shown on the right-hand side of the figures. During the second phase, the wave spectrum decays with time for all
$k$
, whereas there are ranges of both increase and decrease of the NP spectrum in figure 1 for the same range of
$\hat{t}$
. This reflects near energy conservation of the NP component and significant dissipation of the wave one.
Log–log plots of the scaled spherically averaged wave spectrum,
$\hat{E}_{W}/\alpha$
, for (a)
$\hat{\beta }=0$
, (b)
$\hat{\beta }=5$
, (c)
$\hat{\beta }=10$
, (d)
$\hat{\beta }=20$
and (e)
$\hat{\beta }=50$
at times from
$\hat{t}=0$
to
$\hat{t}=3$
in steps of
$0.15$
for (i) and from
$\hat{t}=3$
to
$\hat{t}=19$
in steps of
$2$
for (ii). The dashed lines represent power laws. Arrows representing the direction of time evolution are given.
In figures 4(a), 4(b) and 4(c), at large times,
$\hat{E}_{W}$
becomes roughly independent of
$k$
in the inertial range, which distributes the wave energy across all scales larger than the dissipative ones. This is different from the NP case and also the classical picture of turbulence in the absence of waves, for which, despite the transfer towards smaller ones, the large scales remain the dominant contributors to the energy. Such a strong transfer towards dissipative scales obviously favours dissipation. In figure 4(e), the large-time behaviour is more traditional with the majority of the energy in the large scales, while figure 4(d) represents a case which is somewhat intermediate. Both figures 4(d) and 4(e) show an approximate power law in the inertial range. Such power laws arise at sufficiently large
$\hat{\beta }$
. Once power laws appear, the exponent is increasingly negative as
$\hat{\beta }$
increases. We note that figure 2 of [II] indicates an exponent of about
$-4$
for
$\hat{e}_{W}$
when
$\beta =0.7$
. This corresponds to a spherically averaged exponent of
$-2$
, which may be compared with the value
$-1.5$
of figure 4(e), which is for the largest value of
$\hat{\beta }$
considered here. Given an increasingly negative spectral exponent for increasing
$\hat{\beta }$
, this is in accord with the expectation, noted earlier, that the present results approach those of [II] at large
$| \hat{\beta }|$
.
Figure 5 shows the scaled wave energy flux,
${\varPhi} _{W}(k)/\alpha$
. As for the NP case,
${\varPhi} _{W}(0)=0 ,{} {\varPhi} _{W}(\infty )=0$
and
${\varPhi} _{W}$
is initially zero. It will be observed that the different values of
$\hat{\beta }$
have similar results. Compared to the NP case, a striking feature is that the wave flux is always positive, i.e. there is energy transfer from large to small scales. Furthermore, at large times, the wave energy flux extends through the inertial range before being absorbed in the dissipative range. These results again suggest a forward energy cascade of the wave component.
Plots of the scaled wave energy spectral flux,
${\varPhi} _{W}(k)/\alpha$
, for (a)
$\hat{\beta }=0$
, (b)
$\hat{\beta }=10$
, (c)
$\hat{\beta }=20$
and (d)
$\hat{\beta }=50$
at times from
$\hat{t}=0$
to
$\hat{t}=1$
in steps of
$0.1$
for (i) and from
$\hat{t}=1$
to
$\hat{t}=6$
in steps of
$1$
for (ii). The horizontal axis is logarithmic, while the vertical one is linear. The arrows represent the direction of increasing time.
Finally, we give some results concerning anisotropy. Figure 6 shows contour plots of
$\hat{e}_{\textit{NP}}$
in the
$k_{\bot } , k_{\parallel }$
plane, where
$k_{\bot }=k\sin \theta _{\boldsymbol{k}}$
and
$k_{\parallel }=k\cos \theta _{\boldsymbol{k}}$
are wavenumbers perpendicular and parallel to the rotation/vertical axis. Figure 6(a) is for the initial time and the circular contours reflect the isotropy of the initial conditions. Figure 6(b) is for the largest time,
$\hat{t}=20$
, at which the calculation was performed. The wiggles reflect the fact that these results represent a single realisation, rather than being an ensemble average. Anisotropy favouring the axial direction is apparent. This suggests that NP energy transfer is more effective in that direction than in the transverse direction. These results represent the case
$\beta =1$
, the value appropriate to the present study. Anisotropy of the NP component for other values of
$\beta$
was discussed in § 4.1 of [I]. It was concluded that the NP inertial-range energy distribution is close to isotropic for
$\beta =1.4$
and favours the axial direction for lower values of
$\beta$
(such as the
$\beta =1$
considered here), and the transverse direction for higher values.
Contour plots of
$\hat{e}_{\textit{NP}}$
for (a)
$\hat{t}=0$
and (b)
$\hat{t}=20$
in the
$k_{\bot } , k_{\parallel }$
plane. There are 10 contour values, which are logarithmically spaced from
$10^{-4}$
to
$1$
. The horizontal axis gives
$k_{\bot }$
and the vertical one
$k_{\parallel }$
.
Concerning the wave component, figure 7 shows contour plots of
$\hat{e}_{W}/\alpha$
for
$\hat{t}=20$
and different values of
$\hat{\beta }$
. As for the NP component, anisotropy has a preference for the axial direction. Thus, the total energy density, which is the sum of the NP and wave components, also favours that direction. It further appears from figure 7(b–d) that the energy density decreases near the plane
$k_{\parallel }=0$
. This tendency increases with
$\hat{\beta }$
and reflects the conclusion of [II] that wave dissipation is greatest near the given plane. That this effect strengthens with increasing
$\hat{\beta }$
is to be expected since
$| \hat{\beta }| \rightarrow \infty$
should approach the case considered in [II].
Contour plots of
$\hat{e}_{W}/\alpha$
in the
$k_{\bot } , k_{\parallel }$
plane for
$\hat{t}=20$
and the values (a)
$\hat{\beta }=0$
, (b)
$\hat{\beta }=10$
, (c)
$\hat{\beta }=20$
and (d)
$\hat{\beta }=50$
. There are 10 contour values, which are logarithmically spaced from
$10^{-4}$
to
$1$
. The horizontal axis gives
$k_{\bot }$
and the vertical one
$k_{\parallel }$
.
4. Conclusions
In this article we have considered a particular case, namely small
$\beta -1$
, of weak turbulence (small Rossby or Froude number, the two being equivalent when
$\beta -1$
is small) with rotation and stratification. This case is interesting because it does not permit wave-turbulence analysis, the wave modes being only weakly dispersive. Furthermore, it completes the studies of weak, rotating, stratified turbulence initiated in [I] and carried further in [II].
The flow is expressed as a combination of wave and NP components using modal analysis, as described in § 2, where the weak-turbulence equations for the modal amplitudes were derived. Owing to vanishing of the NP–wave coupling coefficient
$M_{0s^{\prime},-s^{\prime}}(\boldsymbol{k},\boldsymbol{p})$
for
$\beta =1$
and all
$\boldsymbol{k}$
and
$\boldsymbol{p}$
, the NP component is found, as in [I] and [II], to evolve independently of the wave component and according to an amplitude equation which is precisely that of [I] and [II] in the limit
$\beta \rightarrow 1$
. As in [II], the NP component induces coupling between wave modes, though such coupling is no longer constrained by the resonance condition of the earlier study, a condition which resulted from wave-turbulence theory. In the absence of the NP component, nonlinearity has no effect on the wave component, which evolves according to linear dissipation alone.
The parameter
$\hat{\beta }=(\beta -1)/\varepsilon$
appears in the amplitude equations of the wave component, where
$\varepsilon$
is a small parameter (Rossby or Froude number) which measures the weakness of turbulence. Here
$\hat{\beta }$
expresses the departure from
$\beta =1$
, scaled in a manner appropriate to the present work. When
$| \hat{\beta }|$
is sufficiently large, the analysis of [II] should apply, otherwise the present theory is needed. Note that the NP component is independent of
$\hat{\beta }$
.
The time variable appearing in the amplitude evolution equations is
$\hat{t}=\varepsilon t$
. This implies a nonlinear evolution time
$O(\varepsilon ^{-1})$
times the inertial-gravity wave period, which is comparable with the classical large-eddy turnover time for evolution of turbulence without rotation or stratification. This is short compared with the
$O(\varepsilon ^{-2})$
periods typical of wave-turbulence theory and which characterised the wave evolution time in [I] and [II]. Thus, here the wave time scale is
$O(\varepsilon ^{-1})$
, whereas it was
$O(\varepsilon ^{-2})$
before. This suggests that, as
$\beta -1$
decreases from the
$O(1)$
values of [II], wave-component evolution accelerates until the associated time reaches the
$O(\varepsilon ^{-1})$
value of the present work when
$\beta -1=O(\varepsilon )$
.
Energy conservation by nonlinearity was shown to hold for both components individually and hence for their sum, which is the total energy of the flow. Thus, nonlinearity can transfer energy between modes, but not create or destroy it overall. As usual, when dissipation is small at the large scales, which is the usual situation in turbulent flow, nonlinear transfer from large to small scales may lead to a significant energy supply to the smallest ones, which are dissipative. We refer to such a situation as a forward energy cascade. There are, however, examples, such as the NP component considered here, in which nonlinear transfer from small to large scales occurs, but only small energy dissipation results because the energy cascade is towards the large scales (inverse cascade). On the other hand, we found that the dissipation of the wave component was significant, suggesting a forward cascade for that component.
A numerical scheme for solving the mode-amplitude equations was developed. This scheme is similar to that of classical, periodic DNS. The initial conditions ((2.33) and (2.34)) were chosen such that the energy distribution is initially isotropic and concentrated in the large scales. Not surprisingly, it was found that transfer to smaller scales led to the development of inertial and dissipative ranges. This was followed by slower evolution, with rather different consequences for the NP and wave modes.
In the NP case, the result was small energy dissipation and a clear power law for the spherically averaged spectrum in the inertial range. Nonlinear energy transfer from large to small scales plays an important role in the establishment of the inertial range, but the NP component appears to have an inverse, rather than a forward, energy cascade. On the other hand, energy transfer to small scales was much stronger in the wave case, leading to significant dissipation. We conclude that the wave component is subject to a forward energy cascade. At the same time, it should be recalled that this depends on the existence of the NP component, without which coupling of wave modes does not occur. If the NP amplitude is reduced, nonlinear evolution of the wave component takes longer and longer.
In summary, in both [II] and the present work, the NP component acts as a catalyst for nonlinear coupling of wave modes, leading to a forward energy cascade, hence dissipation, of the wave component. Here we use the term catalyst by analogy with its usage in chemistry because the NP component affects the waves, but is not affected in turn. Assuming small dissipation of the large scales and beginning with both components, the wave component is expected to decay, leaving the NP one. Of course, external forcing of the wave component could maintain it.
Numerical results concerning anisotropy were also obtained. The energy distributions of both the NP and wave components show a preference for the axial direction, hence a trend in that direction for the total flow. As
$\hat{\beta }$
was increased from zero, a dip in the intensity of the wave component near the plane perpendicular to the gravity/rotation axis appeared, and intensified at higher
$\hat{\beta }$
. Such behaviour was perhaps to be expected given the conclusion of [II] that wave dissipation is maximal for the given plane.
Finally, it would be interesting to understand why the results for
$-\hat{\beta }$
are nearly the same as those for
$\hat{\beta }$
.
Acknowledgements
Thanks to A. Cadiou for help with the computational aspects of this work and to C. Cambon for helpful discussions and suggestions concerning the manuscript.
Declaration of interests
The author reports no conflict of interest.
Appendix A. Derivation of (2.25)
Equation (2.20) gives
\begin{equation} w_{i,s}^{*}=\int\! c_{s}^{*}(\boldsymbol{p}){v_{i}^{(s)}}^*(\boldsymbol{p})\exp \left(-i\boldsymbol{p}\boldsymbol{.}\boldsymbol{x}\right)\textrm{d}^{3}\boldsymbol{p} .\end{equation}
Thus,
\begin{align} w_{i,s^{\prime}}^{*}w_{j,s^{\prime\prime}}^{*} & =\int\!\! \int c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(\boldsymbol{q}){v_{i}^{(s^{\prime})}}^{*}(\boldsymbol{p}){v_{\!j}^{(s^{\prime\prime})}}^{*}(\boldsymbol{q})\exp (-i(\boldsymbol{p}+\boldsymbol{q})\boldsymbol{.}\boldsymbol{x})\textrm{d}^{3}\boldsymbol{q}\textrm{d}^{3}\boldsymbol{p}\nonumber \\& =\int\!\! \int c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p}){v_{i}^{(s^{\prime})}}^{*}(\boldsymbol{p}){v_{\!j}^{(s^{\prime\prime})}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\exp (i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x})\textrm{d}^{3}\boldsymbol{k}\textrm{d}^{3}\boldsymbol{p}\nonumber \\ & =\int\!\! \int c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p}){v_{i}^{(s^{\prime})}}^{*}(\boldsymbol{p}){v_{\!j}^{(s^{\prime\prime})}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\exp (i\boldsymbol{k}\boldsymbol{.}\boldsymbol{x})\textrm{d}^{3}\boldsymbol{k}. \end{align}
It follows that
Likewise,
Using
$w_{i,-s}=w_{i,s}^{*}$
and
$\zeta _{-s}=\zeta _{s}^{*}$
,
Employing (A3)–(A5), (2.24) gives
\begin{align} & \frac{\partial c_{s}(\boldsymbol{k})}{\partial t}+\big(D_{\textit{ss}}(\boldsymbol{k})+i2^{-3/2}s(\beta -1)\cos 2\theta _{\boldsymbol{k}}\big)c_{s}(\boldsymbol{k})\nonumber\\&\qquad =\sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! N_{{\textit{ss}}^{\prime}{s^{\prime\prime}}}^{*}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}, \end{align}
where
\begin{equation} N_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})=ik_{\!j}{v_{\!j}}^{(s^{\prime\prime})}(-\boldsymbol{k}-\boldsymbol{p})\big(v_{i}^{(s)}(\boldsymbol{k}){v_{i}}^{(s^{\prime})}(\boldsymbol{p})+\eta ^{(s)}(\boldsymbol{k})\eta ^{(s^{\prime})}(\boldsymbol{p})\big) \end{equation}
represents nonlinear coupling between modes. The change of integration variable
$\boldsymbol{q}=-\boldsymbol{k}-\boldsymbol{p}$
gives
\begin{align} & \int\! {N_{{\textit{ss}}^{\prime}s^{\prime\prime}}}^{*}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\nonumber\\&\quad = \int\! {N_{{\textit{ss}}^{\prime}s^{\prime\prime}}}^{*}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{q}\big)c_{s^{\prime}}^{*}\big(-\boldsymbol{k}-\boldsymbol{q}\big)c_{s^{\prime\prime}}^{*}(\boldsymbol{q})\textrm{d}^{3}\boldsymbol{q}\nonumber\\&\quad = \int {N_{{\textit{ss}}^{\prime}s^{\prime\prime}}}^{*}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{p}\big)c_{s^{\prime}}^{*}\big(-\boldsymbol{k}-\boldsymbol{p}\big)c_{s^{\prime\prime}}^{*}(\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}. \end{align}
Thus, permuting the summation indices according to
$s^{\prime}\leftrightarrow s^{\prime\prime}$
,
\begin{align}& \sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! {N_{{\textit{ss}}^{\prime}s^{\prime\prime}}}^{*}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\nonumber\\&\quad = \sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! N_{{\textit{ss}}^{\prime\prime}{s^{\prime}}}^{*}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{p}\big)c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}.\end{align}
It follows that
\begin{align}& \sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! {N_{{\textit{ss}}^{\prime}s^{\prime\prime}}}^{*}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\nonumber\\&\quad = \sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! {M_{{\textit{ss}}^{\prime}s^{\prime\prime}}}^{*}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}, \end{align}
where
\begin{align}&\qquad\quad M_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})=\frac{1}{2}\big(N_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})+N_{{\textit{ss}}^{\prime\prime}s^{\prime}}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{p}\big)\big) . \end{align}
\begin{align} & \frac{\partial c_{s}(\boldsymbol{k})}{\partial t}+\big(D_{\textit{ss}}(\boldsymbol{k})+i2^{-3/2}s (\beta -1)\cos 2\theta _{\boldsymbol{k}}\big)c_{s}(\boldsymbol{k})\nonumber\\&\qquad =\sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! {M_{{\textit{ss}}^{\prime}s^{\prime\prime}}^{*}}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}. \end{align}
Equation (A12) is an evolution equation for the mode amplitudes which is equivalent to (2.20), (2.21) and (2.24).
Adopting Cartesian coordinates with
$x_{3}$
directed along the rotation vector and given that
$v_{i}^{(s)}$
and
$\eta ^{(s)}$
take their values for
$\beta =1$
, (I.2.12)–(I.2.15),
$N={\varOmega} =2^{-1/2}$
and
$\boldsymbol{e}=(0,0,1)$
yield
\begin{equation} \left(\begin{array}{l} v_{1}^{(0)}\\[3pt] v_{2}^{(0)}\\[3pt] v_{3}^{(0)} \end{array}\right)=\frac{i}{k}\left(\begin{array}{l} k_{2}\\ -k_{1}\\ 0 \end{array}\right)\! ,\quad \eta ^{(0)}=\frac{ik_{3}}{k} \end{equation}
and
\begin{equation} \left(\begin{array}{l} v_{1}^{(s)}\\[3pt] v_{2}^{(s)}\\[3pt] v_{3}^{(s)} \end{array}\right)=\frac{1}{2^{1/2}kk_{\bot }}\left(\begin{array}{l} \left(isk_{1}-k_{2}\right)k_{3}\\ \left(isk_{2}+k_{1}\right)k_{3}\\ -isk_{\bot }^{2} \end{array}\right)\! ,\quad \eta ^{(s)}=\frac{k_{\bot }}{2^{1/2}k} \end{equation}
when
$s=\pm 1$
, where
$k=| \boldsymbol{k}|$
and
$k_{\bot }=(k_{1}^{2}+k_{2}^{2})^{1/2}=k\sin \theta _{\boldsymbol{k}}$
. The coefficients,
$M_{{\textit{ss}}^{\prime}s^{\prime\prime}}$
, in (A12) follow from (A7), (A11), (A13) and (A14).
Taking
$s=0$
, (A12) yields
\begin{align} \frac{\partial c_{0}}{\partial t}+D_{\textit{NP}}c_{0}&=\int\! M_{000}^{*}(\boldsymbol{k},\boldsymbol{p})c_{0}^{*}(\boldsymbol{p})c_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\nonumber\\&\quad + \sum _{s^{\prime}=\pm 1}\int\! M_{0s^{\prime},-s^{\prime}}^{*}(\boldsymbol{k},\boldsymbol{p})c_{s^{\prime}}^{*}(\boldsymbol{p})c_{-s^{\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}. \end{align}
\begin{align}& v_{i}^{(0)}(\boldsymbol{k})v_{i}^{(s^{\prime})}(\boldsymbol{p})+\eta ^{(0)}(\boldsymbol{k})\eta ^{(s^{\prime})}(\boldsymbol{p})=k^{-1}Z_{s^{\prime}}(\boldsymbol{k},\boldsymbol{p}) , \end{align}
\begin{align}&\qquad k_{\!j}v_{\!j}^{(-s^{\prime})}(-\boldsymbol{k}-\boldsymbol{p})=s^{\prime}Z_{-s^{\prime}}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{p}\big) , \end{align}
where
$s^{\prime}=\pm 1$
,
\begin{equation} Z_{s^{\prime}}(\boldsymbol{k},\boldsymbol{p})=\frac{i}{2^{1/2}pp_{\bot }}\big(k_{2}p_{3}\big({\textit{is}}^{\prime}p_{1}-p_{2}\big)-k_{1}p_{3}\big({\textit{is}}^{\prime}p_{2}+p_{1}\big)+k_{3}p_{\bot }^{2}\big) ,\end{equation}
$p=| \boldsymbol{p}|$
and
$p_{\bot }=(p_{1}^{2}+p_{2}^{2})^{1/2}$
. Equations (A7), (A16) and (A17) give
\begin{equation} N_{0s^{\prime},-s^{\prime}}(\boldsymbol{k},\boldsymbol{p})={\textit{is}}^{\prime}k^{-1}Z_{s^{\prime}}(\boldsymbol{k},\boldsymbol{p})Z_{-s^{\prime}}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{p}\big) ,\end{equation}
hence
$M_{0s^{\prime},-s^{\prime}}(\boldsymbol{k},\boldsymbol{p})=0$
from (A11). Thus, the second term on the right-hand side of (A15) is zero, while (A11) implies
\begin{align} & \int\! M_{000}^{*}(\boldsymbol{k},\boldsymbol{p})c_{0}^{*}(\boldsymbol{p})c_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\nonumber\\&\quad = \frac{1}{2}\left(\int\! N_{000}^{*}(\boldsymbol{k},\boldsymbol{p})c_{0}^{*}(\boldsymbol{p})c_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\right.\nonumber\\&\qquad + \left.\int\! N_{000}^{*}\big(\boldsymbol{k},-\boldsymbol{k}-\boldsymbol{p}\big)c_{0}^{*}(\boldsymbol{p})c_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\right)\!. \end{align}
\begin{align} & \int\! M_{000}^{*}(\boldsymbol{k},\boldsymbol{p})c_{0}^{*}(\boldsymbol{p})c_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\nonumber\\&\quad = \int\! N_{000}^{*}(\boldsymbol{k},\boldsymbol{p})c_{0}^{*}(\boldsymbol{p})c_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}, \end{align}
hence
Appendix B. Energy conservation by nonlinearity
The scaled version of (A6) is
\begin{align} &\frac{\partial \hat{c}_{s}(\boldsymbol{k})}{\partial \hat{t}}+\big(\hat{D}_{\textit{ss}}(\boldsymbol{k})+i2^{-3/2}s\hat{\beta }\cos 2\theta _{\boldsymbol{k}}\big)\hat{c}_{s}(\boldsymbol{k})\nonumber \\&\quad= \sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! N_{{\textit{ss}}^{\prime}s^{\prime\prime}}^{*}(\boldsymbol{k},\boldsymbol{p}){\hat{c}}_{s^{\prime}}^{*}(\boldsymbol{p}){\hat{c}}_{s^{\prime\prime}}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}. \end{align}
It follows from (2.10) and (2.22) that
\begin{equation} \overline{\hat{c}_{s}^{*}(\boldsymbol{k})\hat{c}_{s}\left(\boldsymbol{k}'\right)}=\hat{A}_{\textit{ss}}(\boldsymbol{k})\delta \left(\boldsymbol{k}-\boldsymbol{k}'\right)\! .\end{equation}
Given statistical homogeneity, the third-order spectral moments have the form
\begin{equation} \overline{\hat{c}_{s}(\boldsymbol{k})\hat{c}_{s^{\prime}}(\boldsymbol{p})\hat{c}_{s^{\prime\prime}} (\boldsymbol{k}'')}={\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})\delta \big(\boldsymbol{k}+\boldsymbol{p}+\boldsymbol{k}''\big) .\end{equation}
Multiplying the complex conjugate of (B1) by
$\hat{c}_{s}(\boldsymbol{k}')$
and taking the ensemble average:
\begin{align} \overline{\frac{\partial \hat{c}_{s}^{*}(\boldsymbol{k})}{\partial \hat{t}}\hat{c}_{s}\big(\boldsymbol{k}'\big)}&=\delta \big(\boldsymbol{k}-\boldsymbol{k}'\big)\big(\tau _{s}(\boldsymbol{k})-\big(\hat{D}_{\textit{ss}}(\boldsymbol{k})-i2^{-3/2}s\hat{\beta }\cos 2\theta _{\boldsymbol{k}}\big)\hat{A}_{\textit{ss}}(\boldsymbol{k})\big) , \end{align}
where
\begin{align} \tau _{s}(\boldsymbol{k})&=\sum _{s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\! N_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p}){\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p} \end{align}
and we have used (B2), (B3) and the fact that
$\hat{D}_{\textit{ss}}$
is real. Given the Dirac function in (B4), we have also replaced
${\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k}',\boldsymbol{p})$
by
${\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})$
in (B5). Taking the conjugate of (B4) and permuting
$\boldsymbol{k}\leftrightarrow \boldsymbol{k}'$
:
\begin{align} \overline{\hat{c}_{s}^{*}(\boldsymbol{k})\frac{\partial \hat{c}_{s}\big(\boldsymbol{k}'\big)}{\partial \hat{t}}} & =\delta \big(\boldsymbol{k}-\boldsymbol{k}'\big)\big(\tau _{s}^{*}\big(\boldsymbol{k}'\big)-\big(\hat{D}_{\textit{ss}}\big(\boldsymbol{k}'\big)+i2^{-3/2}s\hat{\beta }\cos 2\theta _{\boldsymbol{k}'}\big)\hat{A}_{\textit{ss}}\big(\boldsymbol{k}'\big)\big)\nonumber \\ & =\delta \big(\boldsymbol{k}-\boldsymbol{k}'\big)\big(\tau _{s}^{*}(\boldsymbol{k})-\big(\hat{D}_{\textit{ss}}(\boldsymbol{k})+i2^{-3/2}s\hat{\beta }\cos 2\theta _{\boldsymbol{k}}\big)\hat{A}_{\textit{ss}}(\boldsymbol{k})\big), \end{align}
where we have used the fact that
$\hat{D}_{\textit{ss}}$
and
$\hat{A}_{\textit{ss}}$
are real. Taking the sum of (B4) and (B6), the time derivative of (B2) implies
\begin{equation} \frac{\partial \hat{A}_{\textit{ss}}}{\partial \hat{t}}=2\big({\textit{Re}} \big\{\tau _{s}(\boldsymbol{k})\big\}-\hat{D}_{\textit{ss}}(\boldsymbol{k})\hat{A}_{\textit{ss}}(\boldsymbol{k})\big) .\end{equation}
This result shows that the second-order spectral moments
$\hat{A}_{\textit{ss}}$
evolve due to nonlinearity and dissipation. The nonlinear term involves the third-order moments via (B5), which reflects the usual closure problem of turbulence.
The energy density in spectral space is given by (2.13). Hence,
\begin{equation} \frac{{\rm d}}{{\rm d}\hat{t}}\int\! \hat{e}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}={\textit{Re}} \left\{\sum _{s=0,\pm 1}\int\! \tau _{s}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}\right\}-\sum _{s=0,\pm 1}\int\! \hat{D}_{\textit{ss}}(\boldsymbol{k})\hat{A}_{\textit{ss}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} \end{equation}
for the evolution of the total energy. Using (B5),
\begin{equation} \sum _{s=0,\pm 1}\int\! \tau _{s}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=\sum _{s,s^{\prime},s^{\prime\prime}=0,\pm 1}\delta _{s+s^{\prime}+s^{\prime\prime}}\int\!\! \int\! N_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p}){\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\textrm{d}^{3}\boldsymbol{k} .\end{equation}
The identity
$k_{\!j}v_{\!j}^{(s)}(\boldsymbol{k})=0$
follows from the definition, (I.2.12) and (I.2.13), of
$v_{\!j}^{(s)}(\boldsymbol{k})$
. Thus,
$(k_{\!j}+p_{\!j})v_{\!j}^{(s^{\prime\prime})}(-\boldsymbol{k}-\boldsymbol{p})=0$
, hence
$N_{s^{\prime}{\textit{ss}}^{\prime\prime}}(\boldsymbol{p},\boldsymbol{k})=-N_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})$
from (A7). On the other hand,
${\Theta} _{s^{\prime}{\textit{ss}}^{\prime\prime}}(\boldsymbol{p},\boldsymbol{k})={\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})$
according to (B3). It follows that the integrand in (B9) changes sign under the permutation
$s\leftrightarrow s^{\prime} , \boldsymbol{k}\leftrightarrow \boldsymbol{p}$
. Thus,
\begin{equation} \sum _{s=0,\pm 1}\int\! \tau _{s}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=0 ,\end{equation}
hence nonlinearity conserves the total energy, which evolves according to
\begin{equation} \frac{\rm d}{{\rm d}\hat{t}}\int\! \hat{e}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=-\sum _{s=0,\pm 1}\int\! \hat{D}_{\textit{ss}}(\boldsymbol{k})\hat{A}_{\textit{ss}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} .\end{equation}
Next consider the NP component of the energy. As was shown in Appendix A, the second term on the right-hand side of (A15) is zero, while the first term can be expressed using (A21). Thus,
\begin{equation} \frac{\partial \hat{c}_{0}}{\partial \hat{t}}+\hat{D}_{\textit{NP}}\hat{c}_{0}=\int\! N_{000}^{*}(\boldsymbol{k},\boldsymbol{p})\hat{c}_{0}^{*}(\boldsymbol{p})\hat{c}_{0}^{*}(-\boldsymbol{k}-\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p} .\end{equation}
Multiplying the complex conjugate of (B12) by
$\hat{c}_{0}(\boldsymbol{k}')$
and taking the ensemble average:
\begin{align} \overline{\frac{\partial \hat{c}_{0}^{*}(\boldsymbol{k})}{\partial \hat{t}}\hat{c}_{0}\big(\boldsymbol{k}'\big)}&=\delta \big(\boldsymbol{k}-\boldsymbol{k}'\big)\big(\tau (\boldsymbol{k})-\hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{A}_{00}(\boldsymbol{k})\big) , \end{align}
where
\begin{align} \tau (\boldsymbol{k})&=\int\! N_{000}(\boldsymbol{k},\boldsymbol{p}){\Theta} _{000}(\boldsymbol{k},\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p} \end{align}
and we have used (B2), (B3) and the fact that
$\hat{D}_{\textit{NP}}$
is real. Given the Dirac function in (B13), we have also replaced
${\Theta} _{000}(\boldsymbol{k}',\boldsymbol{p})$
by
${\Theta} _{000}(\boldsymbol{k},\boldsymbol{p})$
in (B14). Taking the conjugate of (B13) and permuting
$\boldsymbol{k}\leftrightarrow \boldsymbol{k}'$
:
\begin{align} \overline{\hat{c}_{0}^{*}(\boldsymbol{k})\frac{\partial \hat{c}_{0}(\boldsymbol{k}')}{\partial \hat{t}}} & =\delta \big(\boldsymbol{k}-\boldsymbol{k}'\big)\big(\tau ^{*}\big(\boldsymbol{k}'\big)-\hat{D}_{\textit{NP}}(\boldsymbol{k}')\hat{A}_{00}(\boldsymbol{k}')\big)\nonumber \\ & =\delta \big(\boldsymbol{k}-\boldsymbol{k}'\big)\big(\tau ^{*}(\boldsymbol{k})-\hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{A}_{00}(\boldsymbol{k})\big), \end{align}
where we have used the fact that
$\hat{D}_{\textit{NP}}$
and
$\hat{A}_{00}$
are real. Taking the sum of (B13) and (B15), (B2) implies
\begin{equation} \frac{\partial \hat{A}_{00}}{\partial \hat{t}}=2\left({\textit{Re}} \left\{\tau (\boldsymbol{k})\right\}-\hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{A}_{00}(\boldsymbol{k})\right)\! .\end{equation}
Given the second equation of (2.14),
\begin{equation} \frac{\rm d}{{\rm d}\hat{t}}\int\! \hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}={\textit{Re}} \left\{\int\! \tau (\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}\right\}-2\int\! \hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} \end{equation}
for the evolution of the NP energy. Using (B14),
\begin{equation} \int\! \tau (\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=\int\!\! \int N_{000}(\boldsymbol{k},\boldsymbol{p}){\Theta} _{000}(\boldsymbol{k},\boldsymbol{p})\textrm{d}^{3}\boldsymbol{p}\textrm{d}^{3}\boldsymbol{k} .\end{equation}
It was shown above that
$N_{s^{\prime}{\textit{ss}}^{\prime\prime}}(\boldsymbol{p},\boldsymbol{k})=-N_{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})$
and
${\Theta} _{s^{\prime}{\textit{ss}}^{\prime\prime}}(\boldsymbol{p},\boldsymbol{k})={\Theta} _{{\textit{ss}}^{\prime}s^{\prime\prime}}(\boldsymbol{k},\boldsymbol{p})$
, hence the integrand in (B18) changes sign under
$\boldsymbol{k}\leftrightarrow \boldsymbol{p}$
. Thus, (B18) is zero and (B17) gives
\begin{equation} \frac{{\rm d}}{{\rm d}\hat{t}}\int\! \hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}=-2\int\! \hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} ,\end{equation}
i.e. (2.31).
Since
$\hat{D}_{00}=\hat{D}_{\textit{NP}}$
and
$\hat{D}_{11}=\hat{D}_{-1,-1}=\hat{D}_{W}$
, (2.14) gives
\begin{equation} \sum _{s=0,\pm 1}\int\! \hat{D}_{\textit{ss}}(\boldsymbol{k})\hat{A}_{\textit{ss}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}= 2\int\! \hat{D}_{\textit{NP}}(\boldsymbol{k})\hat{e}_{\textit{NP}}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}+2\int\! \hat{D}_{W}(\boldsymbol{k})\hat{e}_{W}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k}. \end{equation}
Using (B11), (B19), (B20) and
$\hat{e}=\hat{e}_{\textit{NP}}+\hat{e}_{W}$
yields (2.32).
Appendix C. The numerical scheme
In what follows, we use Cartesian coordinates,
$x_{1} , x_{2}, x_{3}$
, where
$x_{3}$
is directed along the rotation vector and the flow is approximated as periodic with respect to
$x_{1} ,\ x_{2}$
and
$x_{3}$
. Thus, (2.27) and (2.28) are replaced by the Fourier series
in which the sums are over
$\boldsymbol{k}$
having components
$k_{1}=2\pi l_{1}/L_{1} ,\ k_{2}=2\pi l_{2}/L_{2} ,\ k_{3}=2\pi l_{3}/L_{3}$
, where
$l_{i}$
take integer values and
$L_{i}$
are the spatial periods. As usual in DNS, the spatial periods,
$L_{i}$
, need to be sufficiently large that periodicity is unimportant. Another way of putting this is that the spacing
$2\pi /L_{i}$
in spectral space should be sufficiently small. The numerical equivalents of (2.29) and (2.30) are
with suitable interpretations of 
and 
in terms of discrete Fourier transforms.
The quantities
$v_{i}^{(s)}$
and
$\eta ^{(s)}$
appearing in (C1)–(C4) are given by (A13) and (A14). There are, however, problems when
$k_{\bot }=0$
, namely zero-by-zero divisions. We first consider the case
$k\neq 0$
, for which there is no such problem for
$v_{i}^{(0)}$
and
$\eta ^{(s)}$
. A difficulty arises for
$s=\pm 1$
because (A14) implies
$v_{i}^{(s)}(\boldsymbol{k})$
has different values in the limit
$k_{\bot }\rightarrow 0$
, depending on the direction from which the limit is approached. This ambiguity for
$k_{\bot }=0$
is resolved by choosing a particular limiting direction. Taking the limit
$k_{1}\searrow 0$
, with
$k_{2}=0$
and
$k_{3}\gt 0$
constant, (A14) gives
\begin{equation} \left(\begin{array}{l} v_{1}^{(s)}\\[3pt] v_{2}^{(s)}\\[3pt] v_{3}^{(s)} \end{array}\right)=2^{-1/2}\left(\begin{array}{l} is\\ 1\\ 0 \end{array}\right)\! ,\quad \eta ^{(s)}=0\quad s=\pm 1 .\end{equation}
To respect the identities
$v_{i}^{(s)}(-\boldsymbol{k})={v_{i}^{(-s)}}^*(\boldsymbol{k})$
and
$\eta ^{(s)}(-\boldsymbol{k})={\eta ^{(-s)}}^*(\boldsymbol{k})$
, (C5) is extended to
$k_{3}\lt 0$
. This provides the values of
$v_{i}^{(s)}(\boldsymbol{k})$
and
$\eta ^{(s)}(\boldsymbol{k})$
used in the numerical implementation for
$k_{\bot }=0 ,\ k_{3}\neq 0$
and
$s=\pm 1$
. Given zero-by-zero division in both (A13) and (A14), a problem remains when
$k=0$
. This is resolved by taking 
. Thus, the
$\boldsymbol{k}=0$
terms in (C1) and (C2) are zero and (C3) and (C4) are not used for zero
$\boldsymbol{k}$
. Note that the right-hand sides of (C3) and (C4) are zero when
$\boldsymbol{k}=0$
, hence the choice 
is consistent with time evolution according to those equations.
In addition to supposing periodicity, the sums in (C1) and (C2) are truncated to
$| l_{i}| \leq K_{i}$
, where
$K_{i}$
are positive integers, which must be sufficiently large that the smallest dynamically significant length scales of the flow are resolved. Evaluating (C1) and (C2) at the discrete points
$x_{1}=p_{1}L_{1}/N_{1} ,\ x_{2}=p_{2}L_{2}/N_{2} ,\ x_{3}=p_{3}L_{3}/N_{3}$
, where
$N_{i}\gt 2K_{i}$
and
$0\leq p_{i}\lt N_{i}$
are integers, the result can be expressed as
\begin{align} \hat{w}_{i,s}&=\sum _{q_{1}=0}^{N_{1}-1}\sum _{q_{2}=0}^{N_{2}-1}\sum _{q_{3}=0}^{N_{3}-1}\xi _{i,s}\left(q_{1},q_{2},q_{3}\right)\exp \left(2\pi i\left(\frac{q_{1}p_{1}}{N_{1}}+\frac{q_{2}p_{2}}{N_{2}}+\frac{q_{3}p_{3}}{N_{3}}\right)\right)\! , \end{align}
\begin{align} \hat{\zeta }_{s}&=\sum _{q_{1}=0}^{N_{1}-1}\sum _{q_{2}=0}^{N_{2}-1}\sum _{q_{3}=0}^{N_{3}-1}\xi _{s}\left(q_{1},q_{2},q_{3}\right)\exp \left(2\pi i\left(\frac{q_{1}p_{1}}{N_{1}}+\frac{q_{2}p_{2}}{N_{2}}+\frac{q_{3}p_{3}}{N_{3}}\right)\right)\! , \end{align}
in which
$\xi _{i,s}=\xi _{s}=0$
unless
$q_{i}\leq K_{i}$
or
$q_{i}\geq N_{i}-K_{i}$
. Otherwise,
where
$\boldsymbol{k}=2\pi (l_{1}/L_{1},l_{2}/L_{2},l_{3}/L_{3}) ,\ l_{i}=q_{i}$
when
$q_{i}\leq K_{i}$
and
$l_{i}=q_{i}-N_{i}$
for
$q_{i}\geq N_{i}-K_{i}$
. The right-hand sides of equations (C6) and (C7) are discrete Fourier transforms which are evaluated using fast Fourier transforms (FFT).
Having determined
$\hat{w}_{i,s}$
and
$\hat{\zeta }_{s}$
in physical space, the products
$\hat{w}_{i,s^{\prime}}\hat{w}_{j,s^{\prime\prime}}$
and
$\hat{\zeta }_{s^{\prime}}\hat{w}_{j,s^{\prime\prime}}$
can be calculated. Performing the discrete Fourier transforms
\begin{align}& \chi _{ijs^{\prime}s^{\prime\prime}}\left(q_{1},q_{2},q_{3}\right)= \frac{1}{N_{1}N_{2}N_{3}}\sum _{p_{1}=0}^{N_{1}-1}\sum _{p_{2}=0}^{N_{2}-1}\sum _{p_{3}=0}^{N_{3}-1}\hat{w}_{i,s^{\prime}}\hat{w}_{j,s^{\prime\prime}}\left(p_{1},p_{2},p_{3}\right)\nonumber\\&\qquad\qquad\qquad\qquad\times \exp \left(-2\pi i\left(\frac{q_{1}p_{1}}{N_{1}}+\frac{q_{2}p_{2}}{N_{2}}+\frac{q_{3}p_{3}}{N_{3}}\right)\right)\!, \end{align}
\begin{align}& \lambda _{js^{\prime}s^{\prime\prime}}\left(q_{1},q_{2},q_{3}\right)= \frac{1}{N_{1}N_{2}N_{3}}\sum _{p_{1}=0}^{N_{1}-1}\sum _{p_{2}=0}^{N_{2}-1}\sum _{p_{3}=0}^{N_{3}-1}\hat{\zeta }_{s^{\prime}}\hat{w}_{j,s^{\prime\prime}}\left(p_{1},p_{2},p_{3}\right)\nonumber\\&\qquad\qquad\qquad\qquad\times \exp \left(-2\pi i\left(\frac{q_{1}p_{1}}{N_{1}}+\frac{q_{2}p_{2}}{N_{2}}+\frac{q_{3}p_{3}}{N_{3}}\right)\right)\!, \end{align}
the nonlinear terms in (C3) and (C4) follow from
where
$q_{i}=l_{i}$
for
$0\leq l_{i}\leq K_{i}$
and
$q_{i}=N_{i}+l_{i}$
for
$-K_{i}\leq l_{i}\lt 0$
. As for (C6) and (C7), FFT is used to evaluate (C10) and (C11).
Writing the numerical mode-amplitude equations, (C3) and (C4), as
a time-stepping scheme is introduced. Time is discretised to the values
$\hat{t}_{n}=n{\varDelta}$
, where
${\varDelta}$
is the time step, and (C13) gives
Given 
is determined using FFT as described above. As a first approximation,
$\gamma _{s}(\hat{t}^{\prime})$
is taken as constant and equal to
$\gamma _{s}(\hat{t}_{n})$
, hence
as an approximate value of 
. To increase the precision to second order, let
$\gamma _{s}^{\dagger }$
be the value obtained using 
, which approximates the value of 
at time
$\hat{t}_{n}+{\varDelta} /2$
. Thus,
completes the time step. The advantage of the above scheme, compared with more traditional ones for integrating differential equations, is that it only assumes that the time scale for significant changes in 
is long compared with
${\varDelta}$
. As
$k$
increases, so does
$\hat{D}_{\textit{ss}}(\boldsymbol{k})$
, and, when
$\hat{D}_{\textit{ss}}(\boldsymbol{k}){\varDelta}$
becomes of
$O(1)$
, dissipative effects are significant within a time step. At this point, the precision of traditional schemes is degraded, whereas it is maintained for the above scheme provided 
does not change too much over a single step.
According to (C1) and (C2) and the scaled versions of (2.7) and (2.19):
Thus,
Given homogeneity, 
when
$\boldsymbol{k}'\neq \boldsymbol{k}$
, hence
Equation (I.2.13) implies that the vectors
$\boldsymbol{e}^{(1)}$
and
$\boldsymbol{e}^{(2)}$
are orthonormal, hence (I.2.12) yields
\begin{equation} {v_{i}^{(s)}}^*(\boldsymbol{k})v_{i}^{(s^{\prime})}(\boldsymbol{k})=\sum _{l=1}^{2}{u_{s}^{(l)}}^*u_{s^{\prime}}^{(l)} .\end{equation}
Using this result in (C20), (I.2.16) gives
Prior to numerical approximation, the scaled version of (2.12) implies
\begin{equation} \frac{1}{2}\big(\overline{\hat{u}_{i}\hat{u}_{i}}+\overline{\hat{\eta }^{2}}\big)=\int\! \hat{e}(\boldsymbol{k})\textrm{d}^{3}\boldsymbol{k} .\end{equation}
Equation (C22) is a numerical approximation of (C23), the integral being represented by a sum. Since the density of discrete
$\boldsymbol{k}$
in spectral space is
$L_{1}L_{2}L_{3}/8\pi ^{3}$
,
approximates the spectral energy density. Comparing the scaled version of (2.13) and (C24),
provides a numerical approximation of
$\hat{A}_{\textit{ss}}=A_{\textit{ss}}/\varepsilon ^{2}$
.
The initial 
must satisfy various constraints. Firstly, 
because the velocity and buoyancy variable of the flow are real. Secondly, they should be of zero mean and, as noted above, homogeneity requires 
when
$\boldsymbol{k}'\neq \boldsymbol{k}$
. Finally,
from (C25). Supposing the initial
$\hat{A}_{\textit{ss}}(\boldsymbol{k})$
are given and are such that
$\hat{A}_{\textit{ss}}(-\boldsymbol{k})=\hat{A}_{-s,-s}(\boldsymbol{k})$
, these conditions are satisfied by taking
where the phases,
$\phi _{s}(\boldsymbol{k})$
, are random variables, uniformly distributed between
$-\pi$
and
$\pi$
and independent apart from the constraint
$\phi _{s}(-\boldsymbol{k})=-\phi _{-s}(\boldsymbol{k})$
. Equation (C27) also implies statistical homogeneity of the initial flow defined by (C1), (C2), another important requirement.
Having described the methods used to initialise and integrate (C3) and (C4), we now want to extract information from the results. The problem with using (C25) to determine
$\hat{A}_{\textit{ss}}(\boldsymbol{k})$
is that it contains an ensemble average. To calculate such an average numerically would require very many runs, which, given the computational cost, is not attempted. We instead extract information from a single run. In so doing, we assume that, prior to discretisation, the flow is statistically axisymmetric. In particular,
$\hat{A}_{\textit{ss}}(\boldsymbol{k})$
is axisymmetric and thus has the form
$\hat{A}_{\textit{ss}}=\hat{A}_{\textit{ss}}(k_{\bot },k_{\parallel })$
, where
$k_{\bot }=k\sin \theta _{\boldsymbol{k}}$
and
$k_{\parallel }=k\cos \theta _{\boldsymbol{k}}$
are wavenumbers perpendicular to and parallel to the rotation vector. The aim here is to obtain results which are representative of
$\hat{A}_{\textit{ss}}(k_{\bot },k_{\parallel })$
without using ensemble averaging.
For simplicity’s sake, from here on we assume
$L_{1}=L_{2}=L_{3}=L$
and
$K_{1}=K_{2}=K_{3}=K$
. Let
\begin{equation} k_{m}=\frac{\pi }{L}\left(2m-1\right) \end{equation}
and define axisymmetric regions
$D(m_{\bot },m_{\parallel })$
in spectral space by
$k_{{m_{\bot }}}\leq k_{\bot }\lt k_{{m_{\bot }}+1}$
and
$k_{{m_{\parallel }}}\leq k_{\parallel }\lt k_{{m_{\parallel }}+1}$
, where
$0\leq m_{\bot }\leq K$
and
$-K\leq m_{\parallel }\leq K$
. Let
where the sum is over discrete wave vectors inside
$D(m_{\bot },m_{\parallel })$
and
$N_{{m_{\bot }}{m_{\parallel }}}$
is the number of such wave vectors (which is non-zero given the construction of
$D(m_{\bot },m_{\parallel })$
). Taking the average of (C29) and using (C25):
\begin{equation} \overline{B_{s}^{m_{\bot }m_{\parallel }}}=\frac{1}{N_{{m_{\bot }}{m_{\parallel }}}}\sum _{\boldsymbol{k}\in D\left(m_{\bot },m_{\parallel }\right)}\hat{A}_{\textit{ss}}(\boldsymbol{k}) .\end{equation}
If the region is small in the
$k_{\bot }$
–
$k_{\parallel }$
plane (large enough
$L$
), variations of
$\hat{A}_{\textit{ss}}(\boldsymbol{k})$
across
$D(m_{\bot },m_{\parallel })$
can be neglected and
$\overline{B_{s}^{m_{\bot }m_{\parallel }}}$
approximates
$\hat{A}_{\textit{ss}}(\boldsymbol{k})$
throughout the given region, in particular at the point
\begin{align} k_{\bot }&=\frac{2\pi }{L}m_{\bot }\quad 0\leq m_{\bot }\leq K , \end{align}
\begin{align} k_{\parallel }&=\frac{2\pi }{L}m_{\parallel }\quad -K\leq m_{\parallel }\leq K . \\[9pt] \nonumber \end{align}
The problem is that, as noted above, it is not feasible to determine the ensemble average using numerical calculations. Thus, we take
$B_{s}^{m_{\bot }m_{\parallel }}$
, rather than its average, as representative of
$\hat{A}_{\textit{ss}}(k_{\bot },k_{\parallel })$
at (C31), (C32).
The spherically averaged energy spectrum for
$k$
given by
\begin{equation} k=\frac{2\pi }{L}m\qquad 0\leq m\leq K \end{equation}
is determined using
where
$S(m)$
denotes the region
$k_{m}\leq k\lt k_{m+1}$
in spectral space. The wave and NP parts of this spectrum are obtained by corresponding restrictions of the sum over
$s$
. The total energy follows from (C22) (without averaging) and wave/NP decomposition results from restricting the sum over
$s$
.
Finally, the spectral fluxes are obtained as follows. The unaveraged wave energy inside the sphere of radius
$k$
is
whose time derivative is
according to (C13). Comparison with (3.4) gives
Similarly,
for the NP radial flux. Equation (C37) and (C38) are applied for the wavenumbers (C33).
The numerical parameters used to obtain the results described in § 3 were
$N_{i}=512 ,\ K_{i}=255 ,\ L_{i}=60 ,\ {\varDelta} =5\times 10^{-3}$
and
$d=1.5\times 10^{-7}$
, where, as noted in the main text, the hyperviscous damping factor is
$\hat{D}_{\textit{ss}}(\boldsymbol{k})=dk^{6}$
.
Open access
