[image] [image] Ressource protégée - Archive ouverte HAL Accéder directement au contenu

*[image] [image]*

• [Documentation](#)

• FR
Français (FR)

• Anglais (EN)
Se connecter

• [HAL science ouverte](#)

• Recherche
Loading...

## Recherche avancée

Information de documents

Titres

○ [Titres](#)

○ [Sous-titre](#)

○ Titre de l'ouvrage

○ Titre du volume (Série)

○ Champ de recherche par défaut (multicritères) + texte intégral des PDF

○ [Résumé](#)

○ Texte intégral indexé des documents PDF

○ [Mots-clés](#)

○ Type de document

○ Sous-type de document

○ Tous les identifiants du document

○ Identifiant HAL du dépôt

○ Langue du document (texte)

○ Pays (Texte)

○ Ville

○ À paraître (true ou false)

Ajouter

Auteur
# Evolution of weak, homogeneous turbulence with rotation and stratification

J.F. Scott, C. Cambon

<figure><img src="image_1.png" /></figure>

## To cite this version:

J.F. Scott, C. Cambon. Evolution of weak, homogeneous turbulence with rotation and stratification. Journal of Fluid Mechanics, 2024, 979, pp.A17. 10.1017/jfm.2023.1046. hal-04389842

HAL Id: hal-04389842

https://hal.science/hal-04389842v1

Submitted on 12 Jan 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
## Auteur (multicritères)

○ Auteur (multicritères)

○ Auteur : Nom complet

○ Auteur : Nom de famille

○ Auteur : Prénom

○ Auteur : Complément de nom, deuxième prénom

○ Auteur : Organisme payeur

○ Auteur : IdHal (chaîne de caractères)

○ Auteur : Fonction

○ Auteur : personID (entier)

○ Auteur : Identifiant de l'organisme payeur

○ Auteur : Identifiant de la structure

○ Directeur de thèse

○ Éditeur

○ Éditeur scientifique

○ Éditeur de la série

## Ajouter

### Structure

### Structure (multicritères)

○ Structure (multicritères)

○ Structure : Acronyme

○ Structure : Nom

○ Structure : Code

○ Structure : Pays

○ Structure : Type

○ Structure : État dans le référentiel

○ Structure : Identifiant HAL de la structure

○ Structure/équipe de recherche : Pays

○ Structure/regroupement d'équipes : Pays

○ Structure/laboratoire : Pays

○ Structure/regroupement de laboratoires : Pays

○ Structure/institution : Pays

○ Structure/regroupement d'institutions : Pays

## Ajouter

## Autres
# Evolution of weak, homogeneous turbulence with rotation and stratification

J. F. Scott, C. Cambon

Laboratoire de Mécanique des Fluides et d'Acoustique, Université de Lyon,
Ecole Centrale de Lyon, 36 avenue Guy de Collongue, 69134 Ecully CEDEX, France

## Abstract

This article concerns the time-evolution, spectral structure and scaling of weak turbulence subject to rotation and stable stratification. The flow is expressed as a combination of particular solutions, referred to as modes, of the linearized governing equations without viscosity or diffusion. Modes are of two types: oscillatory ones which represent inertial-gravity waves and time-independent ones that express a non-propagating (NP) component of the flow. The presence of the NP component, which plays an active role in the dynamics apart from in the case of pure rotation, renders wave-turbulence analysis problematic because the NP mode is nondispersive. Equations are derived for the time evolution of the modal amplitudes, evolution which is due to nonlinearity and visco-diffusion. Subsequent analysis assumes that one or other (or both) of the Rossby and Froude numbers is small (weak turbulence). Given this assumption, the NP component is found to evolve independently of the wave one and a numerical scheme, similar to, though significantly different from classical DNS, is used to determine its time evolution. The treatment of the wave component assumes its amplitude large compared with the NP one, otherwise there are seemingly intractable difficulties of closure in the analysis. Given this further assumption, the wave component decouples from the NP one. Evolution equations for the wave spectra are derived using wave-turbulence analysis and are integrated numerically. As might be expected, these equations indicate that nonlinear coupling of wave modes is dominated by resonances. Results are given for both the NP and wave components.

## 1. Introduction

This paper concerns decaying homogeneous turbulence in a rotating, stably stratified fluid with constant Brunt-Vaisala frequency and assumes that the rotation vector is parallel to gravity. The decay is due to viscous and diffusive effects which are detailed later. We refer the reader to Sagaut and Cambon (2018) for a wide-ranging description of homogeneous turbulence, including, among other cases, rotation and stratification. In the present paper, one or other (or both) of the Rossby and Froude numbers is supposed small (weak turbulence), so nonlinearity is negligible over short enough time spans.

As discussed in Bartello (1995) and Cambon (2001), in the absence of nonlinearity and visco-diffusion, the flow consists of inertial-gravity wave modes, which oscillate in time, and a non-propagating (abbreviated to NP in what follows) mode which is steady. These modes form a complete set and so, even in the presence of nonlinearity and visco-diffusion, the flow can be expressed as a combination of modes. When weak nonlinearity and visco-diffusion are
Champ de recherche par défaut (multicritères)

○ <u>Champ de recherche par défaut (multicritères)</u>

○

○ <u>Revues (multicritères)</u>

○ <u>Revue : Éditeur</u>

○ <u>Revue : Titre abrégé</u>

○ <u>Revue : Titre</u>

○ <u>Revue : Identifiant interne</u>

○ <u>Revue : Date de début de publication</u>

○ <u>Revue : ISSN électronique</u>

○ <u>Revue : ISSN</u>

○ <u>Revue : Éditeur</u>

○ <u>Revue : Couleur dans SHERPA/RoMEO</u>

○ <u>Revue : État dans le référentiel</u>

○

○ <u>Colloque (multicritères)</u>

○ <u>Colloque : Titre</u>

○ <u>Colloque : Organisateur</u>

○ <u>Colloque : date de début (Année)</u>

○ <u>Colloque : date de fin (Année)</u>

○ <u>Conférence invitée (oui/non)</u>

○

○ <u>Projets ANR (multicritères)</u>

○ <u>Projet ANR : Acronyme</u>

○ <u>Projet ANR : Acronyme du programme</u>

○ <u>Projet ANR : Code décision (référence)</u>

○ <u>Projet ANR : Nom</u>

○ <u>Projet ANR : Identifiant interne</u>

○ <u>Projet ANR : État dans le référentiel</u>

○

○ <u>Projets européens (multicritères)</u>

○ <u>Projet européen : Acronyme</u>

○ <u>Projet européen : Identifiant de l'appel à projet</u>

○ <u>Projet européen : Référence</u>

○ <u>Projet européen : Nom</u>

○ <u>Projet européen : Date de fin</u>

○ <u>Projet européen : Financement</u>

○ <u>Projet européen : Date de début</u>

○ <u>Projet européen : État dans le référentiel</u>
allowed for, the mode amplitudes evolve with time, slowly compared to the oscillation period of the waves. The elucidation of the effects of this long-time evolution on the turbulence statistics in the presence of both rotation and stratification is the objective of this paper.

Weakness of turbulence has profound implications for the nonlinear dynamics, as indicated by wave-turbulence analysis, which is the usual approach for weak turbulence and has a long history (see e.g. Benney and Saffman (1966), Benney and Newell (1969), Zhakarov, Lvov and Falkovich (1992), Newell and Rumpf (2011) and Nazarenko (2011)). According to that theory, nonlinear coupling between wave modes is dominated by resonances, non-resonant interactions being suppressed, which reduces the effectiveness of nonlinearity. More recently, studies have been undertaken (see e.g. Deng and Hani (2021)) which aim to place the theory on a rigorous mathematical basis.

Galtier (2003) and Bellet et al. (2006) (henceforth referred to as [B]) used the wave-turbulence approach in the case of pure rotation, but we are unaware of any previous applications of wave-turbulence theory which allow for stratification. The difficulty is that wave-turbulence analysis requires that modes be dispersive, which is not true of the NP ones in the present problem. As a result, wave-turbulence analysis does not allow for the NP component of the flow, which must be analysed separately. In the special case $N = 2\Omega$, where $N$ is the Brunt-Vaisala frequency and $\Omega$ the rotation rate, the wave modes are also nondispersive. To allow the application of wave-turbulence analysis to the wave component, we suppose that $N$ is not too close to $2\Omega$. This condition turns out to be also required by our treatment of the NP component.

Direct numerical simulation (DNS) of the governing equations has often been used (see e.g. Orszag and Patterson (1972), Canuto et al. (1988)) to study homogeneous turbulence and specifically in the present case of rotation and stratification by Coleman, Ferziger and Spalart (1992), Smith and Waleffe (2002) and Liechtenstein, Godeferd and Cambon (2005). In this approach, the infinite flow is approximated as spatially periodic and is represented by Fourier series. However, DNS has difficulties when applied to the present problem if the turbulence is weak. This is because the intervention of nonlinearity requires evolution over many wave periods (according to wave-turbulence theory, the ratio of a typical wave period to the nonlinear time scale behaves as the square of the small parameter characterising the weakness of turbulence). The time step must be small compared with the wave period to resolve the associated oscillations, hence very many steps are needed before nonlinearity intervenes. Furthermore, the Fourier coefficients develop rapid oscillations in spectral space at large times, oscillations which need to be resolved numerically. This requires spatial periods considerably larger than the size of the large scales, otherwise numerical precision is degraded; the weaker the turbulence, the larger the required spatial periods. All of this means that the correct treatment by DNS of weak turbulence places considerable demands on computing time and memory requirements, demands which increase as the turbulence becomes weaker. With this in mind, the analytical and numerical methods developed here directly address the asymptotic limit of weak turbulence.
○ Projet européen : Identifiant interne

O

O Date de production : année

O Date de mise en ligne : année

O Date de publication : année

O Date d'écriture : année

O Date de modification du dépôt : année

O Date de dépôt : année

O Date de publication électronique : année

O

O Collection HAL (multicritères)

O Collection HAL : catégorie

O Collection HAL : Code

O Collection HAL : Nom

O Collection HAL : Identifiant interne

O

O Identifiant interne du contributeur/déposant

O Nom complet du contributeur/déposant

O

O Domaines

O Domaine primaire

O Domaine racine

O Sous-domaine niveau 1

O Sous-domaine niveau 2

O Sous-domaine niveau 3

O

O Statut du document

O Version du document

O Type de dépôt

O Type de document

O ISBN

O Numéro - référence

O Identifiant DOI

O

O Classification

O Audience

O Vulgarisation

O Comité de lecture - texte (oui ou non)

O Actes de colloque

O Référence interne
The spectral closure model EDQNM (Eddy-Damped, Quasi-Normal, Markovian) has also been used to study purely stratified and purely rotating turbulence (see Cambon and Jacquin (1989), Godeferd and Cambon (1994), Cambon, Mansour and Godeferd (1997), Godeferd and Staquet (2003)). EDQNM was introduced by Orszag (1970) for homogeneous, isotropic turbulence and later extended to other cases in which there are dynamical mechanisms, such as stratification or rotation, which render the turbulence anisotropic, even if it is initially isotropic. For the case of pure rotation, the close relationship of EDQNM and wave-turbulence theory is one of the principal subjects of [B], where it is shown that the weak-turbulence limit of a particular version of EDQNM (EDQNM 3) gives the wave-turbulence equations for the spectra. It should, however, be recognised that models like EDQNM are based on ad hoc hypotheses, such as eddy-damping and quasi-normality, whereas wave-turbulence theory follows from detailed asymptotic analysis of weak turbulence (given in appendix B for the present case), analysis which justifies closure in that limit.

The paper is organised as follows. Section 2 concerns the governing equations of the flow, their Fourier transforms, the definition of the modes, the results of modal projection (i.e. the mode amplitudes) and the mode-amplitude evolution equations. It also recognises the random nature of turbulent flow and, using ensemble averaging, defines a spectral matrix, denoted $A$, whose diagonal elements represent the energy distribution in spectral space of the different modes and whose off-diagonal elements express correlations between modes and are less important.

Section 3 introduces the assumption of weak turbulence and presents analyses of the consequences for the NP and wave components of the flow. In section 3.1, the NP component is found to evolve independently of the wave one and a numerical method, closely related to classical DNS of homogeneous turbulence, is proposed. However, although the methods are related, the flow field is here projected onto the NP modes at each time step, which distinguishes the present approach from the classical one. Note that, because we apply DNS only to the NP component, which is non-oscillatory, the problems of the classical approach, discussed above, do not arise.

Section 3.2 and appendix B derive wave-turbulence equations which describe the time evolution of the wave part of the spectral matrix. To avoid the seeming intractability which the presence of the NP component entails if it is of comparable or greater magnitude than the wave one, the derivation requires that the NP component be small compared with the wave one. The end result is a system of equations for the diagonal elements of the wave-component spectral matrix.

Finally, section 4 describes results of numerical calculations for both the NP component using DNS and the wave component using the wave-turbulence equations.

## 2. Formulation

Consider decaying, homogeneous turbulence in a rotating, stably stratified fluid having constant Brunt-Vaisala frequency, $N$, and rotation vector $\Omega$, which is supposed parallel to
○ Financement

○ Collaborations

Ajouter

Lancer la recherche

Recherche experte (SolR)

Recherche experte (SolR)

Lancer la recherche
Vers la recherche avancée

• Déposer

X

X

X

Loading...

X

Le document que vous demandez est en access restreint
jusqu'au 2024-07-11 !

Revenir a la notice

Retourner à l'accueil

[image]

[image]

[image]

[image]

Contact

• Support

Ressources

• Documentation

• FAQ

• API
gravity. It is also supposed that $N$ and $\Omega = |\Omega|$ are not both zero. When $\Omega \neq 0$ an axial direction is defined by the unit vector $\mathbf{e} = \Omega/\Omega$. On the other hand, if $\Omega = 0$, $\mathbf{e}$ is taken in the vertical direction. Henceforth, spatial coordinates, time and velocity are nondimensionalised using $L$, $(N^2 + 4\Omega^2)^{-1/2}$ and $L(N^2 + 4\Omega^2)^{1/2}$, where $L$ is a length scale characterising the initial turbulence.

Using a rotating frame of reference and Cartesian coordinates (which will often be taken such that $\mathbf{e} = (0,0,1)$), as well as the summation convention, the nondimensional Boussinesq equations of motion are

$$ \frac{\partial u_i}{\partial t} + \hat{\Omega} \varepsilon_{ijk} e_j u_k + \hat{N} \eta e_i + \frac{\partial p}{\partial x_j} = - \frac{\partial}{\partial x_j} (u_i u_j) + D_u \frac{\partial^2 u_i}{\partial x_j \partial x_j}, \quad (2.1) $$

$$ \frac{\partial u_i}{\partial x_j} = 0, \quad (2.2) $$

$$ \frac{\partial \eta}{\partial t} - \hat{N} e_i u_i = - \frac{\partial}{\partial x_j} (\eta u_j) + D_{\eta} \frac{\partial^2 \eta}{\partial x_j \partial x_j}, \quad (2.3) $$

where $\varepsilon_{ijk}$ is the alternating tensor, $\hat{\Omega} = 2\Omega(N^2 + 4\Omega^2)^{-1/2}$ and $\hat{N} = N(N^2 + 4\Omega^2)^{-1/2}$. Furthermore, $D_u = \nu(N^2 + 4\Omega^2)^{-1/2} L^{-2}$ and $D_{\eta} = \kappa(N^2 + 4\Omega^2)^{-1/2} L^{-2}$, where $\nu$ is the kinematic viscosity and $\kappa$ the diffusivity associated with the buoyancy variable $\eta$. Note that, since $\hat{N}^2 + \hat{\Omega}^2 = 1$, $\hat{N} = (\beta^2 + 1)^{-1/2}$ and $\hat{\Omega} = \beta(\beta^2 + 1)^{-1/2}$, where $\beta = \hat{\Omega}/\hat{N} = 2\Omega/N$. Thus, the nondimensional governing equations, (2.1)-(2.3), only depend on the parameters $\beta$, $D_u$ and $D_{\eta}$. For simplicity's sake, $\hat{\Omega}$, $\hat{N}$ are denoted $\Omega$, $N$ in what follows. Because we will only be working with nondimensional 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. When $N=0$ (pure rotation), (2.1) and (2.2) decouple from (2.3), so the velocity field can be studied independently of $\eta$, which becomes a passive scalar. In all other cases, there is coupling in both directions between (2.1) and (2.3).

## 2.1 Fourier transforms and modes

Defining the Fourier transforms

$$ \tilde{u}_i(\mathbf{k}) = \frac{1}{8\pi^3} \int u_i(\mathbf{x}) \exp(-i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{x}, \quad (2.4) $$

$$ \tilde{\eta}(\mathbf{k}) = \frac{1}{8\pi^3} \int \eta(\mathbf{x}) \exp(-i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{x}, \quad (2.5) $$
• OAI-PMH

• AuréHAL

## Informations

• À propos

• Données personnelles

• Mentions légales

• Accessibilité

• Conformité RGAA

## Questions juridiques

• Je publie, quels sont mes droits ?

• Loi pour une République numérique

• Stratégie de non-cession des droits

• SHERPA/RoMEO

## Portails

• Portails institutionnels

• HAL

• HAL SHS

• HAL Thèses

• MédiHAL

## CCSD

• CCSD

• Episciences

• Sciencesconf
with similar definitions of $\tilde{p}$, $\tilde{u}_i u_j$ and $\tilde{\eta} u_j$, (2.1)-(2.3) yield

$$ \frac{\partial \tilde{u}_i}{\partial t} + \Omega \varepsilon_{ijk} e_j \tilde{u}_k + N \tilde{\eta} e_i + ik_i \tilde{p} = -ik_j \tilde{u}_i u_j - D_u k^2 \tilde{u}_i, \quad (2.6) $$

$$ k_i \tilde{u}_i = 0, \quad (2.7) $$

$$ \frac{\partial \tilde{\eta}}{\partial t} - N e_i \tilde{u}_i = -ik_j \tilde{\eta} u_j - D_{\eta} k^2 \tilde{\eta}, \quad (2.8) $$

where $k = |\mathbf{k}|$.

Appendix A examines the consequences of (2.6)-(2.8). In the absence of nonlinearity and visco-diffusion, the right-hand sides of (2.6) and (2.8) are zero. There are then three linearly independent solutions, referred to as modes. These solutions have time dependence $\exp(-is\omega(\mathbf{k})t)$, where $s$ takes one of the three values $s = 0, \pm 1$,

$$ \omega(\mathbf{k}) = \frac{(N^2 k_{\perp}^2 + \Omega^2 k_{\parallel}^2)^{1/2}}{k} = (N^2 \sin^2 \theta_k + \Omega^2 \cos^2 \theta_k)^{1/2}, \quad (2.9) $$

$k_{\perp} = |\mathbf{k} \times \mathbf{e}|$ is the transverse wavenumber, $k_{\parallel} = \mathbf{k} \cdot \mathbf{e}$ is the axial wavenumber and $0 \le \theta_k \le \pi$ is the angle between $\mathbf{k}$ and $\mathbf{e}$. The modes with $s = \pm 1$ represent inertial-gravity waves, for which (2.9) is the dispersion relation, while that with $s = 0$ will be referred to as the non-propagating (NP) mode because, when regarded as a wave, it has zero group velocity. The modal solutions of (2.6)-(2.8) (with zero right-hand sides) are

$$ \tilde{u}_i = v_i^{(s)}(\mathbf{k}) \exp(-is\omega(\mathbf{k})t), \quad (2.10) $$

$$ \tilde{\eta} = \eta^{(s)}(\mathbf{k}) \exp(-is\omega(\mathbf{k})t), \quad (2.11) $$

where

$$ v_i^{(s)}(\mathbf{k}) = \sum_{l=1}^{2} u_s^{(l)}(\mathbf{k}) e_i^{(l)}(\mathbf{k}), \quad (2.12) $$

$$ \mathbf{e}^{(1)} = \frac{\mathbf{k} \times \mathbf{e}}{k_{\perp}}, \quad \mathbf{e}^{(2)} = \frac{\mathbf{k} \times \mathbf{e}^{(1)}}{k} \quad (2.13) $$

are unit vectors orthogonal to each other and to $\mathbf{k}$,
$$ \begin{pmatrix} u_s^{(1)} \\ u_s^{(2)} \\ \eta^{(s)} \end{pmatrix} = \frac{1}{2^{1/2} (N^2 k_{\perp}^2 + \Omega^2 k_{\parallel}^2)^{1/2}} \begin{pmatrix} -\Omega k_{\parallel} \\ is (N^2 k_{\perp}^2 + \Omega^2 k_{\parallel}^2)^{1/2} \\ Nk_{\perp} \end{pmatrix}. \quad (2.14) $$

when $s = \pm 1$ and

$$ \begin{pmatrix} u_0^{(1)} \\ u_0^{(2)} \\ \eta^{(0)} \end{pmatrix} = i \frac{1}{(N^2 k_{\perp}^2 + \Omega^2 k_{\parallel}^2)^{1/2}} \begin{pmatrix} Nk_{\perp} \\ 0 \\ \Omega k_{\parallel} \end{pmatrix}. \quad (2.15) $$

It follows from (2.12)-(2.15) and the definitions of $k_{\perp}$ and $k_{\parallel}$ that $v_i^{(s)}(-\mathbf{k}) = v_i^{(-s)*}(\mathbf{k})$ and $\eta^{(s)}(-\mathbf{k}) = \eta^{(-s)*}(\mathbf{k})$, where * denotes complex conjugation. Note the orthonormality relation

$$ \sum_{l=1}^{2} u_{s}^{(l)*} u_{s'}^{(l)} + \eta^{(s)*} \eta^{(s')} = \delta_{ss'}, \quad (2.16) $$

where $\delta_{ss'}$ is the Kronecker delta. Note also that the $s=0$ modes, referred to as NP here, are often described in the geophysical literature as “vortical modes” or “potential vorticity (PV) modes”.

## 2.2 Mode-amplitude equation

The modes form a complete set for $\tilde{u}_i$, $\tilde{\eta}$ satisfying (2.7). Thus,

$$ \tilde{u}_i = \sum_{s=0, \pm 1} b_s(\mathbf{k}, t) v_i^{(s)}(\mathbf{k}), \quad (2.17) $$

$$ \tilde{\eta} = \sum_{s=0, \pm 1} b_s(\mathbf{k}, t) \eta^{(s)}(\mathbf{k}), \quad (2.18) $$

where, as shown in appendix A, the modal coefficients evolve according to

$$ \frac{\partial b_s}{\partial t} + i s \omega b_s = -i k_j \left( v_i^{(s)*} \tilde{u}_i \tilde{u}_j + \eta^{(s)*} \tilde{\eta} u_j \right) - \sum_{\hat{s}=0, \pm 1} D_{s\hat{s}} b_{\hat{s}}, \quad (2.19) $$

in which

$$ D_{s\hat{s}}(\mathbf{k}) = k^2 \left( D_u \sum_{l=1}^{2} u_s^{(l)*}(\mathbf{k}) u_{\hat{s}}^{(l)}(\mathbf{k}) + D_{\eta} \eta^{(s)*}(\mathbf{k}) \eta^{(\hat{s})}(\mathbf{k}) \right) \quad (2.20) $$
is a Hermitian, positive definite matrix expressing viscosity and diffusion. Real $u_i$ and $\eta$ requires that $\tilde{u}_i(-\mathbf{k}) = \tilde{u}_i^*(\mathbf{k})$ and $\tilde{\eta}(-\mathbf{k}) = \tilde{\eta}^*(\mathbf{k})$, conditions which, according to (2.17), (2.18), $v_i^{(s)}(-\mathbf{k}) = v_i^{(-s)*}(\mathbf{k})$ and $\eta^{(s)}(-\mathbf{k}) = \eta^{(-s)*}(\mathbf{k})$, are met provided $b_s(-\mathbf{k}) = b_{-s}^*(\mathbf{k})$.

Using (2.17), (2.18) and the inverse transforms of (2.4) and (2.5), the total flow can be expressed as the sum of wave and NP components:

$$u_i = u_i^W + u_i^{NP}, \quad \eta = \eta^W + \eta^{NP}, \qquad (2.21)$$

where

$$u_i^W = \sum_{s=\pm 1} \int b_s(\mathbf{k}) v_i^{(s)}(\mathbf{k}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k}, \qquad (2.22)$$

$$\eta^W = \sum_{s=\pm 1} \int b_s(\mathbf{k}) \eta^{(s)}(\mathbf{k}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k}, \qquad (2.23)$$

$$u_i^{NP} = \int b_0(\mathbf{k}) v_i^{(0)}(\mathbf{k}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k}, \qquad (2.24)$$

$$\eta^{NP} = \int b_0(\mathbf{k}) \eta^{(0)}(\mathbf{k}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k}. \qquad (2.25)$$

Given $b_s(-\mathbf{k}) = b_{-s}^*(\mathbf{k})$, $v_i^{(s)}(-\mathbf{k}) = v_i^{(-s)*}(\mathbf{k})$ and $\eta^{(s)}(-\mathbf{k}) = \eta^{(-s)*}(\mathbf{k})$, both components are real, the $s = \pm 1$ contributions to the wave component being complex conjugates.

If nonlinearity, viscosity and diffusion were neglected, the solution of (2.19) would have the expected modal form $b_s \propto \exp(-is\omega t)$. The resulting oscillations of $b_{\pm 1}$ due to waves can be suppressed by defining $a_s = b_s \exp(is\omega t)$, which evolves according to

$$\frac{\partial a_s}{\partial t} = -ik_j \left( v_i^{(s)*} \tilde{u}_i u_j + \eta^{(s)*} \tilde{\eta} u_j \right) \exp(is\omega t) - \sum_{\hat{s}=0, \pm 1} D_{s\hat{s}} \exp(i(s-\hat{s})\omega t) a_{\hat{s}}. \quad (2.26)$$

In the absence of nonlinearity and visco-diffusion, the mode amplitudes $a_s$ are time-independent, whereas, when we later specialise to weak turbulence and small visco-diffusion they evolve slowly with time. $b_s(-\mathbf{k}) = b_{-s}^*(\mathbf{k})$ and the definition of $a_s$ imply $a_s(-\mathbf{k}) = a_{-s}^*(\mathbf{k})$.

As shown in appendix A, the nonlinear term in (2.26) can be expressed as
$$-ik_j (v_i^{(s)*} \tilde{u}_i u_j + \eta^{(s)*} \tilde{\eta} u_j) \exp(is\omega t) = \\ \sum_{s_p, s_q=0,\pm 1} \int N_{ss_p s_q}^* (\mathbf{k}, \mathbf{p}) a_{s_p}^* (\mathbf{p}) a_{s_q}^* (-\mathbf{k} - \mathbf{p}) \exp(iF_{ss_p s_q} (\mathbf{k}, \mathbf{p}) t) d^3\mathbf{p}, \quad (2.27)$$

where

$$N_{ss_p s_q} (\mathbf{k}, \mathbf{p}) = ik_j v_j^{(s)} (-\mathbf{k} - \mathbf{p}) (v_i^{(s)} (\mathbf{k}) v_i^{(s_p)} (\mathbf{p}) + \eta^{(s)} (\mathbf{k}) \eta^{(s_p)} (\mathbf{p})) \quad (2.28)$$

represents nonlinear coupling between modes and

$$F_{ss_p s_q} (\mathbf{k}, \mathbf{p}) = s\omega(\mathbf{k}) + s_p\omega(\mathbf{p}) + s_q\omega(-\mathbf{k} - \mathbf{p}). \quad (2.29)$$

A symmetrised version of (2.27), namely

$$-ik_j (v_i^{(s)*} \tilde{u}_i u_j + \eta^{(s)*} \tilde{\eta} u_j) \exp(is\omega t) = \\ \sum_{s_p, s_q=0,\pm 1} \int M_{ss_p s_q}^* (\mathbf{k}, \mathbf{p}) a_{s_p}^* (\mathbf{p}) a_{s_q}^* (-\mathbf{k} - \mathbf{p}) \exp(iF_{ss_p s_q} (\mathbf{k}, \mathbf{p}) t) d^3\mathbf{p}, \quad (2.30)$$

where

$$M_{ss_p s_q} (\mathbf{k}, \mathbf{p}) = \frac{1}{2} (N_{ss_p s_q} (\mathbf{k}, \mathbf{p}) + N_{ss_q s_p} (\mathbf{k}, -\mathbf{k} - \mathbf{p})), \quad (2.31)$$

is also derived in appendix A. Using either (2.27) or (2.30), (2.26) provides an evolution equation for the $a_s$. Both versions are employed in what follows. Note the symmetries $M_{ss_p s_q} (\mathbf{k}, -\mathbf{k} - \mathbf{p}) = M_{ss_q s_p} (\mathbf{k}, \mathbf{p})$ and $F_{ss_p s_q} (\mathbf{k}, -\mathbf{k} - \mathbf{p}) = F_{ss_q s_p} (\mathbf{k}, \mathbf{p})$. Note also that the wave vectors $\mathbf{p}$ and $\mathbf{q} = -\mathbf{k} - \mathbf{p}$, which appear in (2.27) and (2.30), satisfy the usual condition, $\mathbf{k} + \mathbf{p} + \mathbf{q} = 0$, for formation of a triad with $\mathbf{k}$.

## 2.3 The spectral matrix and energy

Here and henceforth, the random nature of turbulent flow is recognized and ensemble averaging introduced. The flow is assumed statistically homogeneous, i.e. its statistical properties are the same at all spatial locations, in particular any one-point average is independent of position. Averaging (2.1) and (2.3) eliminates the nonlinear, pressure and visco-diffusive terms, while the average of (2.2) is automatically satisfied. We have in mind that there is no mean flow, i.e. $\overline{u_i} = 0$, where the overbar denotes ensemble averaging, hence the average of (2.1) gives $\overline{\eta} = 0$. $\overline{u_i} = \overline{\eta} = 0$ leads to $\overline{a_s} = \overline{b_s} = 0$, while (2.22)-(2.25) imply that the wave and NP components are both of zero mean. Those components also inherit the statistical homogeneity of the total flow.
The ensemble-averaged, nondimensional energy density in physical space is the sum of kinetic and potential energies:

$$ \frac{1}{2} \overline{u_i u_i} + \frac{1}{2} \overline{\eta^2} . \qquad (2.32) $$

Given statistical homogeneity,

$$ \overline{\tilde{u}_i^*(\mathbf{k}) \tilde{u}_i(\mathbf{k}')} = \Phi_u(\mathbf{k}) \delta(\mathbf{k} - \mathbf{k}'), \quad \overline{\tilde{\eta}^*(\mathbf{k}) \tilde{\eta}(\mathbf{k}')} = \Phi_\eta(\mathbf{k}) \delta(\mathbf{k} - \mathbf{k}'), \quad (2.33) $$

in which $\delta$ represents the Dirac function. Hence, using the inverse transform of (2.4),

$$ \begin{aligned} \frac{1}{2} \overline{u_i u_i} &= \frac{1}{2} \iint \overline{\tilde{u}_i^*(\mathbf{k}) \tilde{u}_i(\mathbf{k}')} \exp(i(\mathbf{k}' - \mathbf{k}) \cdot \mathbf{x}) d^3\mathbf{k}' d^3\mathbf{k} \\ &= \frac{1}{2} \int \Phi_u(\mathbf{k}) d^3\mathbf{k} \end{aligned} \qquad (2.34) $$

which indicates that the kinetic-energy density in spectral space is $e_K(\mathbf{k}) = \Phi_u / 2$. Similarly, the potential-energy density in spectral space is $e_P(\mathbf{k}) = \Phi_\eta / 2$.

Again using statistical homogeneity,

$$ \overline{a_s^*(\mathbf{k}) a_{s'}(\mathbf{k}')} = A_{ss'}(\mathbf{k}) \delta(\mathbf{k} - \mathbf{k}'), \qquad (2.35) $$

where $A_{ss'}$ is the spectral matrix and $s, s'$ run over the values -1, 0 and +1. $A_{ss'}$ is Hermitian and positive semi-definite. In particular, the diagonal elements of $A_{ss'}$ are real and non-negative. Also, given $a_s(-\mathbf{k}) = a_{-s}^*(\mathbf{k})$, $A_{ss'}(-\mathbf{k}) = A_{-s',-s}(\mathbf{k})$. Using (2.12), (2.17), (2.18), (2.35), $e_i^{(l)}(\mathbf{k}) e_i^{(l)}(\mathbf{k}) = \delta_{ll'}$ and $b_s = a_s \exp(-is\omega t)$ to evaluate the averages in (2.33),

$$ e_K = \frac{1}{2} \sum_{s,s'=0,\pm 1} A_{ss'} \left( u_s^{(1)*} u_{s'}^{(1)} + u_s^{(2)*} u_{s'}^{(2)} \right) \exp(i(s-s')\omega t), \qquad (2.36) $$

$$ e_P = \frac{1}{2} \sum_{s,s'=0,\pm 1} A_{ss'} \eta^{(s)*} \eta^{(s')} \exp(i(s-s')\omega t). \qquad (2.37) $$

Note that (2.16) implies

$$ e = e_K + e_P = \frac{1}{2} \sum_{s=0,\pm 1} A_{ss} \qquad (2.38) $$
for the total energy density. Thus, $A_{ss'}$ determines the energy densities in spectral space, but it contains considerably more statistical information than that: any second-order, two-point moment involving $u_i$ and $\eta$ can be obtained knowing $A_{ss'}$.

The diagonal elements of $A_{ss'}$ can be interpreted as the energy densities of the different modes and will often be referred to as spectra, whereas the off-diagonal ones represent correlations between modes. As usual, it can be shown that the nonlinear term in (2.26) conserves the total energy, which decays due to visco-diffusion according to

$$ \frac{d}{dt} \int e(\mathbf{k}) d^3\mathbf{k} = -\Re \left( \int \sum_{s, \hat{s}=0, \pm 1} D_{ss'}^*(\mathbf{k}) \exp(i(\hat{s}-s)\omega(\mathbf{k})t) A_{\hat{s}s}(\mathbf{k}) d^3\mathbf{k} \right). \quad (2.39) $$

Finally, a spherically-averaged spectrum, $E(k)$, can be defined by averaging $e(k)$ over the sphere $|k|=k$, then multiplying by $4\pi k^2$. When integrated over $k$, $E(k)$ gives the total energy. It represents the energy distribution in spectral space, without regard for anisotropy. Obviously, $E(k)$ contains less information than $e(k)$.

The energy density in spectral space, given by (2.38), can be split into wave and NP contributions as $e = e_W + e_{NP}$, where

$$ e_W = \frac{1}{2} \sum_{s=\pm 1} A_{ss}, \quad e_{NP} = \frac{1}{2} A_{00}. \qquad (2.40) $$

The spherically-averaged spectrum can also be split into wave and NP contributions as $E = E_W + E_{NP}$, where $E_W$ and $E_{NP}$ are obtained from $e_W$ and $e_{NP}$ in the same way that $E$ follows from $e$.

### 3. Weak turbulence

From here on, we suppose that $u_i$ and $\eta$ are small. This means that linear theory applies over time intervals of $O(1)$. Small nonlinearity can nonetheless have significant cumulative effects over longer time scales, effects whose quantification is the aim of this article. However, this requires sufficiently small visco-diffusive dissipation, otherwise it will kill the turbulence before nonlinearity can intervene. Given weak turbulence and small visco-diffusion, we consider large $t$.

As discussed in the introduction, wave-turbulence analysis is the usual approach to describe weak turbulence in the presence of waves. However, it requires that the waves be dispersive. Since the frequency of the NP mode is zero, it has zero group velocity and hence is nondispersive, which rules out direct application of the wave-turbulence approach to the flow as a whole. Furthermore, when $\Omega = N$, $\omega(\mathbf{k}) = \Omega = N$ is constant, in which case the
wave mode is also nondispersive, while, if $|\Omega-N|$ is small, but nonzero, it is only weakly dispersive. Many of the steps in the analysis require $|\Omega-N| \gg t^{-1}$, an assumption which is made from here on. This assumption opens the way for wave-turbulence theory of the wave modes, but care is still needed because of the nondispersive character of the NP mode. Wave-turbulence analysis cannot be applied to the NP component of the flow, which is analysed separately.

The weakness of turbulence is expressed by small parameters, $\epsilon$ and $\epsilon_{NP}$, which respectively measure the amplitudes of the wave and NP components. The parameters need not have the same order of magnitude. Indeed, to close the system of equations for the wave component, we will later suppose that $\epsilon_{NP}$ is small compared with $\epsilon$.

In the following analysis, we suppose that neither $N$ nor $\Omega$ are zero, which excludes pure rotation or stratification. This avoids having to deal with special cases where, rather than being strictly positive, $\omega(\mathbf{k})$ is zero for particular values of $\theta_k$ ($\theta_k = \pi/2$ when $N=0$, $\theta_k = 0, \pi$ when $\Omega=0$). We will nonetheless later give results for the limiting cases $N=0$ and $\Omega=0$.

## 3.1 Evolution of the NP component

Applying (2.26) with $s=0$,

$$ \frac{\partial a_0}{\partial t} = -ik_j \left( v_i^{(0)*} \widetilde{u_i u_j} + \eta^{(0)*} \widetilde{\eta u_j} \right) - \sum_{\hat{s}=0,\pm 1} D_{0\hat{s}} \exp(-i\hat{s}\omega(\mathbf{k})t) a_{\hat{s}}. \quad (3.1) $$

Using the decomposition into wave and NP components, (2.21),

$$ u_i u_j = \underbrace{u_i^{NP} u_j^{NP}}_{NP-NP} + \underbrace{u_i^W u_j^W}_{\text{wave-wave}} + \underbrace{u_i^W u_j^{NP} + u_i^{NP} u_j^W}_{\text{wave-NP}}, \quad (3.2) $$

$$ \eta u_j = \underbrace{\eta^{NP} u_j^{NP}}_{\text{NP-NP}} + \underbrace{\eta^W u_j^W}_{\text{wave-wave}} + \underbrace{\eta^W u_j^{NP} + \eta^{NP} u_j^W}_{\text{wave-NP}}. \quad (3.3) $$

Thus, there are three contributions to the nonlinear term in (3.1): NP-NP, wave-wave and wave-NP.

The right-hand side of (3.1) represents nonlinearity and visco-diffusion. As discussed above, both are negligible over time spans of $O(1)$, but can have significant cumulative effects over the longer time scales considered here. The wave component has oscillations of periods $O(1)$, whereas the NP one is steady on such time scales. Thus, the wave-NP contributions to (3.2) and (3.3) are oscillatory and their cumulative effect following evolution according to (3.1) remains small and is neglected. Furthermore, the $\hat{s}=\pm 1$ contributions to
the visco-diffusive term are also oscillatory and hence negligible. Thus, we obtain the approximation

$$ \frac{\partial a_0}{\partial t} + D_{00}a_0 = -ik_j \left( v_i^{(0)*} \overline{u_i^{NP} u_j^{NP}} + \eta^{(0)*} \overline{\eta^{NP} u_j^{NP}} \right) + f, \quad (3.4) $$

where the first term on the right-hand side expresses nonlinear interactions between NP modes and

$$ f(\mathbf{k}, t) = -ik_j \left( v_i^{(0)*} \overline{u_i^W u_j^W} + \eta^{(0)*} \overline{\eta^W u_j^W} \right) \quad (3.5) $$

represents forcing of the NP mode by the wave component.

Following the procedure which led to (2.30), but without the NP contributions,

$$ f(\mathbf{k}, t) = \sum_{s_p, s_q = \pm 1} \int M_{0s_p s_q}^*(\mathbf{k}, \mathbf{p}) a_{s_p}^*(\mathbf{p}) a_{s_q}^*(-\mathbf{k}-\mathbf{p}) \exp(iF_{0s_p s_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p}. \quad (3.6) $$

Using (2.29), $F_{0s_p s_q}(\mathbf{k}, \mathbf{p}) = s_p \omega(\mathbf{p}) + s_q \omega(-\mathbf{k}-\mathbf{p})$. Since $\omega > 0$, the exponential in (3.6) is oscillatory when $s_p = s_q$ and such terms are therefore neglected. Assuming $|\Omega - N|t$ is large, if $s_p = -s_q$ the exponential is oscillatory with period $O(1)$ for $\mathbf{p}$ away from the surface $\omega(\mathbf{p}) = \omega(-\mathbf{k}-\mathbf{p})$. Neglecting such $\mathbf{p}$, we focus on $\mathbf{p}$ close to the surface, where the exponential oscillates slowly, potentially allowing significant cumulative effects at long times. However, although we have been unable to show it analytically, numerical calculations with different values of $\beta \neq 1$, $\mathbf{k}$ and $\mathbf{p}$ show that $M_{0,s,-s}(\mathbf{k}, \mathbf{p})$ is zero to IEEE double precision when $s \neq 0$ and $\omega(\mathbf{p}) = \omega(-\mathbf{k}-\mathbf{p})$. Assuming this result, which is a priori far from obvious, is exactly true, slow oscillations are suppressed and we suppose the forcing term in (3.4) negligible.

Recall from (2.40) that the energy density in spectral space can be split into wave and NP contributions and that the NP contribution is

$$ e_{NP} = \frac{1}{2} A_{00}. \quad (3.7) $$

Appendix D (supplementary material) shows that the first term on the right-hand side of (3.4), representing nonlinear interactions between NP modes, conserves the total energy of the NP component, which, neglecting the forcing term, decays due to visco-diffusive dissipation according to

$$ \frac{d}{dt} \int e_{NP}(\mathbf{k}) d^3\mathbf{k} = -2 \int D_{00}(\mathbf{k}) e_{NP}(\mathbf{k}) d^3\mathbf{k}. \quad (3.8) $$
Let $\epsilon_{NP}$ be a small parameter measuring the amplitude of the NP component, thus $\hat{u}_i^{NP} = u_i^{NP} / \epsilon_{NP}$ and $\hat{\eta}^{NP} = \eta^{NP} / \epsilon_{NP}$ are $O(1)$. (2.24) and (2.25) with $b_0 = a_0$ give

$$ \hat{u}_i^{NP}(\mathbf{x}) = \int \hat{a}_0(\mathbf{k}) v_i^{(0)}(\mathbf{k}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k}, \qquad (3.9) $$

$$ \hat{\eta}^{NP}(\mathbf{x}) = \int \hat{a}_0(\mathbf{k}) \eta^{(0)}(\mathbf{k}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k}, \qquad (3.10) $$

where $\hat{a}_0 = a_0 / \epsilon_{NP}$. Without the forcing term, (3.4) yields

$$ \frac{\partial \hat{a}_0}{\partial \hat{t}} + \hat{D}_{00} \hat{a}_0 = -ik_j \left( v_i^{(0)*} \overline{\hat{u}_i^{NP} \hat{u}_j^{NP}} + \eta^{(0)*} \overline{\hat{\eta}^{NP} \hat{u}_j^{NP}} \right), \qquad (3.11) $$

where $\hat{t} = \epsilon_{NP}t$ is time, scaled appropriately for evolution of the NP component, and $\hat{D}_{00} = D_{00} / \epsilon_{NP}$. Thus, the NP component evolves according to (3.9)-(3.11). Assuming visco-diffusion sufficiently small that it does not kill the turbulence before nonlinearity intervenes, the time scale for NP evolution is $O(\epsilon_{NP}^{-1})$. This time scale is generally distinct from that, $O(\epsilon^{-2})$, which, as we will see later, characterises the wave component according to wave-turbulence theory.

Before going further, we should discuss the close relationship between the theory of the NP component described here and quasi-geostrophic (QG) theory (see Pedlosky (1987)), which is one of the cornerstones in the study of atmospheric and oceanic flows since its development by Charney (1948, 1971). NP modes represent geostrophic flows (i.e. the Coriolis and pressure-gradient terms in the horizontal momentum equation balance) and the NP component of the real flow can be regarded as its mathematical projection onto such idealised flows. With this in mind, it can be shown that, for the present problem, (3.9)-(3.11) are equivalent to the three-dimensional QG approximation, thus providing support for the present approach.

Embid and Majda (1998) give asymptotic analysis of (2.1)-(2.3) for small Froude number and consider two cases. In the first, the Rossby number is small and they conclude that the flow consists of oscillatory waves and a slowly varying component which evolves according to the quasi-geostrophic approximation, in agreement with the present results. In the second case, the Rossby number is of order one, corresponding to a small value of $\beta$. This makes $N \sim 1$ and $\Omega \sim \beta$, thus the wave frequency, (2.9), is small for small $\theta_k$, i.e. waves having wave vectors near the axis $k_\perp = 0$ are slowly varying. Given such modes, one might question the neglect of the wave-NP contribution to (3.1), which was based on oscillatory waves and led to (3.4). For Rossby numbers of order one, Embid and Majda supposed the slow, horizontal component of the flow to be the sum of two parts (see their equation (3.31)), one of which is the so-called VSHF (vertically sheared horizontal flow), which is independent of $x_1$.
and $x_2$ and can be considered as a combination of wave modes with $k_⊥ = 0$. They derived evolution equations for both parts of the flow (see their equations (3.33) and (3.34)) and found that the VSHF component entered into the equation for the other part. This corresponds to coupling of the wave and NP components and neglect of the wave-NP contribution to (3.1) does not capture such coupling. However, the assumed form of the flow places significant wave energy precisely on the axis, $k_⊥ = 0$, which is not the situation we have in mind. Instead, we envisage a continuous distribution of energy over wave vectors. In that case, we believe that the neglect of the wave-NP contribution to (3.1) is justified at small $β$ by the smallness of the region in $θ_k$ for which the wave frequency is small. Whether or not the flow will eventually evolve towards a state close to that assumed by Embid and Majda is unclear.

Equation (3.11) is integrated numerically. (3.9) and (3.10) yield $\hat{u}_i^{NP}(\mathbf{x})$ and $\hat{\eta}^{NP}(\mathbf{x})$ as inverse Fourier transforms. $\hat{u}_i^{NP}$ and $\hat{\eta}^{NP}$ can then be calculated in physical space, while forward transformation gives $\hat{u}_i^{NP}$ and $\hat{\eta}^{NP}$ for use in (3.11). This resembles classical DNS of homogeneous turbulence and the numerical method employed here is based on that approach and is described in appendix F (supplementary material). However, as noted in the introduction, it differs significantly from classical DNS because of the projection onto the NP modes at each time step, a projection which is implicit in (3.11).

## 3.2 Evolution of the wave component

Terms in the sum of (2.26) with $\hat{s} \neq s$ are oscillatory and hence negligible for the long-time evolution of $a_s$. Adopting this approximation and using (2.30), (2.26) becomes

$$ \frac{\partial a_s(\mathbf{k})}{\partial t} + D_{ss}(\mathbf{k}) a_s(\mathbf{k}) = \sum_{s_p, s_q, s_{p'}, s_{q'}} \int M_{ss_p s_q}^*(\mathbf{k}, \mathbf{p}) a_{s_p}^*(\mathbf{p}) a_{s_q}^*(-\mathbf{k}-\mathbf{p}) \exp(iF_{ss_p s_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p}. \quad (3.12) $$

Appendix B describes wave-turbulence analysis based on (3.12). It is assumed that $|\Omega - N|t$ is large. Furthermore, to obtain closed equations for the wave component, it is supposed that the amplitude of the NP component is small compared with that of the wave component. The result is (B.52) for evolution of the wave elements ($A_{ss'}$, $s, s' \neq 0$) of the spectral matrix. Since the wave-component spectral energy density is given by the first of the equations (2.40), the most interesting application of (B.52) is $s' = s$, hence (B.53), which involves Cauchy principal-value integrals. However, somewhat remarkably, it is found that the sum of such contributions is zero, as is the sum of $s_q = 0$ contributions to the first term on the right-hand side of (B.53). The final result is
$$ \frac{\partial \hat{A}_{ss}(\mathbf{k})}{\partial T} + 2\hat{D}(\mathbf{k}) \hat{A}_{ss}(\mathbf{k}) = 2\pi \sum_{s_p, s_q = \pm 1} \int_{S_{s_p, s_q}(\mathbf{k})} \frac{\hat{A}_{s_p s_p}(\mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k}+\mathbf{p})|} \left( \lambda_{ss_{p} s_q}(\mathbf{k}, \mathbf{p}) \hat{A}_{s_q s_q}(-\mathbf{k}-\mathbf{p}) + \Re(\zeta_{ss_{p} s_q}(\mathbf{k}, \mathbf{p})) \hat{A}_{ss}(\mathbf{k}) \right) d^2\mathbf{p} \quad (3.13) $$

for $s \neq 0$, where $T = \varepsilon^2 t$ is time, scaled appropriately for evolution of the wave component, $\hat{A}_{ss} = A_{ss} / \varepsilon^2$ are the $O(1)$ scaled wave spectra, $\hat{D} = D / \varepsilon^2$,

$$ D(\mathbf{k}) = k^2 \frac{D_u(N^2 k_{\perp}^2 + 2\Omega^2 k_{\|}^2) + D_{\eta} N^2 k_{\perp}^2}{2(N^2 k_{\perp}^2 + \Omega^2 k_{\|}^2)} \quad (3.14) $$

is the damping coefficient of wave modes (i.e. $D_{ss}(\mathbf{k}) = D(\mathbf{k})$ for $s \neq 0$),

$$ \mathbf{c}_g(\mathbf{k}) = \nabla_{\mathbf{k}} \omega(\mathbf{k}) \quad (3.15) $$

gives the group velocity of the wave modes as $s\mathbf{c}_g$,

$$ \lambda_{ss_{p} s_q}(\mathbf{k}, \mathbf{p}) = 2 |M_{ss_{p} s_q}(\mathbf{k}, \mathbf{p})|^2, \quad (3.16) $$

$$ \zeta_{ss_{p} s_q}(\mathbf{k}, \mathbf{p}) = 4 M_{ss_{p} s_q}(\mathbf{k}, \mathbf{p}) M_{s_q s_p}^*(\mathbf{-k} - \mathbf{p}, \mathbf{k}) \quad (3.17) $$

and $S_{ss_{p} s_q}(\mathbf{k})$ is the surface in $\mathbf{p}$-space defined by $F_{ss_{p} s_q}(\mathbf{k}, \mathbf{p}) = 0$. As is apparent from (2.29), the surface $S_{ss_{p} s_q}(\mathbf{k})$ is such that $s\omega(\mathbf{k}) + s_p\omega(\mathbf{p}) + s_q\omega(-\mathbf{k}-\mathbf{p}) = 0$. This condition represents triadic resonances, triadic since it involves three wave vectors, $\mathbf{k}$, $\mathbf{p}$ and $-\mathbf{k}-\mathbf{p}$ which sum to zero, and resonant because it says that the sum of modal frequencies is zero. For this reason, $S_{ss_{p} s_q}(\mathbf{k})$ will be referred to as the resonant surface. If there are no such resonances for given $\mathbf{k}$, $s$, $s_p$ and $s_q$ (i.e. $F_{ss_{p} s_q}(\mathbf{k}, \mathbf{p}) = 0$ has no solutions in $\mathbf{p}$-space), $S_{ss_{p} s_q}(\mathbf{k})$ is the empty set and the surface integral in (3.13) should be interpreted as zero.

It is perhaps interesting to discuss (3.13) in the context of EDQNM. Eddy damping does not appear (for the case of pure rotation, it was shown in [B] that it vanishes in the limit of weak turbulence). Quasi-normality is a consequence of wave-turbulence analysis and is apparent in (3.13) via the products of spectra on the right-hand side which represent fourth-order spectral moments. The fact that, according to (3.13), the spectrum at time $t$ evolves according to its values at the same instant justifies the final, Markovian, closure hypothesis of EDQNM. An important difference between wave-turbulence theory and EDQNM is that the nonlinear term in (3.13) is a surface integral over resonant triads, whereas it is a volume integral over all triads, including non-resonant ones, according to EDQNM. For the case of pure rotation, it was shown in [B] that the volume integral of EDQNM is dominated by
resonant triads in the limit of weak turbulence, hence the wave-turbulence result is approached.

Appendix E (supplementary material) analyses the energetics of the wave component based on (3.13). It is shown that the right-hand side of (3.13) conserves the total wave-component energy, which evolves according to

$$ \frac{d}{dt} \int e_W(\mathbf{k}) d^3\mathbf{k} = -2 \int D(\mathbf{k}) e_W(\mathbf{k}) d^3\mathbf{k}. \qquad (3.18) $$

Taking the sum of (3.8) and (3.18),

$$ \frac{d}{dt} \int e(\mathbf{k}) d^3\mathbf{k} = -2 \int (D_{00}(\mathbf{k}) e_{NP}(\mathbf{k}) + D(\mathbf{k}) e_W(\mathbf{k})) d^3\mathbf{k}. \qquad (3.19) $$

This result may be compared with the exact energy equation (2.39). In the weak-turbulence limit considered here, terms in the sum of (2.39) with $\hat{s} \neq s$ are oscillatory and hence negligible. Dropping these terms gives (3.19). Thus, as regards the total energy, the present approximations agree with the weak-turbulence limit of the exact result.

Appendix C concerns the existence of solutions of $F_{ss_p s_q}(\mathbf{k}, \mathbf{p}) = 0$, where $s$, $s_p$ and $s_q$ take one of the values $\pm 1$. It is shown that the resonant surface does not exist if $s = s_p = s_q$ or $1/2 \le \beta \le 2$ (a result in agreement with Smith and Waleffe (2002), section 6.1). Thus, when $1/2 \le \beta \le 2$, the right-hand side of (3.13) is zero according to the present theory. Nonlinear interactions between wave modes are then absent at the order to which we are working, suggesting the need to go to higher order. However, that lies beyond the scope of the present work. When $\beta < 1/2$ or $\beta > 2$, the resonant surface exists for $s_p = s_q \neq s$ provided

$$ \omega(\mathbf{k}) > 2 \min(\Omega, N) \qquad (3.20) $$

and for $s_p = -s_q$ when

$$ \omega(\mathbf{k}) < |\Omega - N|. \qquad (3.21) $$

If $\beta < 1/3$ or $\beta > 3$, one or other of (3.20) or (3.21) is satisfied for all $\mathbf{k}$, hence resonant surfaces can always be found for such values of $\beta$. However, when $1/3 < \beta < 1/2$ or $2 < \beta < 3$, there are directions of $\mathbf{k}$ for which neither (3.20) nor (3.21) hold, hence no resonant surfaces. Such wave vectors are decoupled from all others according to (3.13), which predicts that these modes simply decay under visco-diffusive dissipation. Thus, spectral space can be divided into two regions. In the first (referred to as A), either (3.20) or (3.21) applies and there is nonlinear coupling between wave modes within that region. In the second (B), such coupling is absent. This distinction expresses different dynamics for the two types of modes, but is only significant for $1/3 < \beta < 1/2$ or $2 < \beta < 3$, otherwise all modes
are of type A. Note that the extent of the region in which inter-mode coupling occurs shrinks down to zero as either $\beta \nearrow 1/2$ or $\beta \searrow 2$, which is the way in which nonlinearity of the order considered here disappears as those boundaries are approached.

Consider a wave vector for which there is one or more resonant surfaces and suppose visco-diffusive dissipation small enough ($D \le O(\varepsilon^2)$) that it does not kill the wave component before nonlinearity intervenes, the time scale for evolution of the wave component implied by (3.13) is $t \sim O(\varepsilon^{-2})$, the usual one for wave turbulence. If there is no resonant surface, we expect the nonlinear evolution time to be longer. As noted earlier, the time scales for evolution of the wave and NP components are generally distinct.

Appendix G (supplementary material) describes the numerical method used to solve (3.13) for the statistically axisymmetric case. The procedure is essentially that of Bellet (2003), who used the wave-turbulence equations for the case of pure rotation. For the applications described in the next section, the initial $A_{ss}(\mathbf{k})$ are not just axisymmetric, but also symmetric under reflection in the plane $k_{\parallel}=0$. This corresponds to statistical symmetry of the underlying flow under reflection in any plane perpendicular to the rotation axis, a symmetry which, like axisymmetry, is preserved by time evolution according to the governing equations. As a result, $A_{ss}(\mathbf{k})$ remains reflection symmetric. Together, $A_{ss}(-\mathbf{k}) = A_{-s,-s}(\mathbf{k})$, axisymmetry and reflection symmetry imply $A_{ss}(\mathbf{k}) = A_{-s,-s}(\mathbf{k})$, hence $A_{11}(\mathbf{k}) = A_{-1,-1}(\mathbf{k}) = e_W(\mathbf{k})$ according to the first of equations (2.40).

Although, for the sake of simplicity in the analysis, the case $N=0$ (pure rotation) was earlier excluded, it may be interesting to consider that case because, to our knowledge, it is the only one for which wave-turbulence theory has previously been applied to the present problem. According to (2.10), (2.12), (2.14) and (2.15), if $N=0$ the velocity field is carried solely by the wave modes, while the scalar field, $\eta$, consists of NP modes alone, i.e. there is a precise correspondence between velocity and wave modes and $\eta$ and NP ones. As discussed towards the beginning of section 2, when $N=0$ the velocity field decouples from the scalar $\eta$, hence a corresponding decoupling of the wave modes from the NP ones. Under these circumstances, it is natural to leave the NP component to one side and consider the wave component alone. Furthermore, since the wave component is exactly decoupled from the NP one, the wave-turbulence equations are closed without the need for the assumption $\varepsilon_{NP} \ll \varepsilon$ made in the general case.

Focusing on the wave component, as $N \to 0$, (3.13) approaches the wave-turbulence equations used by Galtier (2003) and [B] for the case of pure rotation. Numerical integration of those equations by [B] indicate that the wave energy density, $e_W(\mathbf{k})$, develops an infinite, but integrable, singularity at the plane $k_{\parallel}=0$. This is the point of view of wave-turbulence theory, i.e. $\varepsilon=0$. On the other hand, for small, but nonzero $\varepsilon$, there is energy transfer
towards the plane $k_\parallel = 0$ (see Waleffe (1991)), leading to large, but bounded, values of $e_\omega$ for small $k_\parallel$. Thus, we see how the wave-turbulence result is approached in the limit $\epsilon \to 0$. That the singularity of $e_\omega(\mathbf{k})$ is integrable in that limit is significant because it implies that the contribution of small $k_\parallel$ to the total energy is small, rather than dominant, as suggested by some authors (e.g. Hossain (1994)).

## 4. Numerical results

Here, we describe results of numerical solutions of (3.9)-(3.11) for the NP component and (3.13) for the wave one for different values of $\beta = \Omega/N$. It should be borne in mind that the derivations of (3.9)-(3.11) and (3.13) suppose that $|\Omega - N|t$ is large, hence the avoidance of $\Omega = N$ ($\beta = 1$) in what follows. As usual in theoretical studies of turbulence, we aim for dissipation to have as little effect as possible, but it cannot be entirely removed in a numerical study because the minimum numerically resolvable length scale in physical space is nonzero. If there is an energy cascade toward smaller scales, it requires mopping up by dissipation, otherwise it accumulates at the smallest scales. Furthermore, the time scale for evolution of small scales tends to decrease as they become smaller, yielding potential problems of numerical instability if the minimum scale is insufficiently limited by dissipation.

The initial ($t = 0$) distribution of energy of the different modes is chosen to be

$$ \hat{A}_{ss}(\mathbf{k}) = k^2 \exp(-k^2), \qquad (4.1) $$

where $\hat{A}_{00} = A_{00} / \epsilon_{NP}^2$ is the $O(1)$ scaled NP spectra, while, as before, $\hat{A}_{ss} = A_{ss} / \epsilon^2$ for $s = \pm 1$ are the scaled wave spectra. Note that the spherically averaged spectra resulting from (4.1) have the form $k^4 \exp(-k^2)$, a form often used in theoretical studies of turbulence. Because visco-diffusion is small, we expect time evolution to result in the creation of an inertial range (by which we mean a range of large $k$, resulting from nonlinear energy transfer from small to large $k$ and in which nonlinearity dominates dissipation, but which is not necessarily associated with an energy cascade), followed by the establishment of a dissipative range which limits the spectral extent.

The above considerations suggest the use of a hyperviscosity to extend the inertial range (see e.g. Haugen and Brandeburg (2004)). Based on numerical experimentation, the numerical visco-diffusive coefficient used in the NP calculations was chosen as $\hat{D}_{00}(\mathbf{k}) = d_0 k^6$, while the wave-component ones employ $\hat{D}(\mathbf{k}) = dk^4$, the lesser exponent for the wave component being the result of the higher maximum wavenumber attainable for that component (due to the logarithmic wavenumber distribution for the wave-turbulence equations, described in appendix G (supplementary material), but which is unavailable in DNS). The values of $d_0$ and $d$, as well as those of the other numerical parameters, are given at the end of appendices F and G (supplementary material).
Results are first presented for the NP component, then for the wave one. When considering these results, it should be borne in mind that the total flow is the sum of both components. Each component alone only gives a partial view of the end result.

## 4.1 DNS results for the NP component

Let us begin with figure 1, which gives log-log plots of the scaled NP spectra, $\hat{A}_{00} = A_{00}(k, \theta_k)/\varepsilon_{NP}^2$, where $\theta_k$ is the angle between **k** and **e**, as a function of $k$ for different values of $\beta$ and $\theta_k$. Perhaps the most striking feature is the wiggles in the curves. These are due to two traditional limitations of DNS. Firstly, the results are based on a single run of DNS, not the many runs which would be needed to accurately compute the ensemble average required by the definition of $A_{00}$, runs which would render the calculations too computationally expensive. Thus, there are random fluctuations about the smooth curves which would presumably result from averaging. Secondly, DNS discretizes the wave vector **k**, leading to the considerable jumping around apparent for smaller values of $k$. These problems decrease in importance as $k$ increases and we focus on $k$ having larger values.

Figures 1a-h represent the scaled, NP spectra for $\hat{t} = 8$ (recall that $\hat{t} = \varepsilon_{NP}t$ is the scaled time for the NP component) for the different values of $\beta$ given in the figure caption. Each figure shows results for eight equally spaced values of $\theta_k$ between $\pi/32$ and $15\pi/32$. The dashed line represents the power law $k^{-5}$ and has the same location in all figures. For those readers more used to plots of spherically averaged spectra, we recall that the $k^{-5}$ for the spectral density shown here corresponds to behaviour like $k^{-3}$ of $E_{NP}(k)$. The symmetry $A_{00}(k, \theta_k) = A_{00}(k, \pi - \theta_k)$ allows the results for values of $\theta_k$ above $\pi/2$ to be deduced from those given.

In most cases, an inertial range with power-law exponent close to -5 is apparent. In particular, this applies for the values of $\beta$ represented by figures 1b-g for all $\theta_k$ and for figure 1a (pure stratification), apart from $\theta_k = \pi/32$. Such generic near constancy of the exponent is remarkable. However, for $\beta = 10$ (near to pure rotation) the situation is rather different and harder to interpret, since the power law observed in the other cases is no longer apparent. Indeed, it is hard to discern large-$k$ power laws from these results and the limit of small $N$ merits further study of the NP component. Furthermore, the lack of a $k^{-5}$ range in figure 1a for $\theta_k = \pi/32$ suggests that the case of small $\Omega$ and small $\theta_k$ is also worth further investigation.

In addition to the inertial-range exponent, figure 1 gives information on anisotropy of the small scales. For the $\beta$ of figures 1a-1d, the large-$k$ spectra decrease as $\theta_k$ increases, thus the small-scale NP energy is concentrated towards the poles, $\theta_k = 0$ and $\theta_k = \pi$. When $\beta = 1.4$ (figure 1e) the energy distribution is almost isotropic, whereas for larger $\beta$ the
energy is concentrated towards the equator, $\theta_k = \pi/2$. This is in agreement with the general consensus that dominant stratification favours energy transfer towards the poles, whereas dominant rotation sends it towards the equator. However, it should again be recalled that we are only considering one component of the flow here.

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>
<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /><figcaption>Figure 1: Log-log plots of the scaled NP spectra for $\hat{t}=8$ and (a) $\beta=0$, (b) $\beta=0.25$, (c) $\beta=0.5$, (d) $\beta=0.7$, (e) $\beta=1.4$, (f) $\beta=2$, (g) $\beta=4$ and (h) $\beta=10$. Each figure shows curves for eight equally spaced values of $\theta_k$ between $\pi/32$ and $15\pi/32$. The dashed line</figcaption></figure>
corresponds to the power law $\hat{A}_{00} \propto k^{-5}$. If required, $N$ and $\Omega$ follow from $\beta$ using
$N = (\beta^2 + 1)^{-1/2}$ and $\Omega = \beta(\beta^2 + 1)^{-1/2}$.

<figure><img src="image_11.png" /><figcaption>Figure 2: Contour plots of $\hat{A}_{00}$ as functions of $k_{\perp}$ and $k_{\parallel}$ for $\hat{t} = 8$ and (a) $\beta = 0.1$, (b) $\beta = 10$. There are ten contours, whose heights are logarithmically spaced from $10^{-8}$ to 1.</figcaption></figure>
<figure><img src="image_12.png" /></figure>

Figure 3: Log-log plots of the scaled, spherically-averaged NP spectra for $\beta = 0.5$ and times $\hat{t} = 0, 0.25, 0.5, \dots, 3.75, 4$. The arrow indicates the direction in which the large-$k$ spectral curves move with increasing time in the early stages.

Figure 2 provides another view of anisotropy in which contours of $\hat{A}_{00}$ are plotted in the $k_{\perp}-k_{\parallel}$ plane for the cases $\beta = 0.1$ and $\beta = 10$. We once again see that close to pure stratification the NP energy density tends towards the axial direction, whereas it concentrates near the equator as the case of pure rotation is approached.

As regards time evolution, figure 3 shows a representative example ($\beta = 0.5$) of the spherically-averaged NP energy spectrum $E_{NP}(k)$, defined earlier, at 17 equally spaced times from $\hat{t} = 0$ to $\hat{t} = 4$. It will be seen that, at early times, there is energy transfer from large to small scales, thus forming the inertial and dissipative ranges. At later times, the spectral peak moves slowly towards smaller $k$. The inertial range has approximate power-law behaviour close to $k^{-3}$, corresponding to the $k^{-5}$ of figure 1c and agreeing with the spectral power law proposed by Charney (1941).

The total NP energy decreases with time by an amount which depends on the value of $\beta$. The decrease from $\hat{t} = 0$ to $\hat{t} = 8$ is 0.7% for $\beta = 0.7$ and $\beta = 1.4$, 1% when $\beta = 0.5$ and $\beta = 2$, 3% for $\beta = 0.25$ and $\beta = 4$, 11% when $\beta = 0.1$, 13% for $\beta = \infty$, 16% when $\beta = 10$ and 33% for $\beta = 0$. Considering that these small to moderate decreases correspond to a time, $\hat{t} = 8$, much greater than the establishment time of the dissipative range, it appears that there is no energy NP energy cascade. This conjecture is reinforced by log-log plots (not shown here) of the total NP energy as a function of time which give no indication of the
power laws which might be expected if there were an energy cascade. The absence of a cascade is in agreement with Charney (1941).

In conclusion, our results for the NP component are consistent with the theoretical predictions of Charney (1941), both in terms of power laws and the lack of a cascade. They also provide detailed information on spectral anisotropy and its variation with $\beta$. It has been suggested (see e.g. Herbert, Pouquet and Marino (2014)) that the absence of a cascade of energy to large $k$ is because the cascade is towards small $k$ (an inverse cascade). This is related to the existence of two QG invariants, namely energy and potential enstrophy, in the absence of visco-diffusion. The theoretical basis of the inverse cascade for 3D QG turbulence is analogous to that of the classical 2D case. Several DNS studies (e.g. Marino, Mininni, Rosenberg and Pouquet (2013)) have claimed that the inverse cascade is strongest in the range $1/2 < \beta < 2$, though why this should be, given that this condition appears to refer to the wave component, is unclear.

## 4.2 Results for the wave component

Results for $1/2 \le \beta \le 2$ are not given, because, as discussed earlier, there are no nonlinear effects on the evolution of the wave spectra according to the present theory, hence pure linear dissipation, leading to negligible time evolution given the small dissipation coefficient used here. Two quantities, not previously introduced, are the total wave energy

$$E_{W}^{tot}(T) = \int e_{W} d^{3}k \qquad (4.2)$$

and the angular spectrum

$$E_{ang}(\theta_k, T) = 2\pi \sin \theta_k \int_{0}^{\infty} k^2 e_W dk, \qquad (4.3)$$

which represents the angular distribution of wave energy, the total energy being the integral of $E_{ang}(\theta_k, T)$ over $0 \le \theta_k \le \pi$. As discussed earlier, when $1/3 < \beta < 1/2$ or $2 < \beta < 3$, there are two distinct types of mode having decoupled and different dynamics. For this reason, the integral in (4.2) is decomposed as $E_{W}^{tot} = E_{W}^{A} + E_{W}^{B}$, where, as defined earlier, A denotes modes which are coupled by nonlinearity and B represents those which are not. Whereas $E_{W}^{B}$ is nearly constant, thanks to small modal damping of the large scales, A-modes undergo significant dissipation. Recall that, unless $1/3 < \beta < 1/2$ or $2 < \beta < 3$, all modes are of type A, hence $E_{W}^{B} = 0$ and $E_{W}^{A} = E_{W}^{tot}$, there being no need for the distinction between $E_{W}^{A}$ and $E_{W}^{tot}$. In what follows, it should also be recalled that $T = \varepsilon^2 t$, the scaled time variable appropriate to the wave component.

Beginning with the case of pure rotation, $N=0$ ($\beta = \infty$), allows comparisons with [B]. Figure 4 shows a log-log plot of the total wave energy, $E_{W}^{tot}$, as a function of time. In the early stages, before the dissipative range is established, the energy is very nearly constant.
Subsequently, a power law close to $T^{-0.7}$ is evident, suggesting an energy cascade. The exponent is not far from the value, -0.8, found by [B]. Figure 5 shows the energy spectra as a function of $k$ for $T=1$ and the same angles as figure 1. An inertial ranges is apparent and there are approximate power laws close to $k^{-4}$, the exponent identified in [B] for $\theta_k$ near $\pi/2$. The figure also illustrates the expected concentration of energy density near the equator.

Figure 6 shows a plot of $E_{ang}(\theta_k)$ in which the infinite singularity at $\theta_k = \pi/2$, identified by [B], is apparent. Of course, the existence of this singularity means that the precision of numerical results for angles close to $\theta_k = \pi/2$ is likely to be poor. To what extent this affects results away from the equator is unclear. Since the total energy is the integral of $E_{ang}(\theta_k)$, i.e. the area under the curve in figure 6, it is evident that the contribution of the singularity to the total energy is small, as concluded in [B].

<figure><img src="image_13.png" /></figure>

Figure 4: Log-log plot of the total scaled wave energy as a function of time for the case of pure rotation ($\beta = \infty$). The dashed line represents the power law $E_W^{tot} \propto T^{-0.7}$.
<figure><img src="image_14.png" /><figcaption>Figure 5: Log-log plots of the scaled wave energy spectra for $N=0$, $T=1$ and the same angles as figure 1. The dashed line represents the power law $e_w \propto k^{-4}$, which corresponds to $k^{-2}$ for the spherically averaged spectrum.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Figure 6: Scaled angular spectrum for $N=0$, $T=0.25$.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>Figure 7: Log-log plots of the total scaled energy of A-modes as a function of time for different values of (a) $\beta > 2$ and (b) $\beta < 1/2$. The straight-line asymptotes of the curves for large $T$ imply power laws of the form $T^{-\alpha}$. (a) from left to right (decreasing $\beta$), $\alpha$ has values close to 0.7, 0.3, 0.5 and 0.5, (b) $\alpha$ is close to 0.4 for $\beta = 0$ and to 0.5 for the other values of $\beta$.</figcaption></figure>
For cases other than $N=0$, the assumption $\epsilon_{NP} \ll \epsilon$ is required for closure of the wave-
turbulence equations. Furthermore, if $1/3 < \beta < 1/2$ or $2 < \beta < 3$, the distinction, discussed
earlier, is made between modes which undergo nonlinear coupling (component A) and others
(component B) which do not. Figure 7 shows log-log plots of the total wave energy of
component A ($E_w^A$) as a function of time for different values of $\beta$, including the case of pure
rotation, already covered by figure 4, but which is reproduced here for comparison with the
other cases. The exclusion of decoupled modes (component B) due to their different
dynamics results in the lower initial value when $\beta = 2.5$ and $\beta = 0.4$. For all values of $\beta$,
there is a range of times in which the total energy is very nearly constant, during which an
inertial range is created by transfer from small to large wavenumber, followed by a transition
to a power law, $T^{-\alpha}$, once the dissipative range is established. This behaviour suggests the
existence of an energy cascade for the wave component. The scaled time required for creation
of the dissipative range increases as $\beta$ decreases for $\beta > 2$, whereas it does not vary greatly
for $\beta < 1/2$.

<figure><img src="image_17.png" /></figure>
<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>
<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /></figure>
<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>
<figure><img src="image_24.png" /><figcaption>Figure 8: Angular spectra for (a) $β = 0$, (b) $β = 0.1$, (c) $β = 0.2$, (d) $β = 0.3$, (e) $β = 0.4$, (f) $β = 2.5$, (g) $β = 3$, (h) $β = 5$, (i) $β = 10$. Figures (a)-(e), (h) and (i) are for $T = 2$, (f) is for $T = 30$ and (g) for $T = 10$. Figures (b)-(i) show the range of $θ_k$ for which (3.20) applies, while (e) and (f) also indicate the range, $B$, of decoupled modes, for which neither (3.20) nor (3.21) hold. Concerning (3.20) and (3.21), it may help to recall that $N = (β^2 + 1)^{-1/2}$, $Ω = β(β^2 + 1)^{-1/2}$ and $ω(k) = (N^2 \sin^2 θ_k + Ω^2 \cos^2 θ_k)^{1/2}$.</figcaption></figure>
Figure 8 shows plots of the angular spectrum for different values of $\beta$ and, for each $\beta$, a $T$ sufficiently large that the dissipative range is established. Perhaps the most striking feature is the discontinuity in figures 8b-i at the boundaries of the range of $\theta_k$ defined by (3.20). This occurs because the resonant surface with $s_p = s_q = -s$ is present when (3.20) holds, but is absent once the boundary is crossed. This discontinuity would no doubt be revealed to be a region of rapid spectral variation by asymptotic analysis of the limit $\epsilon \to 0$ close to the boundary, but this is not attempted here. Also note the range of $\theta_k$ in figures 8e and 8f indicated by B. This corresponds to the decoupled modes discussed earlier. Given small dissipation, such modes maintain their initial value, $E_{ang}(\theta_k) = 3\pi^{3/2}\epsilon^2 \sin\theta_k / 4$.

Figure 8a (pure stratification) shows the expected large energy density near the poles, though, recalling that the total wave energy is the integral of $E_{ang}$, those regions do not dominate the total energy. Figures 8b-8i indicate that the energy density is considerably lower when (3.20) is satisfied, leading to the expected concentration near the poles when $\beta < 1/2$ and near the equator when $\beta > 2$. Taken together, figures 8a-8i and 6 illustrate the evolution of the wave component as $\beta$ increases.

<figure><img src="image_25.png" /></figure>

Figure 9: Log-log plot of the scaled spectra as a function of $k$ for $\beta = 3$, $T = 10$ and the same angles as figure 1. The dashed straight line represents the power law $e_w \propto k^{-4}$.
<figure><img src="image_26.png" /><figcaption>Figure 10: Log-log plot of the scaled spectra as a function of $k$ for $\beta = 0$, $T = 2$ and the same angles as figure 1. The dashed straight line represents the power law $e_W \propto k^{-4}$.</figcaption></figure>

Figure 9 shows a representative log-log plot of the wave spectrum as a function of $k$ for different values of $\theta_k$. The density of dashes on the curves increases with $\theta_k$, allowing the identification of particular $\theta_k$, while the dashed straight line represents the power law $e_W \propto k^{-4}$. Similar results are obtained for all $\beta > 2$, but the picture, as regards power laws with respect to $k$, is less clear when $\beta < 1/2$, as illustrated by figure 10. Here, inertial-range power laws as a function of $k$ are hard to distinguish and are, at best, rough approximations. For some values of $\beta < 1/2$, there are clearer power laws, but, overall, the results indicate that the wave spectra cannot be relied on to even follow approximate inertial-range power laws with respect to $k$ when $\beta < 1/2$.

In summary, temporal power laws for the total wave energy are quite clear and suggest an energy cascade, while the angular spectra show that the existence/nonexistence of the resonant surface with $s_p = s_q = -s$ leads to polar concentration of energy density when $\beta < 1/2$ and equatorial concentration for $\beta > 2$. However, inertial-range power laws for the $k$-spectra are, at best, approximate.

## 5. Conclusions

To our knowledge, this is the first analytically based study of weak rotating/stratified turbulence, other than those of Galtier (2003) and [B] for the case of pure rotation. The results depend strongly on the ratio, $\beta$, of twice the rotation rate to the Brunt-Vaisala frequency. The theory in section 2 uses Fourier analysis and modal decomposition to express the flow as
a sum of wave ($s=\pm1$) and NP ($s=0$) components, resulting in the mode amplitude equation (2.26), which has terms representing nonlinear interactions between modes and visco-diffusive dissipation. The nonlinear term can be expressed in terms of the modal amplitudes, $a_s$, using either (2.27) or (2.30), equivalent expressions which are both employed at different points in the analysis, as is the original expression in (2.26).

Section 2 also introduces the spectral matrix, $A_{ss'}(\mathbf{k})$, where $\mathbf{k}$ is the wave vector resulting from Fourier transformation. The diagonal elements of this matrix represent the energy density of the different modes in spectral space, whereas the off-diagonal ones correspond to modal correlations and are of lesser importance. The wave and NP spectral energy densities are given by (2.40).

The consequences of weakness of the turbulence are examined in section 3. To avoid difficulties in the analysis, we suppose that $\beta$ is not close to 1. It is found that the NP component evolves according to equations (3.9)-(3.11) independently of the wave component. These equations are equivalent to the three-dimensional, quasi-geostrophic approximation, which is one of the cornerstones in the study of atmospheric and oceanic flows. The wave component is treated using wave-turbulence analysis in appendix B. To rid the wave equations of the NP spectra, and hence close the equations for the wave component, the NP amplitude, $\varepsilon_{NP}$, is assumed small compared to that of the waves, $\varepsilon$. The result is equation (3.13).

Section 4 gives results of numerical solution of the evolution equations derived in section 3, namely (3.9)-(3.11) for the NP component and (3.13) for the wave one. Both (3.11) and (3.13) involve dissipation coefficients, $\hat{D}_{00}(\mathbf{k})$ for (3.11) and $\hat{D}(\mathbf{k})$ for (3.13). As usual in studies of turbulence, dissipation is chosen small for the large scales ($k=|\mathbf{k}|=O(1)$), but increases with $k$. Hyperviscous dissipation is used in an attempt to extend the expected inertial ranges. The initial conditions are expressed by (4.1), where $\hat{A}_{00} = A_{00} / \varepsilon_{NP}^2$ and $\hat{A}_{ss} = A_{ss} / \varepsilon^2$ for $s=\pm1$ are the NP and wave spectra, scaled to be of order 1. Such scaling is also used for time, the scaled time being $\hat{t} = \varepsilon_{NP}\hat{t}$ for the NP component and $T = \varepsilon^2\hat{t}$ for the wave one. From a numerical point of view, the advantage of these scaled variables is that the amplitudes, $\varepsilon_{NP}$ and $\varepsilon$, are eliminated from the problem.

Figures 1-3 show NP results. Figure 1 gives results for the NP spectra as a function of $k$ at sufficiently large scaled time that the dissipative range is established and different values of $\beta$ and $\theta_k$, the angle between the wave vector and the axis of rotation. With a few exceptions, noted in the discussion of the figure, it indicates an inertial range close to $k^{-5}$ for the NP energy density. The exceptions include the case $\beta=10$, near to the case of pure rotation, a limit which may merit further study of the NP component. Figure 1 also indicates that dominant stratification favours energy transfer towards the poles in $\mathbf{k}$-space, whereas dominant rotation sends it towards the equator. Figure 2 gives a different view of the same
trend, while figure 3 illustrates the formation of the inertial and dissipative ranges by energy transfer from large to small scales.

Perhaps the most significant conclusion concerning the NP component arises from the observation that the overall energy dissipation is generally small, even when the dissipation range is very well developed. This suggests that there is no energy cascade associated with this component. If true, this implies it can persist to much longer times than might otherwise be thought given the time scaling $t \sim \epsilon_{NP}^{-1}$ for its evolution. Both the lack of a cascade and the appearance of $k^{-5}$ power laws are in agreement with the theoretical results of Charney (1971) for the three-dimensional, quasi-geostrophic approximation.

Section 4.2 gives results for the wave component, the validity of which depend on the assumption $\epsilon_{NP} \ll \epsilon$ unless $N=0$ ($\beta = \infty$, pure rotation). In the latter case, our results are in general agreement with those of [B]. In particular, Figure 4 shows the total wave energy as a function of time. In the early stages, before the dissipative range is established, the energy is very nearly constant. Subsequently, a power law close to $T^{-0.7}$ is evident, suggesting an energy cascade. The exponent is not far from the value, -0.8, found by [B]. Figure 5 indicates inertial-range behaviour with a power law close to $k^{-4}$, the exponent identified in [B] for $\theta_k$ near $\pi/2$. The figure also shows the expected concentration of energy density near the equator. This concentration is even more evident in figure 6, which plots the angular spectrum of the wave mode, $E_{\text{ang}}(\theta_k)$, defined by (4.3) and whose integral over $0 \le \theta_k \le \pi$ gives the total wave energy. The figure illustrates the infinite singularity found by [B] and also the conclusion that it does not dominate the total energy (the area under the curve).

Turning to cases other than pure rotation, it is important to recognise that the nonlinear terms in the wave-turbulence equation (3.13), which governs the wave component according to the present theory, are integrals over the resonant surfaces, $S_{ss_p s_q}$, defined in section 2. When $1/2 \le \beta \le 2$, there are no such surfaces (as shown previously by Smith and Waleffe (2002)), hence pure linear dissipation. Given small dissipation of the large scales, time evolution is negligible and results for such values of $\beta$ are not shown. The vanishing of the nonlinear terms for this range of $\beta$ suggests going to higher order in $\epsilon$, but this is beyond the scope of this article. If $1/3 < \beta < 1/2$ or $2 < \beta < 3$, there are values of $\theta_k$ for which resonant surfaces do not exist. Thus, there are modes which are nonlinearly coupled to others, type A, while the remainder, type B, are decoupled. Type B modes undergo pure linear dissipation, hence negligible spectral time evolution.

Figure 7 shows the time evolution of the total energy of component A for different values of $\beta$. Component B is excluded, in recognition of the different dynamics of components A and B, the energy of component B being nearly constant. This exclusion only matters for $\beta = 2.5$ and $\beta = 0.4$, since all modes are of type A for the other cases covered by the figure. In all cases, there is a range of times in which the energy is very nearly constant and during which an inertial range is created by transfer from small to large wavenumber, followed by a
transition to a temporal power law once the dissipative range is established. This behaviour
suggests the existence of an energy cascade for the wave component. The time required for
creation of the dissipative range increases significantly as $\beta$ decreases for $\beta > 2$, whereas it
does not vary greatly for $\beta < 1/2$.

Figure 8 shows angular spectra for different $\beta$ at scaled times sufficiently large that the inertial and dissipative ranges are established. Modes of type B exist for $\beta = 2.5$ and $\beta = 0.4$ and, as expected, the spectra are very nearly equal to their initial values in the corresponding ranges of $\theta_k$ (see figures 8e and 8f). For $\beta = 0$ (pure stratification), the spectrum is largest near the poles, though the contribution to the total energy from other angles is comparable. When $\beta \neq 0$, there are spectral discontinuities at the boundaries of the range defined by (3.20), the range of $\theta_k$ in which the resonant surface $s_p = s_q = -s$ exists. Exiting the given range, that resonant surface abruptly ceases to exist, leading to the discontinuity. Of course, we do not expect a real spectral discontinuity, but rather a narrow range of $\theta_k$, whose width tends to zero as $\varepsilon \to 0$ (analogous to a shock wave in a compressible fluid as the dissipation goes to zero). One might expect similar behaviour to arise from the resonant surfaces $s_p = -s_q$ at the boundaries of the range (3.21). However, it can be shown that the surface $s_p = s_q = -s$ disappears by going to $|\mathbf{p}| \to \infty$, whereas, when $s_p = -s_q$, the surfaces shrink down to a point. The latter behaviour makes the contributions to (3.13) tend to zero as the boundaries of (3.21) are approached, hence no discontinuity arises, whereas the former allows the integral in (3.13) with $s_p = s_q = -s$ to tend to a nonzero limit at the boundaries of (3.20), leading to a discontinuity.

As apparent from figure 8, for $\beta \neq 0$, the energy density is considerably greater outside the range defined by (3.20). This leads to the expected higher density near the poles for $\beta < 1/2$ and near the equator when $\beta > 2$. However, the cause identified here, namely the appearance/disappearance of the resonant surface with $s_p = s_q = -s$, is perhaps surprising. It appears that the presence of that surface leads to energy transfer away from the corresponding $\mathbf{k}$-

Concerning the behaviour of the developed wave spectra as a function of $k$, there are
rough power laws $e_w \sim k^{-4}$ when $\beta > 2$, as illustrated by figure 9. However, such power laws
are hard to identify for $\beta < 1/2$ (see e.g. figure 10). In conclusion, inertial-range power laws
for the wave $k$-spectra are approximate at best.

Finally, there remain at least two open questions. Firstly, what happens to the wave
component of weak turbulence if the condition $\varepsilon_{NP} \ll \varepsilon$, used for closure of the wave-
turbulence equations (unless $N=0$), does not apply? Given the apparent lack of an energy
cascade for the NP component and sufficiently small dissipation, it maintains its energy,
while, according to the present results, the wave component decays due to a cascade. Starting
from initial conditions such that $\epsilon_{NP} \ll \epsilon$, it appears that $\epsilon$ will eventually become comparable with $\epsilon_{NP}$, at which point the treatment of the wave component used here no longer applies. The fate of the wave component is also unclear when $\epsilon_{NP} \sim \epsilon$ or $\epsilon_{NP} \gg \epsilon$ initially. However, the results concerning the NP component are unaffected because they do not depend on the assumption $\epsilon_{NP} \ll \epsilon$.

Another question concerns the behaviour of the wave component at or near the case of pure rotation (large $\beta$). As discussed earlier, in the case of pure rotation and according to wave-turbulence theory, there is an infinite singularity of the wave energy density at $k_\parallel = 0$. This suggests that the limit of small Rossby number is singular for $k_\parallel / k_\perp \ll 1$ without stratification. When small stratification is included, it seems likely that the infinite singularity is removed, but the detailed asymptotics in the double limit in which the Rossby number goes to zero and $\beta \to \infty$ remain an open question. In fact, to our knowledge, even the small Rossby number limit with pure rotation has yet to be clarified analytically.

**Declaration of interests:** The authors report no conflict of interest.

## Appendix A: Modal projection

Since the unit vectors $e^{(1)}$ and $e^{(2)}$, given by (2.13), are orthogonal to each other and to $\mathbf{k}$, (2.7) implies

$$ \tilde{u}_i = \sum_{l=1}^{2} \tilde{u}^{(l)} e_i^{(l)}, \qquad (A.1) $$

where $\tilde{u}^{(l)} = e_i^{(l)} \tilde{u}_i$. Multiplying (2.6) by $e_i^{(l)}$ ($l=1,2$) and using $\tilde{u}^{(l)} = e_i^{(l)} \tilde{u}_i$ and $e_i^{(l)} k_i = 0$ gives

$$ \frac{\partial \tilde{u}^{(l)}}{\partial t} + \Omega e_i^{(l)} \epsilon_{ijk} e_j \tilde{u}_k + N \tilde{\eta} e_i^{(l)} e_i = R^{(l)}, \qquad (A.2) $$

where

$$ R^{(l)} = -ik_j e_i^{(l)} \tilde{u}_i u_j - D_u k^2 \tilde{u}^{(l)}. \qquad (A.3) $$

The rotational term in (A.2) can be simplified using (2.7), (2.13) and the vector identity $(\mathbf{k} \times \mathbf{a}) \cdot (\mathbf{e} \times \tilde{\mathbf{u}}) = (\mathbf{k} \cdot \mathbf{e})(\mathbf{a} \cdot \tilde{\mathbf{u}}) - (\mathbf{k} \cdot \tilde{\mathbf{u}})(\mathbf{a} \cdot \mathbf{e})$. Thus, we obtain

$$ e_i^{(l)} \epsilon_{ijk} e_j \tilde{u}_k = \frac{k_\parallel}{k_\perp} e_i \tilde{u}_i, \qquad (A.4) $$
$$e_i^{(2)} \varepsilon_{ijk} e_j \tilde{u}_k = \frac{k_{\parallel}}{k} e_i^{(1)} \tilde{u}_i = \frac{k_{\parallel}}{k} \tilde{u}^{(1)}. \qquad (A.5)$$

(2.13) implies $e^{(2)} = (kk_{\perp})^{-1} k \times (k \times e)$. Using the vector identity $a \times (b \times c) = (a \cdot c)b - (a \cdot c)c$ gives $e^{(2)} = (kk_{\perp})^{-1}(k_{\parallel}k - k^2e)$, hence $e_i e_i^{(2)} = (kk_{\perp})^{-1}(k_{\parallel}^2 - k^2) = -k_{\perp}/k$, where we have used $k^2 = k_{\perp}^2 + k_{\parallel}^2$. Employing (A.1) and $e_i e_i^{(1)} = 0$ (which follows from the first of equations (2.13)), we find $e_i \tilde{u}_i = -k_{\perp} \tilde{u}^{(2)}/k$ so (A.4) gives

$$e_i^{(1)} \varepsilon_{ijk} e_j \tilde{u}_k = -\frac{k_{\parallel}}{k} \tilde{u}^{(2)}. \qquad (A.6)$$

Thus, the rotational term in (A.2) can be expressed using (A.5) and (A.6). $e_i^{(1)}e_i = 0$ and $e_i^{(2)}e_i = -k_{\perp}/k$ also allow the evaluation of the term in (A.2) which represents stratification. The result is

$$\frac{\partial \tilde{u}^{(1)}}{\partial t} - \frac{\Omega k_{\parallel}}{k} \tilde{u}^{(2)} = R^{(1)}, \qquad (A.7)$$

$$\frac{\partial \tilde{u}^{(2)}}{\partial t} + \frac{\Omega k_{\parallel}}{k} \tilde{u}^{(1)} - \frac{N k_{\perp}}{k} \tilde{\eta} = R^{(2)}. \qquad (A.8)$$

Finally, $e_i \tilde{u}_i = -k_{\perp} \tilde{u}^{(2)}/k$ and (2.8) yield

$$\frac{\partial \tilde{\eta}}{\partial t} + \frac{N k_{\perp}}{k} \tilde{u}^{(2)} = R^{(3)}, \qquad (A.9)$$

where

$$R^{(3)} = -ik_j \tilde{\eta} u_j - D_{\eta} k^2 \tilde{\eta}. \qquad (A.10)$$

The quantities $\tilde{u}^{(1)}$ and $\tilde{\eta}$ evolve according to (A.7)-(A.9).

In the absence of nonlinearity, viscosity and diffusion, the right-hand sides of (A.7)-(A.9) are zero. Looking for solutions of the form $\exp(-i\sigma t)$ gives the eigenvalue problem

$$\left( \begin{array}{ccc} 0 & i\Omega k_{\parallel} / k & 0 \\ -i\Omega k_{\parallel} / k & 0 & iNk_{\perp} / k \\ 0 & -iNk_{\perp} / k & 0 \end{array} \right) \left( \begin{array}{c} \tilde{u}^{(1)} \\ \tilde{u}^{(2)} \\ \tilde{\eta} \end{array} \right) = \sigma \left( \begin{array}{c} \tilde{u}^{(1)} \\ \tilde{u}^{(2)} \\ \tilde{\eta} \end{array} \right). \qquad (A.11)$$
Given that the matrix is Hermitian, the eigenvalues are real and there is a complete set of orthogonal eigenvectors. The eigenvalues are $\sigma = s\omega(\mathbf{k})$, where s takes one of the three values $s = 0, \pm 1$ and $\omega(\mathbf{k})$ is given by (2.9). The normalised eigenvectors are given by (2.14) and (2.15), while (2.16) reflects orthogonality and normalisation. Using (A.1), we obtain (2.10) and (2.11), hence the modal solutions described in the main text.

Completeness of the eigenvectors of (A.11) implies

$$ \begin{pmatrix} \tilde{u}^{(1)} \\ \tilde{u}^{(2)} \\ \tilde{\eta} \end{pmatrix} = \sum_{s=0,\pm 1} b_s(\mathbf{k},t) \begin{pmatrix} u_s^{(1)} \\ u_s^{(2)} \\ \eta^{(s)} \end{pmatrix}, \qquad (A.12) $$

hence (2.17) and (2.18) according to (A.1). Using (A.11) and (A.12), (A.7)-(A.9) imply

$$ \sum_{s'=0,\pm 1} \left( \frac{\partial b_{s'}}{\partial t} + i s' \omega b_s \right) \begin{pmatrix} u_{s'}^{(1)} \\ u_{s'}^{(2)} \\ \eta^{(s')} \end{pmatrix} = \begin{pmatrix} R^{(1)} \\ R^{(2)} \\ R^{(3)} \end{pmatrix}. \qquad (A.13) $$

Left-multiplying (A.13) by the row vector $(u_s^{(1)*}, u_s^{(2)*}, \eta^{(s)*})$ and using the orthonormality relation (2.16),

$$ \frac{\partial b_s}{\partial t} + i s \omega b_s = \sum_{l=1}^{2} u_s^{(l)*} R^{(l)} + \eta^{(s)*} R^{(3)}. \qquad (A.14) $$

Employing (A.3), (A.10) and (A.12) gives (2.19).

It remains to derive (2.27) and (2.30). The inverse transform of (2.4) and real $u_i$ imply

$$ u_i(\mathbf{x}) = \int \tilde{u}_i^*(\mathbf{p}) \exp(-i\mathbf{p}\cdot\mathbf{x}) d^3\mathbf{p}, \qquad (A.15) $$

hence

$$ \begin{aligned} u_i u_j &= \iint \tilde{u}_i^*(\mathbf{p}) \tilde{u}_j^*(\mathbf{q}) \exp(-i(\mathbf{p}+\mathbf{q})\cdot\mathbf{x}) d^3\mathbf{q} d^3\mathbf{p} \\ &= \iint \tilde{u}_i^*(\mathbf{p}) \tilde{u}_j^*(-\mathbf{k}-\mathbf{p}) \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k} d^3\mathbf{p} . \qquad (A.16) \\ &= \iint \tilde{u}_i^*(\mathbf{p}) \tilde{u}_j^*(-\mathbf{k}-\mathbf{p}) d^3\mathbf{p} \exp(i\mathbf{k}\cdot\mathbf{x}) d^3\mathbf{k} \end{aligned} $$

It follows that

$$ \tilde{u}_i u_j = \int \tilde{u}_i^*(\mathbf{p}) \tilde{u}_j^*(-\mathbf{k}-\mathbf{p}) d^3\mathbf{p} \qquad (A.17) $$
and, similarly,

$$ \tilde{\eta u_j} = \int \tilde{\eta}^*(\mathbf{p}) \tilde{u}_j^*(\mathbf{-k-p}) d^3\mathbf{p}. \qquad (A.18) $$

Employing (2.17) and (2.18), (A.17), (A.18) and the definition of $a_s$ give (2.27).

Changing the integration variable in (2.27) to $-\mathbf{k}-\mathbf{p}$, using $F_{ss_p s_q}(\mathbf{k}, -\mathbf{k}-\mathbf{p}) = F_{ss_q s_p}(\mathbf{k}, \mathbf{p})$ and permuting the summation indices $s_p \leftrightarrow s_q$,

$$ -ik_j (v_i^{(s)} \tilde{u}_i \tilde{u}_j + \eta^{(s)} \tilde{\eta} u_j) \exp(is\omega t) = \qquad (A.19) $$

$$ \sum_{s_p, s_q = 0, \pm 1} N_{ss_p s_q}^*(\mathbf{k}, -\mathbf{k}-\mathbf{p}) a_{s_p}^*(\mathbf{p}) a_{s_q}^*(\mathbf{-k}-\mathbf{p}) \exp(iF_{ss_p s_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p}. $$

Taking the sum of (2.27) and (A.19) and dividing by two gives (2.30).

## Appendix B: Wave-turbulence analysis

Statistical homogeneity implies

$$ \overline{a_s(\mathbf{k}) a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}''')} = \Theta_{ss_p s_q}(\mathbf{k}, \mathbf{p}) \delta(\mathbf{k} + \mathbf{p} + \mathbf{k}''). \qquad (B.1) $$

The ensemble average of the complex conjugate of (3.12) multiplied by $a_s(\mathbf{k}')$ gives

$$ \overline{\frac{\partial a_s^*(\mathbf{k})}{\partial t}} a_{s'}(\mathbf{k}') = (\tau_{ss'}(\mathbf{k}) - D_{ss}(\mathbf{k}) A_{ss'}(\mathbf{k})) \delta(\mathbf{k} - \mathbf{k}'), \qquad (B.2) $$

where

$$ \tau_{ss'}(\mathbf{k}) = \sum_{s_p, s_q = 0, \pm 1} \int M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) \Theta_{s' s_p s_q}(\mathbf{k}, \mathbf{p}) \exp(-iF_{ss_p s_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p}. \qquad (B.3) $$

In deriving (B.2), (2.35), (B.1) and the fact that $D_{ss}$ is real have been used, while, given the Dirac function in (B.2), $\Theta_{s' s_p s_q}(\mathbf{k}', \mathbf{p})$ has been replaced by $\Theta_{s' s_p s_q}(\mathbf{k}, \mathbf{p})$ in (B.3). Permuting $s \leftrightarrow s'$ and $\mathbf{k} \leftrightarrow \mathbf{k}'$, the complex conjugate of (B.2) yields

$$ \begin{aligned} \overline{a_s^*(\mathbf{k}) \frac{\partial a_{s'}(\mathbf{k}')}{\partial t}} &= (\tau_{s's}(\mathbf{k}') - D_{s's}(\mathbf{k}') A_{s's}(\mathbf{k}')) \delta(\mathbf{k} - \mathbf{k}'), & (B.4) \\ &= (\tau_{s's}(\mathbf{k}) - D_{s's}(\mathbf{k}) A_{s's}(\mathbf{k})) \delta(\mathbf{k} - \mathbf{k}') \end{aligned} $$

where we have employed the Hermitian character of $A_{ss'}$ and reality of $D_{s's}$. Taking the time derivative of (2.35), (B.2) and (B.4) imply
$$ \frac{\partial A_{ss'}(\mathbf{k})}{\partial t} + (D_{ss}(\mathbf{k}) + D_{s's'}(\mathbf{k})) A_{ss'}(\mathbf{k}) = \tau_{ss'}(\mathbf{k}) + \tau_{s's}^*(\mathbf{k}). \quad (B.5) $$

(B.5) governs the evolution of the spectral matrix, evolution which we expect to be slow for weak turbulence with small visco-diffusion. Thus, we consider large $t$. The appearance of the third-order spectral moments, represented by $\Theta$ in (B.3), is a consequence of the usual turbulence closure problem.

The time derivative of (B.1) yields

$$ \delta(\mathbf{k} + \mathbf{p} + \mathbf{k}'') \frac{\partial \Theta_{s's_p s_q}(\mathbf{k}, \mathbf{p})}{\partial t} = \frac{\overline{\partial a_{s'}(\mathbf{k})}}{\partial t} a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'') + \frac{\overline{\partial a_{s_p}(\mathbf{p})}}{\partial t} a_{s'}(\mathbf{k}) a_{s_q}(\mathbf{k}'') + \frac{\overline{\partial a_{s_q}(\mathbf{k}'')}}{\partial t} a_{s'}(\mathbf{k}) a_{s_p}(\mathbf{p}) \quad (B.6) $$

Using (3.12) and (B.1),

$$ \overline{\frac{\partial a_{s'}(\mathbf{k})}{\partial t} a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'')} = \sum_{s_p, s_q, s_p' = 0, \pm 1} \int \overline{M_{s's_p's_q}^*(\mathbf{k}, \mathbf{p}') \exp(iF_{s's_p's_q}(\mathbf{k}, \mathbf{p}')t)} \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'')} d^3\mathbf{p}' \\ - \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') D_{s's'}(\mathbf{k}) \Theta_{s's_p s_q}(\mathbf{k}, \mathbf{p}) \quad (B.7) $$

The fourth-order moments in (B.7) are treated as follows.

Suppose, as is often done in theoretical studies of homogeneous turbulence, that turbulent quantities are statistically independent at sufficiently large separations and consider fourth-order, four-point moments of $u_i(\mathbf{x})$ and $\eta(\mathbf{x})$. Let $r$ denote the largest of the distances between two points. The four-point moments go to zero as $r \to \infty$, unless the points form two pairs, each pair having bounded separation as the distance between the pairs goes to infinity. In that case, the four-point moments are the product of two-point moments, one from each pair, hence a quasi-normal limit applies at large $r$. The four-point moments can thus be written as the sum of their quasi-normal values and another component, often referred to as the cumulant correction, which tends to zero as $r \to \infty$. In spectral terms, the corresponding result is

$$ \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'')} = \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'')} + \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_q}(\mathbf{k}'')} + \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}(\mathbf{k}'')} a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_p}(\mathbf{p}) + K_{s_p s_q s_p s_q}(\mathbf{k}, \mathbf{p}, \mathbf{p}') \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') \quad (B.8) $$
where $K_{s'_p s'_q s_p s_q}$ represents the cumulant correction and the remainder of the right-hand side is the quasi-normal contribution. Using (2.35) and $a_s(-\mathbf{k}) = a_{-s}(\mathbf{k})$,

$$ \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} \overline{a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'')} = A_{s_p', -s_q'}(\mathbf{p}') A_{-s_p, s_q}(-\mathbf{p}) \delta(\mathbf{k}) \delta(\mathbf{p}+\mathbf{k}') \quad (B.9) $$

$$ \overline{a_{s_p}^*(\mathbf{p}') a_{s_p}(\mathbf{p})} \overline{a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_q}(\mathbf{k}'')} = A_{s_p', s_p}(\mathbf{p}') A_{s_q', s_q}(-\mathbf{k}-\mathbf{p}') \delta(\mathbf{p}'-\mathbf{p}) \delta(\mathbf{k}+\mathbf{k}''+\mathbf{p}') \quad (B.10) $$

$$ \overline{a_{s_p}^*(\mathbf{p}') a_{s_q}(\mathbf{k}'')} \overline{a_{s_q}^*(-\mathbf{k}-\mathbf{p}')} a_{s_p}(\mathbf{p}) = A_{s_p', s_q}(\mathbf{p}') A_{s_q', s_p}(-\mathbf{k}-\mathbf{p}') \delta(\mathbf{p}'-\mathbf{k}'') \delta(\mathbf{k}+\mathbf{p}+\mathbf{p}') \quad (B.11) $$

are the quasi-normal contributions to (B.8).

(B.8)-(B.11) are used in (B.7). Given the factor $\delta(\mathbf{k})$ in (B.9), the corresponding contribution to (B.7) depends on the value of $M_{s'_p s'_q}(\mathbf{k}, \mathbf{p}')$ at $\mathbf{k}=0$, which is zero according to (2.28) and (2.31). (B.10) gives

$$ \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') \sum_{s_p, s'_p, s_q, s'_q = 0, \pm 1} M_{s'_p s'_q s_q, s_p s_p}^*(\mathbf{k}, \mathbf{p}) \exp(iF_{s'_p s'_q, s_q}(\mathbf{k}, \mathbf{p})t) A_{s_p s_q}(\mathbf{p}) A_{s'_q s_p}(\mathbf{k}-\mathbf{p}) \quad (B.12) $$

while (B.11) yields

$$ \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') \sum_{s_p, s'_p, s_q, s'_q = 0, \pm 1} M_{s'_p s'_q s_q, s_p s_p}^*(\mathbf{k}, \mathbf{k}'') \exp(iF_{s'_p s'_q, s_q}(\mathbf{k}, \mathbf{k}'')t) A_{s_p s_q}(\mathbf{k}'') A_{s'_q s_p}(\mathbf{k}-\mathbf{k}'') \quad (B.13) $$

Given the Dirac function, $\mathbf{k}''$ is replaced by $-\mathbf{k}-\mathbf{p}$ inside the sum of (B.13). Using the relations $M_{s'_p s'_q}(\mathbf{k}, -\mathbf{k}-\mathbf{p}) = M_{s'_q s'_p}(\mathbf{k}, \mathbf{p})$ and $F_{s'_p s'_q}(\mathbf{k}, -\mathbf{k}-\mathbf{p}) = F_{s'_q s'_p}(\mathbf{k}, \mathbf{p})$, permutation of the summation indices $s'_p$ and $s'_q$ shows that (B.13) equals (B.12), hence the total quasi-normal contribution to (B.7) is $\delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') \Xi_{s'_p s'_q}^{QN}(\mathbf{k}, \mathbf{p})$, where

$$ \Xi_{s'_p s'_q}^{QN}(\mathbf{k}, \mathbf{p}) = 2 \sum_{s'_p, s'_q = 0, \pm 1} M_{s'_p s'_q}^*(\mathbf{k}, \mathbf{p}) \exp(iF_{s'_p s'_q}(\mathbf{k}, \mathbf{p})t) A_{s'_p s'_q}(\mathbf{p}) A_{s'_q s'_p}(-\mathbf{k}-\mathbf{p}) \quad (B.14) $$

Including the cumulant and dissipative contributions, (B.7) gives
$$ \overline{\frac{\partial a_{s'}(\mathbf{k})}{\partial t} a_{s_p}(\mathbf{p}) a_{s_q}(\mathbf{k}'')} = \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') (\Xi_{s's_p s_q}(\mathbf{k}, \mathbf{p}) - D_{s's'}(\mathbf{k}) \Theta_{s's_p s_q}(\mathbf{k}, \mathbf{p})) \quad (B.15) $$

where

$$ \Xi_{s's_p s_q}(\mathbf{k}, \mathbf{p}) = \Xi_{s's_p s_q}^{QN}(\mathbf{k}, \mathbf{p}) + \sum_{s_p', s_q' = 0, \pm 1} \int M_{s's_p's_q'}^*(\mathbf{k}, \mathbf{p}') \exp(iF_{s's_p's_q'}(\mathbf{k}, \mathbf{p}')t) K_{s's'_q s_p s_q}(\mathbf{k}, \mathbf{p}, \mathbf{p}') d^3\mathbf{p}'. \quad (B.16) $$

(B.15) provides the first term on the right-hand side of (B.6). The other two also follow from (B.15) as

$$ \overline{\frac{\partial a_{s_p}(\mathbf{p})}{\partial t} a_{s'}(\mathbf{k}) a_{s_q}(\mathbf{k}'')} = \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') (\Xi_{s's_p s_q}(\mathbf{p}, \mathbf{k}) - D_{s_p s_p}(\mathbf{p}) \Theta_{s_p s'_q}(\mathbf{p}, \mathbf{k})) \quad (B.17) $$

$$ \overline{\frac{\partial a_{s_q}(\mathbf{k}'')}{\partial t} a_{s'}(\mathbf{k}) a_{s_p}(\mathbf{p})} = \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') (\Xi_{s_q s'_q}(\mathbf{k}'', \mathbf{k}) - D_{s_q s_q}(\mathbf{k}''')) \Theta_{s'_q s'_p}(\mathbf{k}'', \mathbf{k}) \quad (B.18) $$

Given the Dirac function in (B.18),

$$ \overline{\frac{\partial a_{s_q}(\mathbf{k}'')}{\partial t} a_{s'}(\mathbf{k}) a_{s_p}(\mathbf{p})} = \delta(\mathbf{k}+\mathbf{p}+\mathbf{k}'') (\Xi_{s_q s'_q}(\mathbf{-k}-\mathbf{p}, \mathbf{k}) - D_{s_q s_q}(\mathbf{-k}-\mathbf{p}) \Theta_{s'_q s'_p}(\mathbf{-k}-\mathbf{p}, \mathbf{k})) \quad (B.19) $$

According to (B.1), $\Theta_{s_p s'_q}(\mathbf{p}, \mathbf{k}) = \Theta_{s'_q s'_p}(\mathbf{-k}-\mathbf{p}, \mathbf{k}) = \Theta_{s's_p s_q}(\mathbf{k}, \mathbf{p})$, hence (B.6), (B.15), (B.17) and (B.19) yield

$$ \overline{\frac{\partial \Theta_{s's_p s_q}(\mathbf{k}, \mathbf{p})}{\partial t} + D_{s's_p s_q}^{\Theta}(\mathbf{k}, \mathbf{p}) \Theta_{s's_p s_q}(\mathbf{k}, \mathbf{p})} = \Xi_{s's_p s_q}(\mathbf{k}, \mathbf{p}) + \Xi_{s_p s'_q}(\mathbf{p}, \mathbf{k}) + \Xi_{s'_q s'_p}(\mathbf{-k}-\mathbf{p}, \mathbf{k}) \quad (B.20) $$

where

$$ D_{s's_p s_q}^{\Theta}(\mathbf{k}, \mathbf{p}) = D_{s's'}(\mathbf{k}) + D_{s_p s_p}(\mathbf{p}) + D_{s_q s_q}(\mathbf{-k}-\mathbf{p}). \quad (B.21) $$

The solution of (B.20) is
$$ \Theta_{s's_p's_q}(\mathbf{k}, \mathbf{p}, t) = \exp(-D_{s's_p's_q}^{\Theta}(\mathbf{k}, \mathbf{p})t) \Theta_{s's_p's_q}(\mathbf{k}, \mathbf{p}, 0) + \\ \Gamma_{s's_p's_q}(\mathbf{k}, \mathbf{p}, t) + \Gamma_{s_p's's_q}(\mathbf{p}, \mathbf{k}, t) + \Gamma_{s_q's's_p}(-\mathbf{k} - \mathbf{p}, \mathbf{k}, t), \qquad (B.22) $$

where

$$ \Gamma_{s's_p's_q}(\mathbf{k}, \mathbf{p}, t) = \int_0^t \Xi_{s's_p's_q}(\mathbf{k}, \mathbf{p}, t') \exp(D_{s's_p's_q}^{\Theta}(\mathbf{k}, \mathbf{p})(t'-t)) dt' \qquad (B.23) $$

and we have used the fact that $D_{s's_p's_q}^{\Theta}(\mathbf{p}, \mathbf{k}) = D_{s's_p's_p}^{\Theta}(-\mathbf{k} - \mathbf{p}, \mathbf{k}) = D_{s's_p's_q}^{\Theta}(\mathbf{k}, \mathbf{p})$, which follows from (B.21).

Using (B.3), the contribution of the final term in (B.22) to $\tau_{ss'}(\mathbf{k})$ is

$$ \sum_{s_p, s_q = 0, \pm 1} \int M_{ss_p's_q}(\mathbf{k}, \mathbf{p}) \Gamma_{s_q's's_p}(-\mathbf{k} - \mathbf{p}, \mathbf{k}, t) \exp(-iF_{ss_p's_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p}. \qquad (B.24) $$

Employing $M_{ss_p's_q}(\mathbf{k}, \mathbf{p}) = M_{ss_q's_p}(\mathbf{k}, -\mathbf{k} - \mathbf{p})$ and $F_{ss_p's_q}(\mathbf{k}, \mathbf{p}) = F_{ss_q's_p}(\mathbf{k}, -\mathbf{k} - \mathbf{p})$, changing the integration variable to $-\mathbf{k} - \mathbf{p}$ and permuting the summation indices $s_p$ and $s_q$ shows that the contributions to $\tau_{ss'}(\mathbf{k}, t)$ of the final two terms of (B.22) are equal. Thus, (B.3) and (B.22) give

$$ \begin{aligned} \tau_{ss'}(\mathbf{k}) = & \sum_{s_p, s_q = 0, \pm 1} \int M_{ss_p's_q}(\mathbf{k}, \mathbf{p}) \Theta_{s's_p's_q}(\mathbf{k}, \mathbf{p}, 0) \\ & \exp\left(-iF_{ss_p's_q}(\mathbf{k}, \mathbf{p}) + D_{s's_p's_q}^{\Theta}(\mathbf{k}, \mathbf{p})\right)t \right) d^3\mathbf{p} + \\ & \sum_{s_p, s_q = 0, \pm 1} \int M_{ss_p's_q}(\mathbf{k}, \mathbf{p}) \Gamma_{s's_p's_q}(\mathbf{k}, \mathbf{p}, t) \exp\left(-iF_{ss_p's_q}(\mathbf{k}, \mathbf{p})\right)t \right) d^3\mathbf{p} + \\ & 2 \sum_{s_p, s_q = 0, \pm 1} \int M_{ss_p's_q}(\mathbf{k}, \mathbf{p}) \Gamma_{s's_p's_q}(\mathbf{p}, \mathbf{k}, t) \exp\left(-iF_{ss_p's_q}(\mathbf{k}, \mathbf{p})\right)t \right) d^3\mathbf{p} \end{aligned} \qquad (B.25) $$

We henceforth suppose that $s$ and $s'$ are nonzero, i.e. we specialise to the wave component, so (B.5) gives

$$ \frac{\partial A_{ss'}(\mathbf{k})}{\partial t} + 2D(\mathbf{k}) A_{ss'}(\mathbf{k}) = \tau_{ss'}(\mathbf{k}) + \tau_{s's}(\mathbf{k}), \qquad (B.26) $$

where $D(\mathbf{k})$, given by (3.14), is the damping factor of wave modes. Up to now, the only approximation used was the replacement of the dissipative term in (2.26) to obtain (3.12). From here on, we use large $t$ to derive the wave-turbulence equations and, to avoid killing the turbulence by visco-diffusion prior to the given time, the modal damping coefficient is assumed comparable to or smaller than $t^{-1}$ for the large scales of turbulence. Furthermore, to obtain closed wave-turbulence equations from (B.26), we suppose the amplitude, $\epsilon_{NP}$, of the
NP component small compared with that of the wave component, $\epsilon$, and consider the effect of each of the terms on the right-hand side of (B.25) on the long-term evolution of $A_{ss'}$ according to (B.26).

Given large $t$ and (2.29), the integrand of the first term on the right-hand side of (B.25) has rapid, self-cancelling oscillations with respect to $\mathbf{p}$ unless $s_p = s_q = 0$. Neglecting nonzero $s_p$ and $s_q$ leaves

$$ \exp(-is\omega(\mathbf{k})t) \int M_{s00}(\mathbf{k}, \mathbf{p}) \Theta_{s'00}(\mathbf{k}, \mathbf{p}, 0) \exp(-D_{s'00}^{\Theta}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p} \quad (B.27) $$

as the contribution to $\tau_{ss'}$. Since $s \neq 0$, (B.27) has temporal oscillations of period $O(1)$ and hence does not contribute significantly to the long-term evolution of $A_{ss'}$ according to (B.26). Thus, the first term in (B.25) is negligible.

Concerning the cumulant contribution, the integrand in (B.16) has rapid, self-cancelling oscillations with respect to $\mathbf{p}'$ unless $s'_p = s'_q = 0$, hence other terms are dropped. Using (B.23), the resulting contribution to (B.25) is

$$ \begin{aligned} & \sum_{s_p, s_q = 0, \pm 1} \int_0^t \int_0^t M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) M_{s'00}^*(\mathbf{k}, \mathbf{p}') K_{00s_p s_q}(\mathbf{k}, \mathbf{p}, \mathbf{p}', t') \\ & \exp(is'\omega(\mathbf{k})t' + D_{s's_p s_q}^{\Theta}(\mathbf{k}, \mathbf{p})(t'-t) - iF_{ss_p s_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p} d^3\mathbf{p}' dt' + \\ & 2 \sum_{s_p, s_q = 0, \pm 1} \int_0^t \int_0^t M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) M_{s_p 00}^*(\mathbf{p}, \mathbf{p}') K_{00s's_q}(\mathbf{p}, \mathbf{k}, \mathbf{p}', t') \\ & \exp(is_p \omega(\mathbf{p})t' + D_{s_p s_q}^{\Theta}(\mathbf{p}, \mathbf{k})(t'-t) - iF_{ss_p s_q}(\mathbf{k}, \mathbf{p})t) d^3\mathbf{p} d^3\mathbf{p}' dt' \end{aligned} \quad (B.28) $$

Avoidance of rapid oscillations with respect to $\mathbf{p}$ leads to

$$ \begin{aligned} & \exp(-is\omega(\mathbf{k})t) \int_0^t \int_0^t M_{s00}(\mathbf{k}, \mathbf{p}) M_{s'00}^*(\mathbf{k}, \mathbf{p}') K_{0000}(\mathbf{k}, \mathbf{p}, \mathbf{p}', t') \\ & \exp(is'\omega(\mathbf{k})t' + D_{s'00}^{\Theta}(\mathbf{k}, \mathbf{p})(t'-t)) d^3\mathbf{p} d^3\mathbf{p}' dt' + \\ & 2 \exp(-is\omega(\mathbf{k})t) \sum_{s_p=0, \pm 1} \int_0^t \int_0^t M_{ss_p 0}(\mathbf{k}, \mathbf{p}) M_{s_p 00}^*(\mathbf{p}, \mathbf{p}') K_{00s'0}(\mathbf{p}, \mathbf{k}, \mathbf{p}', t') \\ & \exp(is_p \omega(\mathbf{p})(t'-t) + D_{s_p s'0}^{\Theta}(\mathbf{p}, \mathbf{k})(t'-t)) d^3\mathbf{p} d^3\mathbf{p}' dt' \end{aligned} \quad (B.29) $$

Since $s' \neq 0$, the integrand of the first term of (B.29) is oscillatory with respect to $t'$, so that term does not grow with $t$ and has order of magnitude determined by $K_{0000} = O(\epsilon_{NP}^4)$. The resulting contribution to (B.26) can, at most, induce a change of $O(\epsilon_{NP}^4 t)$ in $A_{ss'}$ over the time span $t$. However, as we shall see, the appropriate time scale for evolution of $A_{ss'}$ is $t \sim \epsilon^{-2}$,
giving a change in $A_{ss'}$ due to the first term in (B.29) of at most $O(\epsilon_{NP}^4\epsilon^{-2})$, which is negligible since $A_{ss'} = O(\epsilon^2)$ and $\epsilon_{NP} \ll \epsilon$.

The $s_p \neq 0$ components of the second term in (B.29) have integrands which have rapid oscillations with respect to $\mathbf{p}$ unless $t-t' = O(1)$. Neglecting other $t'$, the order of magnitude of the $s_p \neq 0$ contributions is determined by $K_{00s'0} = O(\epsilon_{NP}^3\epsilon)$, implying negligible effects on the evolution of $A_{ss'}$ when $t \sim \epsilon^{-2}$ because $\epsilon_{NP} \ll \epsilon$. Finally, the $s_p = 0$ component of the second term in (B.29) has oscillations of period $O(1)$ due to the term $\exp(-is\omega(\mathbf{k})t)$. The amplitude of these oscillations is $O(\epsilon_{NP}^3\epsilon t)$. The solution of (B.26) inherits these oscillations, but the result is negligible when $t \sim \epsilon^{-2}$ because $\epsilon_{NP} \ll \epsilon$. In summary, thanks to $\epsilon_{NP} \ll \epsilon$, the cumulant contribution is negligible, but there remains the quasi-normal one, which we now investigate.

Employing (B.14) and (B.23), the quasi-normal contribution to the second term on the right-hand side of (B.25) is

$$2 \sum_{\substack{s_p, s_q = 0, \pm 1 \\ s'_p, s'_q = 0, \pm 1}} \int_0^t \int M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) M_{s's'_p s'_q}^*(\mathbf{k}, \mathbf{p}) \\ \exp\left(i(F_{s's'_p s'_q}(\mathbf{k}, \mathbf{p})t' - F_{ss_p s_q}(\mathbf{k}, \mathbf{p})t) + D_{s's_p s_q}^{\Theta}(\mathbf{k}, \mathbf{p})(t'-t)\right). \quad (B.30) \\ A_{s'_p s_p}(\mathbf{p}, t') A_{s'_q s_q}(-\mathbf{k}-\mathbf{p}, t') d^3\mathbf{p} dt'$$

Given

$$F_{s's'_p s'_q}(\mathbf{k}, \mathbf{p})t' - F_{ss_p s_q}(\mathbf{k}, \mathbf{p})t = \\ (s't' - st)\omega(\mathbf{k}) + (s'_p t' - s_p t)\omega(\mathbf{p}) + (s'_q t' - s_q t)\omega(-\mathbf{k}-\mathbf{p}), \quad (B.31)$$

there are rapid, self-cancelling oscillations with respect to $\mathbf{p}$ of the integrand in (B.30) unless $s'_p = s_p, s'_q = s_q$. The case $s'_p = s'_q = s_p = s_q = 0$ gives the contribution

$$2 \exp(-is\omega(\mathbf{k})t) \int_0^t \int M_{s00}(\mathbf{k}, \mathbf{p}) M_{s'00}^*(\mathbf{k}, \mathbf{p}) A_{00}(\mathbf{p}, t') A_{00}(-\mathbf{k}-\mathbf{p}, t') \\ \exp(is'\omega(\mathbf{k})t' + D_{s'00}^{\Theta}(\mathbf{k}, \mathbf{p})(t'-t)) d^3\mathbf{p} dt' \quad (B.32)$$

Since $s' \neq 0$, the integrand is oscillatory with respect to $t'$ and, given $A_{00} = O(\epsilon_{NP}^2)$, (B.32) is of order $O(\epsilon_{NP}^4)$, hence inducing negligible effects on the evolution of $A_{ss'}$ according to (B.26). On the other hand, when $s'_p = s_p, s'_q = s_q$ and one or other of $s_p$ and $s_q$ is nonzero,
there are rapid oscillations with respect to **p** of the integrand of (B.30) unless $t-t'=O(1)$. The latter condition has two consequences. The first is that, given the small visco-diffusive dissipation needed to stop it killing the turbulence before nonlinearity is effective, the visco-diffusive term in the exponential of (B.30) can be neglected. The second follows from the slow evolution of the spectral matrix A, which allows the approximation $A(t') = A(t)$. Thus, (B.30) becomes

$$2 \sum_{s_p, s_q = 0, \pm 1} \int_0^t \int_0^t M_{ss_{p}s_q}(\mathbf{k}, \mathbf{p}) M_{s'_{p}s_q}^*(\mathbf{k}, \mathbf{p}) \\ \exp\left(i(F_{s'_{p}s_q}(\mathbf{k}, \mathbf{p})t' - F_{ss_{p}s_q}(\mathbf{k}, \mathbf{p})t)\right) A_{s_p s_p}(\mathbf{p}, t) A_{s_q s_q}(-\mathbf{k} - \mathbf{p}, t) d^3\mathbf{p} dt', \quad (B.33)$$

where the sum excludes the term $s_p = s_q = 0$. The integral over $t'$ can now be evaluated giving

$$2 \exp(i(s'-s)\omega(\mathbf{k})t) \sum_{s_p, s_q = 0, \pm 1} \int \Delta_{s'_{p}s_q}(\mathbf{k}, \mathbf{p}, t) M_{ss_{p}s_q}(\mathbf{k}, \mathbf{p}) M_{s'_{p}s_q}^*(\mathbf{k}, \mathbf{p}) \quad (B.34) \\ A_{s_p s_p}(\mathbf{p}, t) A_{s_q s_q}(-\mathbf{k} - \mathbf{p}, t) d^3\mathbf{p}$$

for the second-term on the right-hand side of (B.25), where again $s_p = s_q = 0$ is excluded from the sum,

$$\Delta_{s'_{p}s_q}(\mathbf{k}, \mathbf{p}, t) = \frac{1 - \exp(-iF_{s'_{p}s_q}(\mathbf{k}, \mathbf{p})t)}{iF_{s'_{p}s_q}(\mathbf{k}, \mathbf{p})} \quad (B.35)$$

and $F_{ss_{p}s_q}(\mathbf{k}, \mathbf{p}) = F_{s'_{p}s_q}(\mathbf{k}, \mathbf{p}) + (s-s')\omega(\mathbf{k})$ has been used.

Turning attention to the quasi-normal contribution to the final term of (B.25), we employ similar reasoning to that used above. The equivalent of (B.30) is

$$4 \sum_{\substack{s_p, s_q = 0, \pm 1 \\ s'_p, s'_q = 0, \pm 1}} \int_0^t \int_0^t M_{ss_{p}s_q}(\mathbf{k}, \mathbf{p}) M_{s'_{p}s'_q}^*(\mathbf{p}, \mathbf{k}) \\ \exp\left(i(F_{s'_{p}s_p s'_q}(\mathbf{k}, \mathbf{p})t' - F_{ss_{p}s_q}(\mathbf{k}, \mathbf{p})t) + D_{s_p s'_q}^{\oplus}(\mathbf{p}, \mathbf{k})(t'-t)\right), \quad (B.36) \\ A_{s_p s'_q}(\mathbf{k}, t') A_{s'_q s_q}(-\mathbf{k} - \mathbf{p}, t') d^3\mathbf{p} dt'$$

where we have used $F_{s_p s'_p s'_q}(\mathbf{p}, \mathbf{k}) = F_{s'_p s_p s'_q}(\mathbf{k}, \mathbf{p})$. Avoidance of rapid oscillations with respect to **p** requires $s'_q = s_q$. The case $s_p = s'_q = s_q = 0$ gives the contribution
$$4\exp(-is\omega(\mathbf{k})t)$$

$$\sum_{s'_p = 0, \pm 1} \int_0^t \int_{00} M_{s00}(\mathbf{k}, \mathbf{p}) M_{0s'_p 0}^*(\mathbf{p}, \mathbf{k}) A_{s'_p s'}(\mathbf{k}, t') A_{00}(-\mathbf{k}-\mathbf{p}, t'). \quad (B.37)$$
$$\exp(is'_p \omega(\mathbf{k})t' + D_{0s'_p 0}^{\Theta}(\mathbf{p}, \mathbf{k})(t'-t)) d^3\mathbf{p} dt'$$

When $s'_p \neq 0$, the integrand is oscillatory and the result has order of magnitude $\epsilon_{NP}^2 \epsilon^2$, hence negligible effects on the evolution of $A_{ss'}$ for $t \sim \epsilon^{-2}$. For $s'_p = 0$, there are oscillations of period $O(1)$ and amplitude $O(\epsilon_{NP}^3 \epsilon t)$. These oscillations are inherited by $A_{ss'}$ via (B.26), but are negligible because $\epsilon_{NP} \ll \epsilon$.

Finally, when $s'_q = s_q$ and one or other of $s_p$ and $s_q$ is nonzero, there are rapid oscillations with respect to $\mathbf{p}$ of the integrand of (B.36) unless $t-t' = O(1)$, a condition which makes the visco-diffusive term in (B.36) negligible and allows the approximation $A(t') = A(t)$. Evaluating the integral over $t'$ and using $F_{ss_p s_q}(\mathbf{k}, \mathbf{p}) = F_{s'_p s'_p s_q}(\mathbf{k}, \mathbf{p}) + (s - s'_p)\omega(\mathbf{k})$, gives

$$4 \sum_{s_p, s_q, s'_p = 0, \pm 1} \exp(i(s'_p - s)\omega(\mathbf{k})t) A_{s'_p s'}(\mathbf{k}, t) \int \Delta_{s'_p s'_p s_q}(\mathbf{k}, \mathbf{p}, t) M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) M_{s'_p s'_p s_q}^*(\mathbf{p}, \mathbf{k}) A_{s_q s_q}(-\mathbf{k}-\mathbf{p}, t) d^3\mathbf{p}, \quad (B.38)$$

where the sum excludes $s_p = s_q = 0$. $\Delta_{s'_p s'_p s_q}(\mathbf{k}, \mathbf{p}, t) = \Delta_{s'_p s'_p s_q}(\mathbf{k}, -\mathbf{k}-\mathbf{p}, t)$ (which follows from (B.35) and $F_{s'_p s'_p s_q}(\mathbf{k}, \mathbf{p}) = F_{s'_p s'_p (-\mathbf{k}-\mathbf{p})}(\mathbf{k}, \mathbf{p})$) and $M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) = M_{s'_p (-\mathbf{k}-\mathbf{p})}(\mathbf{k}, \mathbf{p})$ are then used, followed by a change of integration variable to $-\mathbf{k}-\mathbf{p}$ and permutation of the summation indices $s_p$ and $s_q$. Thus,

$$4 \sum_{s_p, s_q, s'_p = 0, \pm 1} \exp(i(s'_p - s)\omega(\mathbf{k})t) A_{s'_p s'}(\mathbf{k}, t) \int \Delta_{s'_p s'_p s_q}(\mathbf{k}, \mathbf{p}, t) M_{ss_p s_q}(\mathbf{k}, \mathbf{p}) M_{s'_p s'_p s_q}^*(-\mathbf{k}-\mathbf{p}, \mathbf{k}) A_{s_q s_q}(\mathbf{p}, t) d^3\mathbf{p} \quad (B.39)$$

expresses the final term in (B.25), where again $s_p = s_q = 0$ is excluded from the sum.

The integrals in (B.34) and (B.39) have the form

$$\int \Delta_{s'_p s'_p s_q}(\mathbf{k}, \mathbf{p}, t) \Phi(\mathbf{p}) d^3\mathbf{p}, \quad (B.40)$$

where $s''$ and $\Phi$ depend on which integral is considered. Such integrals always arise in wave-turbulence analysis and need to be evaluated in the limit of large $t$. To this end, (B.35) is decomposed into real and imaginary parts:
$$ \Delta_{s^*s_p s_q}(\mathbf{k}, \mathbf{p}, t) = \frac{\sin(F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})t)}{F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})} - i \frac{1 - \cos(F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})t)}{F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})}. \quad (B.41) $$

Without going into the details (see e.g. Benney and Newell (1969)), the first term approaches a Dirac function, $\pi\delta(F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p}))$, as $t \to \infty$. Using $\omega(-\mathbf{k}-\mathbf{p}) = \omega(\mathbf{k}+\mathbf{p})$, (2.29) and (3.15), $\nabla_{\mathbf{p}} F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p}) = s_p c_g(\mathbf{p}) + s_q c_g(\mathbf{k}+\mathbf{p})$. Thus, converting the volume integral with the Dirac function into a surface integral, the contribution to (B.40) is

$$ \pi \int_{S_{s^*s_p s_q}(\mathbf{k})} \frac{\Phi(\mathbf{p})}{|s_p c_g(\mathbf{p}) + s_q c_g(\mathbf{k}+\mathbf{p})|} d^3\mathbf{p}, \quad (B.42) $$

where $S_{s^*s_p s_q}(\mathbf{k})$ is the surface in $\mathbf{p}$-space defined by $F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})=0$. Of course, $S_{s^*s_p s_q}(\mathbf{k})$ may be the empty set, in which case integrals over $S_{s^*s_p s_q}(\mathbf{k})$ should be interpreted as zero. The cosine in the second term of (B.41) turns out to be negligible provided that the result is interpreted as a Cauchy principal value, so the contribution to (B.40) is

$$ -iP \int \frac{\Phi(\mathbf{p})}{F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p}, \quad (B.43) $$

where the $P$ before the integral sign indicates a principal value. Finally, the sum of (B.42) and (B.43) yields

$$ \int \Delta_{s^*s_p s_q}(\mathbf{k}, \mathbf{p}, t) \Phi(\mathbf{p}) d^3\mathbf{p} \rightarrow \pi \int_{S_{s^*s_p s_q}(\mathbf{k})} \frac{\Phi(\mathbf{p})}{|s_p c_g(\mathbf{p}) + s_q c_g(\mathbf{k}+\mathbf{p})|} d^3\mathbf{p} - iP \int \frac{\Phi(\mathbf{p})}{F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p} \quad (B.44) $$

as the $t \to \infty$ limit of (B.40).

The Cauchy principal value has been introduced above without defining what is meant. If there is no $\mathbf{p}$ for which $F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})=0$ (i.e. $S_{s^*s_p s_q}(\mathbf{k})$ is the empty set), the integrand in (B.43) is non-singular and the principal value is just a normal integral. On the other hand, $F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})=0$ gives an infinite singularity of the integrand and the integral needs more careful interpretation. To remove the singularity, a small region in $\mathbf{p}$-space is excluded from the integral (see figure B.1). This region can be defined by $|F_{s^*s_p s_q}(\mathbf{k}, \mathbf{p})|<\delta$, where $\delta>0$ is small. Taking the limit $\delta \searrow 0$ gives the principal value.
<figure><img src="image_27.png" /></figure>

Figure B.1: Illustration of the volume in p-space which is excluded from the integral in (B.43) when defining the Cauchy principal value. The excluded volume is bounded by the dashed lines.

Using (B.44) in (B.34) and (B.39), there are two types of terms: those which are oscillatory over time spans of $O(1)$ and those which are not. Given that we are looking for contributions that result in cumulative evolution of $A_{ss'}$ according to (B.26) over long time spans, we neglect oscillatory terms, i.e. $s' \neq s$ for (B.34) and $s'_p \neq s$ for (B.39). Taking the sum of (B.34) and (B.39) without these terms,

$$ \tau_{ss'}(\mathbf{k}) = \sum_{s_p, s_q = 0, \pm 1} \int \Delta_{ss_p s_q}(\mathbf{k}, \mathbf{p}, t) \Phi_{ss' p s'_q}(\mathbf{p}, \mathbf{k}) d^3\mathbf{p}, \quad (B.45) $$

where the sum is missing the term $s_p = s_q = 0$,

$$ \Phi_{ss's_p s_q}(\mathbf{p}, \mathbf{k}) = A_{s_p s_p}(\mathbf{p}) (\delta_{ss'} \lambda_{ss_p s_q}(\mathbf{k}, \mathbf{p}) A_{s_q s_q}(-\mathbf{k} - \mathbf{p}) + \zeta_{ss_p s_q}(\mathbf{k}, \mathbf{p}) A_{ss'}(\mathbf{k})) \quad (B.46) $$

and $\lambda_{ss_p s_q}(\mathbf{k}, \mathbf{p})$ and $\zeta_{ss_p s_q}(\mathbf{k}, \mathbf{p})$ are given by (3.16) and (3.17). Using (B.44), (B.45) and (B.46) imply

$$ \tau_{ss'}(\mathbf{k}) = \mu_s(\mathbf{k}) A_{ss'}(\mathbf{k}) + \delta_{ss'} \sum_{s_p, s_q = 0, \pm 1} \left( \pi \int_{S_{ss'p'sq}} \frac{\lambda_{ss_p s_q}(\mathbf{k}, \mathbf{p}) A_{s_p s_p}(\mathbf{p}) A_{s_q s_q}(-\mathbf{k} - \mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k} + \mathbf{p})|} d^3\mathbf{p} - iP \int \frac{\lambda_{ss_p s_q}(\mathbf{k}, \mathbf{p}) A_{s_p s_p}(\mathbf{p}) A_{s_q s_q}(-\mathbf{k} - \mathbf{p})}{F_{ss_p s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p} \right) \quad (B.47) $$

where
$$ \mu_s(\mathbf{k}) = \sum_{s_p, s_q = 0, \pm 1} \left( \pi \int_{S_{ss,p,q}(\mathbf{k})} \frac{\zeta_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k} + \mathbf{p})|} d^2\mathbf{p} - iP \int \frac{\zeta_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p})}{F_{ss,p,s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p} \right) \quad (B.48) $$

and the term $s_p = s_q = 0$ is excluded from the sums in (B.47) and (B.48). Given the assumed smallness of the NP component, the contributions of $A_{00}$ to (B.47) and (B.48) are negligible, thus

$$ \begin{aligned} \tau_{ss'}(\mathbf{k}) ={}& \mu_s(\mathbf{k}) A_{ss'}(\mathbf{k}) + \\ & \delta_{ss'} \sum_{s_p, s_q = \pm 1} \left( \pi \int_{S_{ss,p,q}(\mathbf{k})} \frac{\lambda_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p}) A_{s_q, s_q}(-\mathbf{k} - \mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k} + \mathbf{p})|} d^2\mathbf{p} - iP \int \frac{\lambda_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p}) A_{s_q, s_q}(-\mathbf{k} - \mathbf{p})}{F_{ss,p,s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p} \right) \end{aligned} \quad (B.49) $$

$$ \mu_s(\mathbf{k}) = \sum_{\substack{s_p=\pm 1 \\ s_q=0, \pm 1}} \left( \pi \int_{S_{ss,p,q}(\mathbf{k})} \frac{\zeta_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k} + \mathbf{p})|} d^2\mathbf{p} - iP \int \frac{\zeta_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p})}{F_{ss,p,s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p} \right) \quad (B.50) $$

Both integrals in (B.49) are real, hence, given the Hermitian character of $A_{ss'}$ and the factor $\delta_{ss'}$,

$$ \begin{aligned} \tau_{ss'}(\mathbf{k}) + \tau_{s's}^*(\mathbf{k}) ={}& (\mu_s(\mathbf{k}) + \mu_{s'}^*(\mathbf{k})) A_{ss'}(\mathbf{k}) + \\ & 2\pi \delta_{ss'} \sum_{s_p, s_q = \pm 1} \int_{S_{ss,p,q}(\mathbf{k})} \frac{\lambda_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p}) A_{s_q, s_q}(-\mathbf{k} - \mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k} + \mathbf{p})|} d^2\mathbf{p} \end{aligned} \quad (B.51) $$

is the final asymptotic approximation of the right-hand side of (B.26), giving the wave-turbulence equation

$$ \begin{aligned} & \frac{\partial A_{ss'}(\mathbf{k})}{\partial t} + 2D(\mathbf{k}) A_{ss'}(\mathbf{k}) = (\mu_s(\mathbf{k}) + \mu_{s'}^*(\mathbf{k})) A_{ss'}(\mathbf{k}) + \\ & 2\pi \delta_{ss'} \sum_{s_p, s_q = \pm 1} \int_{S_{ss,p,q}(\mathbf{k})} \frac{\lambda_{ss,p,s_q}(\mathbf{k}, \mathbf{p}) A_{s_p, s_p}(\mathbf{p}) A_{s_q, s_q}(-\mathbf{k} - \mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k} + \mathbf{p})|} d^2\mathbf{p} \end{aligned} \quad (B.52) $$
for $s, s' \neq 0$. Given $A_{ss'} = O(\epsilon^2)$, the evolution time of the wave component can be estimated using (B.50) and (B.52) as $O(\epsilon^{-2})$, a result noted earlier and which is typical of wave turbulence.

Since the wave-component spectral energy density is given by the first of the equations (2.40), the most interesting application of (B.52) is $s' = s$. Using (B.50) and reality of the diagonal components of the spectral matrix $A$, (B.52) yields

$$ \frac{\partial A_{ss}(\mathbf{k})}{\partial t} + 2D(\mathbf{k}) A_{ss}(\mathbf{k}) = 2\pi \sum_{\substack{s_p=\pm 1 \\ s_q=0, \pm 1}} \int_{S_{ss_p s_q}(\mathbf{k})} \frac{A_{s_p s_p}(\mathbf{p})}{|s_p \mathbf{c}_g(\mathbf{p}) + s_q \mathbf{c}_g(\mathbf{k}+\mathbf{p})|} \left( (1-\delta_{s_q 0}) \lambda_{ss_p s_q}(\mathbf{k}, \mathbf{p}) A_{s_q s_q}(-\mathbf{k}-\mathbf{p}) + \Re(\zeta_{ss_p s_q}(\mathbf{k}, \mathbf{p})) A_{ss}(\mathbf{k}) \right) d^3\mathbf{p} + 2A_{ss}(\mathbf{k}) \sum_{\substack{s_p=\pm 1 \\ s_q=0, \pm 1}} P \int \frac{\Im(\zeta_{ss_p s_q}(\mathbf{k}, \mathbf{p})) A_{s_p s_p}(\mathbf{p})}{F_{ss_p s_q}(\mathbf{k}, \mathbf{p})} d^3\mathbf{p} \quad (B.53) $$

Consider the principal-value contribution to (B.53), which can be written

$$ 2A_{ss}(\mathbf{k}) \sum_{s_p=\pm 1} P \int \sum_{s_q=0, \pm 1} \frac{\Im(\zeta_{ss_p s_q}(\mathbf{k}, \mathbf{p}))}{F_{ss_p s_q}(\mathbf{k}, \mathbf{p})} A_{s_p s_p}(\mathbf{p}) d^3\mathbf{p}. \quad (B.54) $$

As noted following its definition, (2.35), the spectral matrix is such that $A_{ss'}(-\mathbf{k}) = A_{-s', -s}(\mathbf{k})$, hence, switching the sign of the summation index $s_p$ and changing the integration variable to $\mathbf{q} = -\mathbf{p}$, (B.54) becomes

$$ 2A_{ss}(\mathbf{k}) \sum_{s_p=\pm 1} P \int \sum_{s_q=0, \pm 1} \frac{\Im(\zeta_{s_p s_q}(\mathbf{k}, -\mathbf{q}))}{F_{s_p s_q}(\mathbf{k}, -\mathbf{q})} A_{s_p s_q}(\mathbf{q}) d^3\mathbf{q}. \quad (B.55) $$

Since (B.54) and (B.55) are equal, they are also equal to half their sum and hence to

$$ A_{ss}(\mathbf{k}) \sum_{s_p=\pm 1} P \int \sum_{s_q=0, \pm 1} \left( \frac{\Im(\zeta_{ss_p s_q}(\mathbf{k}, \mathbf{p}))}{F_{ss_p s_q}(\mathbf{k}, \mathbf{p})} + \frac{\Im(\zeta_{s_p s_q}(\mathbf{k}, -\mathbf{p}))}{F_{s_p s_q}(\mathbf{k}, -\mathbf{p})} \right) A_{s_p s_p}(\mathbf{p}) d^3\mathbf{p}. \quad (B.56) $$

Although we have been unable to show it analytically, numerical evaluation of
$$ \sum_{s_q=0,\pm 1} \left( \frac{\Im(\zeta_{ss_p,s_q}(\mathbf{k},\mathbf{p}))}{F_{ss_p,s_q}(\mathbf{k},\mathbf{p})} + \frac{\Im(\zeta_{s,-s_p,s_q}(\mathbf{k},-\mathbf{p}))}{F_{s,-s_p,s_q}(\mathbf{k},-\mathbf{p})} \right) \qquad (B.57) $$

for nonzero $s$ and $s_p$ and different values of $\beta$, $\mathbf{k}$ and $\mathbf{p}$ show it to be zero to IEEE double precision. Assuming this result is exactly true, the principal-value term in (B.53) is zero and is therefore dropped.

Next consider the $s_q = 0$ contribution to the first term on the right-hand side of (B.53), which is

$$ 2\pi A_{ss}(\mathbf{k}) \sum_{s_p=\pm 1} \int_{S_{s,s_p,0}(\mathbf{k})} \frac{\Re(\zeta_{s,s_p,0}(\mathbf{k},\mathbf{p})) A_{s_p,s_p}(\mathbf{p})}{|\mathbf{c}_g(\mathbf{p})|} d^2\mathbf{p}. \qquad (B.58) $$

If $s_p = s$, $F_{s,s_p,0}(\mathbf{k},\mathbf{p}) = s(\omega(\mathbf{k}) + \omega(\mathbf{p})) = 0$ has no solution because $\omega > 0$, hence the contribution to (B.58) is zero. On the other hand, if $s_p = -s$, $S_{s,s_p,0}(\mathbf{k})$ consists of the double cone $\omega(\mathbf{p}) = \omega(\mathbf{k})$. Once again, although we have been unable to prove it analytically, numerical calculations show that $\zeta_{s,-s,0}(\mathbf{k},\mathbf{p})$ is zero to IEEE double precision when $s \neq 0$ and $\omega(\mathbf{p}) = \omega(\mathbf{k})$. Assuming this is exactly true, (B.58) is zero and (B.53) becomes (3.13) when the scaled variables $T = \varepsilon^2 t$ and $\hat{A}_{ss} = A_{ss} / \varepsilon^2$ are used.

Given $\varepsilon_{NP} \ll \varepsilon$, it might have been tempting to neglect the effects of the NP component on the wave one from the start, i.e. to drop the terms in (3.12) with $s_p = 0$ or $s_q = 0$. However, it is then found that the principal-value terms do not cancel, i.e. they persist in the final result. This is because cancellation requires the $s_q = 0$ contribution to (B.57).

## Appendix C: Existence of the resonant surface

This appendix derives conditions for the existence of solutions of $F_{ss_p,s_q}(\mathbf{k},\mathbf{p}) = 0$, where $s$, $s_p$ and $s_q$ take one of the values $\pm 1$, for a given $\mathbf{k}$. To this end, we determine the maximum and minimum values of $F_{ss_p,s_q}(\mathbf{k},\mathbf{p})$ as a function of $\mathbf{p}$. If these values straddle zero the resonant surface exists. Since $\mathbf{k}$ is fixed, looking for extrema of $F_{ss_p,s_q}(\mathbf{k},\mathbf{p})$, given by (2.29), is equivalent to searching for those of

$$ s_p s_q \omega(\mathbf{p}) + \omega(\mathbf{k} + \mathbf{p}), \qquad (C.1) $$

where we have used $\omega(-\mathbf{k} - \mathbf{p}) = \omega(\mathbf{k} + \mathbf{p})$.
Define spherical polar coordinates, $p = |\mathbf{p}|$, $0 \le \theta \le \pi$ and $-\pi < \psi \le \pi$, in **p**-space, where $\theta$ is the angle between the vectors **p** and **e** and **k** lies in the plane $\psi = 0$. We first look for extrema of (C.1) at constant **p** and $\theta$, i.e. only $\psi$ is varied. Given fixed $\theta$, $\omega(\mathbf{p})$ is constant, so we consider extrema of $\omega(\mathbf{k}+\mathbf{p})$. Using (2.9),

$$ 
\begin{aligned}
\omega(\mathbf{k}+\mathbf{p}) &= \left( \frac{N^2 q_{\perp}^2 + \Omega^2 q_{\parallel}^2}{q_{\perp}^2 + q_{\parallel}^2} \right)^{1/2} \\
&= \left( N^2 + (\Omega^2 - N^2) \frac{q_{\parallel}^2}{q_{\perp}^2 + q_{\parallel}^2} \right)^{1/2},
\end{aligned} 
\qquad (C.2) $$

where $\mathbf{q} = \mathbf{k} + \mathbf{p}$. Employing the spherical coordinates defined above,

$$ q_{\parallel} = k_{\parallel} + p \cos \theta, \qquad (C.3) $$

$$ 
\begin{aligned}
q_{\perp}^2 &= (k_{\perp} + p \sin \theta \cos \psi)^2 + p^2 \sin^2 \theta \sin^2 \psi \\
&= k_{\perp}^2 + p^2 \sin^2 \theta + 2k_{\perp} p \sin \theta \cos \psi.
\end{aligned} 
\qquad (C.4) $$

Since $p$ and $\theta$ are constant, so is $q_{\parallel}$, while the factor $\cos \psi$ induces variations of $q_{\perp}^2$. It is apparent from (C.2)-(C.4) that the extrema we are looking for arise from those of $\cos \psi$, which occur at $\psi = 0$ and $\psi = \pi$. Thus, we can restrict attention to **p** lying in the plane defined by **k** and the rotation axis (see figure C.1).

<figure><img src="image_28.png" /><figcaption>Figure C.1: The plane defined by **k** and the rotation axis, showing the vectors **e**, **p**, $-\mathbf{k}$ and $\mathbf{k}+\mathbf{p}$, as well as the angle $\theta_{\mathbf{k}+\mathbf{p}}$ and the line $L$.</figcaption></figure>
Consider the effect of letting **p** move along the diagonal dashed line, *L*, in figure E.1. This keeps $\omega(\mathbf{p})$ constant in (C.1). As a result, we look for extrema of

$$ \omega(\mathbf{k}+\mathbf{p}) = (N^2 \cos^2 \theta_{\mathbf{k}+\mathbf{p}} + \Omega^2 \sin^2 \theta_{\mathbf{k}+\mathbf{p}})^{1/2}, \quad (C.5) $$

where, as illustrated by the figure, $\theta_{\mathbf{k}+\mathbf{p}}$ is the angle between the vectors $\mathbf{k}+\mathbf{p}$ and $\mathbf{e}$. As $\mathbf{p}$ moves along the line *L*, $\theta_{\mathbf{k}+\mathbf{p}}$ varies and $\omega(\mathbf{k}+\mathbf{p})$, given by (C.5), takes on all values between $\Omega$ and $N$. Thus, the minimum and maximum values of (C.1) for the line *L* are

$$ s_p s_q \omega(\mathbf{p}) + \min(\Omega, N), \quad s_p s_q \omega(\mathbf{p}) + \max(\Omega, N). \quad (C.6) $$

We next consider changes of the line *L*, first for the case $s_p = s_q$. Thus, (C.6) gives

$$ \omega(\mathbf{p}) + \min(\Omega, N), \quad \omega(\mathbf{p}) + \max(\Omega, N). \quad (C.7) $$

As the line changes, $\omega(\mathbf{p})$ takes on every value between $\Omega$ and $N$, hence

$$ 2 \min(\Omega, N), \quad 2 \max(\Omega, N) \quad (C.8) $$

are the minimum and maximum of (C.1) over all $\mathbf{p}$. It follows from (2.29) that the extremal values of $F_{ss_p s_q}(\mathbf{k}, \mathbf{p})$ are

$$ s\omega(\mathbf{k}) + 2s_p \min(\Omega, N), \quad s\omega(\mathbf{k}) + 2s_p \max(\Omega, N). \quad (C.9) $$

As noted earlier, the existence of the resonant surface requires that the two values in (C.9) straddle zero. If $s_p = s$, this is impossible because each of $\omega(\mathbf{k})$, $\Omega$ and $N$ is positive. This means that the resonant surface does not exist when $s = s_p = s_q$, a result which is to be expected given (2.29) and $\omega > 0$. On the other hand, when $s_p = -s$, zero straddling of (C.9) requires

$$ 2 \min(\Omega, N) < \omega(\mathbf{k}) < 2 \max(\Omega, N). \quad (C.10) $$

Given that $\omega(\mathbf{k}) \le \max(\Omega, N)$ according to (2.29), the second inequality in (C.10) is automatically satisfied, while the first requires $2\min(\Omega, N) < \max(\Omega, N)$. The latter inequality is only satisfied if either $\beta < 1/2$ or $\beta > 2$. For such values of $\beta$, the resonant surface with $s_p = s_q \ne s$ exists provided (3.20) holds.

Turning attention to the case, $s_p = -s_q$, (C.6) yields
$$-\omega(\mathbf{p}) + \min(\Omega, N), \quad -\omega(\mathbf{p}) + \max(\Omega, N) \qquad (C.11)$$

as the minimum and maximum values of (C.1) for the line $L$. As the line, $L$, changes, $\omega(\mathbf{p})$ takes on all values between $\Omega$ and $N$, hence $\pm|\Omega - N|$ for the global extremal values of (C.1). It follows from (2.29) that the extremal values of $F_{ss_p s_q}(\mathbf{k}, \mathbf{p})$ are $s\omega(\mathbf{k}) \pm s_p|\Omega - N|$. The requirement that these values straddle zero is $\omega(\mathbf{k}) < |\Omega - N|$. Since $\omega(\mathbf{k}) \ge \min(\Omega, N)$, this implies $\min(\Omega, N) < |\Omega - N|$, which, as when $s_p = s_q$, requires $\beta < 1/2$ or $\beta > 2$.

In summary, the resonant surface does not exist if $s = s_p = s_q$ or $1/2 \le \beta \le 2$. Otherwise it exists for $s_p = s_q$ provided (3.20) holds and for $s_p = -s_q$ when (3.21) applies.

## References

* Bartello, P. 1995 Geostrophic adjustment and inverse cascades in rotating stratified turbulence. J. Atmos. Sci., 52, 4410–4428.

* Benney, D.J., Saffman, P.G. 1966 Nonlinear interactions of random waves in a dispersive medium. Proc. Roy. Soc. A 289, 301–320.

* Benney, D.J., Newell, A.C. 1969 Random wave closures. Stud. Appl. Math. 48, 29–53.

* Bellet, F. 2003 Etude asymptotique de la turbulence d'ondes en rotation. PhD thesis, Ecole Centrale de Lyon.

* Bellet, F., Godeferd, F.S., Scott, J.F., Cambon, C. 2006 Wave turbulence in rapidly rotating flows. J. Fluid Mech. 562, 83–121.

* Cambon, C. 2001 Turbulence and vortex structures in rotating and stratified flows. Eur. J. Mech. B (Fluids) 20, 489–510.

* Cambon, C., Jacquin, L. 1989 Spectral approach to non-isotropic turbulence subjected to rotation. J. Fluid Mech. 202, 295–317.

* Cambon, C., Mansour, N.N., Godeferd, F.S. 1997 Energy transfer in rotating turbulence. J. Fluid Mech. 337, 303–332.

* Canuto, C., Hussaini, M.Y., Quarteroni, A., Zang, T.A. 1988 Spectral Methods in Fluid Dynamics. Springer-Verlag.

* Charney, J.G. 1948 On the scale of atmospheric motions. Geofys. Publ. Oslo, 17(2), 1–17.

* Charney, J.G. 1971 Geostrophic Turbulence. J. Atmos. Sci., 28, 1087–1095.
Coleman, G.N., Ferziger, J.H., Spalart, P.R. 1992 Direct simulation of the stably stratified turbulent Ekman layer. J. Fluid Mech. 244, 677-712. (Corrigendum: J. Fluid Mech. 252, 721.)

Deng, Y., Hani, Z. 2021 On the derivation of the wave kinetic equation for NLS. Forum of Mathematics, Pi 9:e6, 1-37.

Embid, P.F., Madja, A.J. 1998 Low Froude number limiting dynamics for stably stratified flow with small or finite Rossby numbers. Geophys. Astrophys. Fluid Dynamics, 87, 1-50.

Galtier, S. 2003 Weak inertial-wave turbulence theory. Phys. Rev. E 68, 015301.

Godeferd, F.S., Cambon, C. 1994 Detailed investigation of energy transfers in homogeneous stratified turbulence. Phys. Fluids 6, 2084-2100.

Godeferd, F.S., Staquet, C. 2003 Statistical modelling and direct numerical simulations of decaying stably stratified turbulence. Part 2. Large-scale and small-scale anisotropy. J. Fluid Mech. 486, 115-159.

Haugen, N.E.L, Brandenburg, A. 2004 Inertial range scaling in numerical turbulence with hyperviscosity. Phys. Rev. E 70, 026405.

Herbert, C., Pouquet, A., Marino, R. 2014 Restricted equilibrium and the energy cascade in rotating and stratified flows. J. Fluid Mech., 758, 374-406.

Hossain, M. 1994 Reduction in the dimensionality of turbulence due to a strong rotation. Phys. Fluids 6, 1077-1080.

Liechtenstein, L., Godeferd, F.S., Cambon, C. 2005 Nonlinear formation of structures in rotating stratified turbulence. J. Turbul. 6, 1-18.

Marino, R., Mininni, P.D., Rosenberg, D.L., Pouquet, A. 2013 Inverse cascades in rotating stratified turbulence: Fast growth of large scales. Europhys. Lett., 102, 44006.

Nazarenko, S. 2011 Wave Turbulence. Springer.

Newell, A.C., Rumpf, B. 2011 Wave turbulence. Annu. Rev. Fluid Mech. 43, 59-78.

Orszag, S.A. 1970 Analytical theories of turbulence. J. Fluid Mech., 363-386.

Orszag, S.A., Patterson, G.S. 1972 Numerical simulation of three-dimensional homogeneous isotropic turbulence. Phys. Rev. Lett. 28, 76-79.

Pedlosky, J. 1987. Geophysical Fluid Dynamics, Springer-Verlag, New York.

Sagaut, P., Cambon, C. 2018 Homogeneous Turbulence Dynamics (second edition). Springer.
Smith, L.M., Waleffe, F. 2002 Generation of slow large scales in forced rotating stratified turbulence. J. Fluid Mech. 451, 145-168.

Waleffe, F. 1991 Non-linear interactions in homogeneous turbulence with and without background rotation. Centre for Turbulence Research, Annual Research Briefs, 1991.

Zakharov, V.E., Lvov, V.S., Falkovich, G.E. 1992 Kolmogorov Spectra of Turbulence I – Wave Turbulence. Springer-Verlag.