The dynamics of the Earth's atmosphere is characterized by a wide spectrum of oscillations, ranging from hourly to interdecadal and beyond. The low-frequency component of the atmospheric variability cannot be understood solely in terms of linear atmospheric waves that have shorter timescales. A newly proposed mechanism, the precession resonance mechanism, is a regime of highly efficient energy transfer in the spectral space in turbulent systems. Here, we investigate the role of the precession resonance, and the alignment of dynamical phases, in the generation of low-frequency oscillations and the redistribution of energy/enstrophy in the spectral space using the barotropic vorticity equation. First, the mechanism and its ability to generate low-frequency oscillations are demonstrated in low-order models consisting of four and five nonlinearly interacting Rossby–Haurwitz waves. The precession resonance onset is also investigated in the full barotropic vorticity equation, and the results are in agreement with the reduced models. Efficiency peaks in the energy/enstrophy transfer also correspond to regimes of strong excitation of low-frequency oscillations. The results suggest that the organization of the dynamical phases plays a key role in the redistribution of energy in the spectral space, as well as the generation of low frequencies in the barotropic vorticity equation.

## I. INTRODUCTION

The dynamics of a barotropic nondivergent flow in a rotating frame has been extensively studied in the fluid dynamics literature due to its applicability as a prototype model for large- and planetary-scale motions of the atmosphere and ocean.^{1,2} According to the seminal scale analysis of large-scale motions in mid-latitudes performed by Charney,^{3} these large-scale atmospheric disturbances exhibiting a high spatial coherence, and a typical timescale ranging from a few days to a week, are characterized by a quasi-horizontal and quasi-solenoidal flow. Consequently, the atmospheric flow can be thought of as a “noisy sea” of fast propagating gravity and acoustic waves embedded in such large-scale disturbances. The linear wave theory of the so-called barotropic nondivergent model was introduced by Rossby^{4} with the *β*-plane approximation and extended by Haurwitz^{5} for the spherical geometry context. This model was the one utilized by Charney *et al.*^{6} in the first known successful numerical weather prediction and still to date constitutes an essential piece of intermediate steps in the development of general circulation models.^{7} Such wave motions have been labeled as Rossby–Haurwitz waves, and the linear wave theory has been further extended some decades later with the inclusion of wave forcings^{8} and/or more realistic background flows.^{9–13}

Apart from the linear wave theory, a significant research effort has also been made on the role of nonlinearity on the barotropic nondivergent model. Platzman^{14} described the spectral representation of the nonlinear version of the barotropic nondivergent model, with the explicit computation of the coupling coefficients among each spherical harmonic triplet, and discussed some issues regarding the truncation influence on the integral constraints arising from the original model equations. Platzman^{15} described the time evolution of the nonlinear spectral equations by analytically solving a highly truncated version of the system for a single triad of interacting harmonics in terms of Jacobian elliptic functions. Longuet-Higgins *et al.*^{16} demonstrated the possibility of triads of Rossby waves to match the so-called “resonance relations,” which impose the relative linear phase of the wave triplet not to change in time. The wave triplets satisfying these resonance relations are called resonant triads and are believed to be associated with the leading-order nonlinear effects for a weakly nonlinear regime of any dispersive wave system having quadratic nonlinearities. A background review of the three-wave resonant interaction can be found in Craik^{17} for many fluid dynamics problems, or in Pedlosky^{18} and Lynch^{19} for the Rossby wave context. Longuet-Higgins *et al.*^{16} also computed the energy transfer rate associated with the resonant Rossby wave triads as a result of the integral constraints imposed by energy and enstrophy conservations. The findings of Longuet-Higgins *et al.*^{16} and Platzman^{15} verified the constraints of the modal energy transfer due to the so-called Fjørtoft's theorem for a two-dimensional nondivergent flow.^{20} This theorem states that, as a consequence of total energy and total enstrophy conservations, the mode with the intermediate total wavenumber either receives energy from or supplies energy to the remaining modes of an interacting triad. For a resonant triad, this mode has also the highest absolute eigenfrequency and is usually referred to as “pump mode” from the plasma physics literature.^{21}

Consequently, resonant triads have been useful to understand the stability properties of a single Rossby mode to small-amplitude perturbations, for either a *β*-plane^{22} or the spherical geometry^{23,24} contexts. Due to the aforementioned energy constraints for a triad interaction resulting from energy and enstrophy invariants, a single wave mode can only be unstable to small-amplitude perturbations^{25} if (i) it is the pump mode of a triad interaction and (ii) the mismatch among the linear eigenfrequencies of the modes of this triad be smaller than a critical value that depends on the initial amplitude of the pump wave, as will be shown in Sec. II A. Thus, in the weakly nonlinear regime, resonant triads are crucial for triggering the wave instability. Lynch^{24} performed numerical simulations with the barotropic nondivergent model with a vorticity forcing that resonantly excites a specific Rossby–Haurwitz wave mode. He showed that, after being resonantly excited by the forcing, this mode preferentially loses its energy to modes that constitute a nearly resonant triad with it.

Nevertheless, although resonant triads have shown to be crucial in triggering the instability of the primary wave field, the spectral broadening of the solution due to energy cascade processes that characterize the transition to turbulence appears to rely on nonresonant triad interactions. Reznik *et al.*^{26} studied the nonlinear interaction of Rossby–Haurwitz waves and showed that relaxing the resonance condition for the mode eigenfrequencies to consider near resonances increases the spectral energy redistribution. In fact, in the weakly nonlinear limit, as the wave triads are assumed to satisfy the resonance condition for the linear eigenfrequencies, the possible effect of the wave phases on the time evolution of the wave amplitudes is neglected, leading to the so-called random phase hypothesis.^{27,28} On the other hand, by relaxing the weakly nonlinear limit to consider off-resonant triads, some studies have pointed out the important role of the nonlinear dynamical phases on the energy flow throughout the modal space.^{29–31} Chian *et al.*^{31} showed that an amplitude-phase synchronization associated with multiscale interactions in a one-dimensional model of regularized long-wave equation is responsible for an on–off intermittency at the onset of spatiotemporal chaos. This mechanism has also been demonstrated in other turbulence contexts, such as magnetized three-dimensional Keplerian shear flows.^{32}

Another mechanism that has shown to be responsible for a significant energy flow throughout the whole spectral space refers to the so-called precession resonance mechanism. The name of this mechanism was chosen due to the precession of the nonlinear dynamical phases, and it has a loose analogy to the precession of celestial bodies, but should not be confused with it. The mechanism has been proposed by Bustamante *et al.*^{33} to maximize the energy transfer between different interacting wave triplets and to enhance the cascade process, in the context of drift waves. This resonance takes place when there is a balance between the characteristic precession frequency of the combined phase of one triad and the nonlinear frequency of amplitude modulation of another wave triplet. According to Bustamante *et al.*,^{33} for a moderate wave amplitude regime, the characteristic precession frequency of the combined phase of an interacting triad receives significant contribution from the mismatch among its linear eigenfrequencies. This mechanism has also been analyzed in the context of Burgers' equation in Murray and Bustamante,^{34} who showed that bursts of strong energy transfer throughout the spectral space are interchanged with periods of weak interactions; this interchange between strong and weak energy fluxes throughout the spectral space was found to be associated with synchronized and unsynchronized regimes, respectively, between amplitudes and dynamical phases of the wave triplets. See also Arguedas-Leiva *et al.*,^{35} where the relation between phase synchronization across spatial scales and turbulent intermittency was demonstrated. Furthermore, the deep connection between precession resonance and the region of convergence of the wave-turbulence normal form transformation in the barotropic vorticity equation was studied in Walsh and Bustamante.^{36}

Another important aspect of the precession resonance refers to its ability in generating low-frequency fluctuations in a wave system, as demonstrated by Raphaldini *et al.*^{37} in the context of magnetohydrodynamic (MHD) Rossby waves in the solar tachocline. According to Raphaldini *et al.*,^{37} the inter-triad energy transfers associated with the precession resonance occur in a longer timescale than the intra-triad energy exchanges, thus yielding longer-period modulations on the energy exchanges within a resonant triad. This aspect of the precession resonance should also be useful in the atmospheric dynamics context. In fact, in the atmospheric flow, the low-frequency variability is defined as oscillations with periods beyond 10 days and therefore is responsible for modulating the weather systems and consequently yielding climate anomalies. To date, the mechanisms behind the generation of the atmospheric low-frequency variability still constitute a topic of current research in the atmospheric sciences, as these timescales are longer than the characteristic periods of the atmospheric normal modes.^{1} Nonlinear triad interactions may play some roles in the dynamics of these low-frequency modes. This has been suggested by Kartashova and L'vov^{38} and Raupp and Silva Dias^{39} for the case of intraseasonal variations, an important component of the atmospheric low-frequency variability.

Here, we analyze the precession resonance mechanism in the barotropic nondivergent model, with parameters relevant for the study of atmospheric flows. The main purpose is to study the role of this mechanism in the redistribution of energy in spectral space and the generation of low-frequency oscillations. We first investigate the mechanism in the reduced spectral equations of four (two interacting triads connected by two modes) and five (two interacting triads connected by a single mode) wave modes, selecting triads involving a prevalent atmospheric mode that has been shown to be unstable in previous studies.

This mode is the pump wave in many triad interactions, making the excitation of neighboring modes in spectral space efficient. We vary the initial conditions and show that the maximal efficiency of energy/enstrophy transfer between two different triads is associated with the precession resonance regime. Furthermore, numerical simulations with the full spectral equations (i.e., considering a regular truncation), with the initial condition projected only onto a resonant triad, demonstrate that the precession resonance associated with the four-wave problem is responsible for the energy leakage from the primary resonant triad and, therefore, the initiation of the cascade process. The results of both the reduced and the full spectral models show that the precession resonance yields low-frequency modulations in the energy exchanges associated with the primary resonant interaction. Therefore, based on our theoretical and numerical results, we conjecture that this mechanism may be responsible for generating considerable low-frequency variability in the atmosphere and ocean.

## II. BASIC EQUATIONS

The dynamics of a barotropic nondivergent flow in a rotating frame is governed by the so-called barotropic nondivergent model. This model is described by the conservation of absolute vorticity, which can be written in the conservative form as follows:

In the equation above, $v\u2192(\varphi ,\theta )$ is the solenoidal velocity field, and $\zeta =(\u2207\xd7v\u2192)\xb7r\u2192$ is the relative vorticity, with $r\u2192$ representing a radial unit vector with respect to the sphere; $f=2\Omega E\u2009sin\u2009(\theta )$ is the planetary vorticity, where Ω_{E} is the Earth's rotation rate and $(\varphi ,\theta )$ refers to the longitude–latitude coordinate system ($\varphi \u2208[0,2\pi ]$ and $\theta \u2208[\u2212\pi /2,\pi /2]$). Using the streamfunction relations

and

Eq. (1) now reads

where $J(\xb7,\xb7)$ is the Jacobian operator in spherical coordinates

for any differentiable functions *h* and *g*, and *a* denotes the Earth's radius.

From (1), it follows that the mean relative vorticity is conserved. However, the conserved quantities that are relevant for the wave interactions are those having a quadratic dependence in terms of the disturbance variables.^{16,20,40} The quadratic conserved quantities of (1) refer to total energy *E* and total enstrophy $E$, defined, respectively, by

where *S* represents the surface of the sphere.

Equation (2) can be solved according to the following expansion in terms of the eigensolutions of its linearized version:

where $i=\u22121$ refers to the imaginary unit. In (6), $An,m(t)$ denotes the complex-valued spectral amplitudes, *m* is the zonal wavenumber, and *n* is the meridional structure index; $Ynm(\varphi ,\theta )$ are the spherical harmonics, which are eigenfunctions of the spherical Laplace's operator

$Pnm(sin\u2009\theta )$ indicates the associated Legendre polynomials, and the eigenfrequencies $\omega n,m$ satisfy the dispersion relation of Rossby–Haurwitz waves

Expansion (6) is an exact solution of the barotropic vorticity equation (2) provided the spectral amplitudes satisfy

In the equation above, we have adopted a simplified notation in which each spherical harmonic is specified by a single index, viz., $(nj,mj)\u2192j,\u2009(nk,mk)\u2192k$ and $(nl,ml)\u2192l$. Thus, $\delta k,lj=\omega k+\omega l\u2212\omega j$ is the mismatch among the mode eigenfrequencies of a wave triplet {*j*, *k*, *l*}, and $Ck,lj\u2208\mathbb{R}$ refers to their nonlinear coupling coefficients, given by

where $z=sin\u2009\theta $. Notice that $Ck,lj=Cl,kj$ in general. In order for the nonlinear coupling coefficients to be nonzero, the respective mode indices must satisfy the so-called Ellsaesser conditions^{41}

The coupling constants $Ck,lj$ contain all the information on the nonlinearity of the original equation (2). In addition, as a consequence of the integral constraints (4) and (5), for any triad of wave modes satisfying conditions (10), the nonlinear coupling coefficients must satisfy the relations (see the Appendix A for details)

Equation (11) refers to the Fjørtoft's theorem for a two-dimensional nondivergent flow.^{20} It implies that

Thus, the mode with intermediate total wavenumber will always receive energy from or supply energy to the remaining triad components. This defines the pump mode. For resonant triads ($\delta k,lj=0$), the pump mode has also the highest absolute eigenfrequency. Condition (11) also implies that total energy and total enstrophy are conserved for an arbitrary truncation of the spectral equations (8). Condition (11) has also been demonstrated by Longuet-Higgins *et al.*^{16} for both the barotropic vorticity equation and the quasi-geostrophic potential vorticity equation in which a finite divergence is considered in the *β*-plane context. Ripa^{40} presented a more general framework valid for nonlinear systems having two quadratic conserved quantities and obtained a similar result, but with the total wavenumber squared in (11) being replaced by the slowness index defined as the inverse of the wave phase speed. However, as can be noted from (7), for Rossby waves the slowness index is proportional to the total wavenumber squared so that (11) is consistent with Ripa's theoretical results.

Before describing the precession resonance mechanism in more detail, it is suitable to review some basic concepts of a single wave triad dynamics.

### A. Three-wave interaction and wave stability

Let us truncate expansion (6) to consider only a set of three spherical harmonics $(n1,m1),(n2,m2),(n3,m3)$ satisfying conditions (10). In this context, Eqs. (8) now read

^{24}that the coefficients $C2,31,C3,12$, and $C1,23$ are real and have the same sign for any given triad satisfying the Elsaesser conditions (10). For simplicity of presentation, we will assume in this subsection that the mode

*A*

_{3}is a pump mode, namely, that $C1,23$ is the coefficient with the largest size, as per Eq. (12). The equations above constitute the most elementary truncation of (8) that exhibits nonlinearity. Following Teruya and Raupp,

^{42}to analyze the stability of a single wave mode in a triad interaction, it is suitable to eliminate the oscillatory exponential factor in the RHS of (13) by making the transformation $A2=A\u03032ei\delta t$. In this context, the above system can be written as

where we have omitted the “∼” in the equations above for simplicity. Also, assuming that the pump mode holds most of the initial energy, $|A3(0)|\u226b|A1(2)(0)|$, the system above may be approximated by its linearized version

Thus, instability occurs if

and

The first condition is always satisfied as per Lynch's result previously mentioned. For obvious reasons, this instability is called pump wave instability.^{21} The second condition implies that off-resonant triads require a higher amplitude level for pump wave instability to be possible. Consequently, for an initial amplitude level compatible with a weakly nonlinear regime, resonant and nearly resonant triads play a crucial role to trigger wave instability.

The above theoretical result on the pump wave instability in an interacting triad is consistent with numerical studies on Rossby wave instability. The stability of a Rossby wave field has been a topic of considerable research effort due to its role in initiating the energy cascade process responsible for the transition to a turbulent regime and therefore the breakdown of weather predictability. Lorenz^{43} studied the barotropic instability of a flow associated with a Rossby wave mode, and Hoskins^{44} performed numerical simulations to analyze the stability of a Rossby–Haurwitz wave field superimposed to a zonal flow. Gill^{22} and Baines^{23} were perhaps the firsts to study the Rossby wave instability on the optics of nonlinear triad interaction. Gill^{22} showed that a Rossby wave mode is always unstable irrespective of its amplitude, but for small amplitudes, the destabilizing disturbances always constitute a resonant triad with the primary wave. However, Gill^{22} considered an infinite $\beta \u2212$ plane in which the continuum spectrum allows a much broader possibility of resonant triads. By contrast, in the spherical geometry, Baines^{23} found a narrower range of possible unstable wave modes than Gill,^{22} with some modes requiring a high amplitude to be unstable, possibly due to the more restricted set of resonant triads in the spherical geometry as a consequence of the discrete spectrum of wave modes.

Moreover, comparing the studies of Hoskins^{44} and Baines^{23} in the spherical geometry context, Hoskins^{44} found a more restricted set of unstable modes, showing that only modes with *n *>5 are unstable; however, his spectral model considered only zonal wavenumbers that are multiples of 4, therefore excluding many possibilities of interacting triads. On the other hand, Baines^{23} considered a larger number of destabilizing disturbances (30 modes) and obtained that traveling wave modes having wavenumbers *n *=* *3 and 4 can be unstable, with required amplitudes being smaller than those proposed by Hoskins.^{44} This difference between the findings of Hoskins^{44} and Baines^{23} can be attributed to the higher possibilities of resonant triad interactions in the latter's spectral equations. In fact, Lynch^{24} demonstrated that the instability of the wave mode $(n,m)=(5,4)$, which had also been verified in the numerical simulations of Thuburn and Li,^{45} is due to the nearly resonant triad interaction with modes (3, 1) and (7, 3). These wave modes will be considered in our numerical simulations displayed in Secs. III and IV.

The nonlinear system (13) has the following conserved quantities: the Manley–Rowe invariants

and the Hamiltonian

With these quantities, the triad equations may be integrated by quadrature, reducing it to a single equation. The solution is periodic and may be represented in terms of Jacobian elliptic functions, as will be shown in Sec. II B. The single triad system with a constant forcing, which is relevant in many contexts,^{39} was also shown to be integrable.^{46}

### B. The precession resonance

The precession resonance was introduced by Bustamante *et al.*^{33} as a highly efficient energy transfer mechanism throughout the modal space that takes place when multiple wave triads are connected. In order to introduce this concept in more detail, it is convenient to use the polar representation of the spectral coefficients

where $Bj=|Aj|$ and $\phi j=arg(Aj)$. In this representation, the system of three coupled complex valued differential equations (13) is reduced to a system of four real-valued equations, three for the moduli of the spectral amplitudes and one for the combined phases

where

is the so-called dynamical phase, and $\delta :=\omega 1+\omega 2\u2212\omega 3$ is the mismatch among the linear eigenfrequencies. In the case *δ* = 0, the solution for the moduli *B _{j}* can be expressed as (e.g., Refs. 15 and 47) follows:

where

with $sn,\u2009cn$, and $dn$ representing the Jacobi elliptic functions (sine, cosine, and delta amplitude, respectively) and $K(\mu )$ the elliptic integral of the first kind

where the parameter *μ* is given by

with $\rho =I23/I13$ and $\nu \u2208[0,\pi ]$ being an angle such that

In addition, in Eqs. (20), the constants *C*_{1}, *C*_{2}, and *C*_{3} refer to the nonlinear coupling coefficients, $C1=C2,31,\u2009C2=C1,32$, and $C3=C1,23$. The solutions for the amplitudes and dynamical phase presented above are periodic in time, with period *T* being defined as follows:

In the case $\delta \u22600$, the solutions have a similar form. In particular, the system is still integrable, and the amplitudes $B1,B2,B3$ oscillate with a period *T*, satisfying a formula similar to Eq. (25). See Harris,^{48} Chap. 3.2, for details, including the solution for the dynamical phase. Crucially, a remarkable difference is noted in the behavior of the dynamical phase $\Phi $. In the case *δ* = 0, the dynamical phase $\Phi $ is periodic in time, with period *T*, but when $\delta \u22600$ is large enough, the dynamical phase $\Phi $ has phase slips that occur periodically with period *T*, leading to an overall “drift” or “precession” when considering the dynamical phase as being defined on the torus $T$. Hence, in general one can define two characteristic frequencies associated with the dynamics of an interacting wave triad: $\Gamma =2\pi /T$, the frequency of energy exchanges among the waves, and Ω, the so-called precession frequency of triad (1, 2, 3), given by

with the dynamical phase $\Phi $ being defined in (19) and $\u27e8\xb7\u27e9$ indicating the time average. A schematic view of the meaning of the precession frequency will be shown in Sec. III A, Fig. 3.

As an integrable system, the triad equations are completely characterized by these two frequencies (Ω and Γ). In general, when these frequencies are incommensurable ($\Omega /\Gamma \u2208\mathbb{R}$\$\mathbb{Q}$), the system will be quasi-periodic and the orbits will live in a two-dimensional torus $T2$. By contrast, for special cases in which $\Omega /\Gamma $ is a rational number, the solution is periodic in time. Thus, by manipulating the initial conditions, one can change the values of the constants of motion and consequently these characteristic frequencies of the solution. An easy way to do this is to re-scale the initial amplitudes as $Bj(0)\u2192\alpha Bj(0)$, *j *=* *1, 2, 3, with *α* being a nondimensional real number. From Eq. (25), it is evident that the period of energy exchange, *T*, scales as the inverse of the initial amplitudes ($T\u221d1/I\u221dBj\u22121$), where *I* is a Manley–Rowe constant. This type of manipulation will be performed throughout this study in cases of nonintegrable systems of nonlinearly interacting wave modes as a way to explore different balances between the different frequencies of the system.

When two triads, *a* and *b*, are connected by either four- or five-wave configurations, although the system is no longer integrable, it is still useful to characterize each of the triads by these two frequencies, $(\Gamma a,\Omega a)$ and $(\Gamma b,\Omega b)$. Assuming that the energy of triad *a* is initially much higher than that of triad *b*, the coupling terms between the two triads will make the primary triad, triad *a*, to act as a forcing for the secondary one, triad *b*. In this context, resonance occurs when the precession frequency of triad *b* matches one of the Fourier harmonic frequencies of the amplitude modulation of triad *a*, viz.,

where $p\u2208\mathbb{Z}$. The equation above defines the so-called precession resonance between triads *a* and *b*.^{33} In practical applications, the primary triad that contains most of the initial energy is nearly resonant ($\delta a\u22480$), as discussed in Sec. II A. In addition, in the secondary wave triplet, the amplitude of one (two) of the modes is negligible for the four (five)-wave system. As we will see in Sec. II C, in the five-wave system, a small perturbation in the energy of the secondary modes of triad *b* is required for these modes to be excited; for the four-wave system (Sec. II D), however, as the triads are connected by two modes, the secondary mode of triad *b* can be excited even if it has zero energy initially. Consequently, the precession frequency of the secondary triad is initially dominated by the mismatch among the linear eigenfrequencies, $\Omega b\u2248\delta b$. However, as time evolves and the secondary modes of triad *b* are sufficiently excited, the nonlinear terms involving the mode amplitudes in (18d) will become important as well so that the synchronization defined by (27) will involve the fully nonlinear dynamical phases of the triads, as will be shown in Secs. III and IV.

Bustamante *et al.*^{33} showed that, in the precession resonance regime, there is a peak in the efficiency of energy exchange between the triads. This process gets more intricate in systems involving many triad interactions, such as in full partial differential equation models. The effects of the precession resonance are still noticeable in these full spectral models.^{33,34} Here, we first investigate the onset of the precession resonance mechanism in the reduced models of five- and four-wave modes shown, respectively, in Secs. II C and II D. We demonstrate that, in the precession resonance regime, the solutions display long-term modulations as a result of inter-triad energy exchanges, an aspect that was not explored in Bustamante *et al.*,^{33} but may have important consequences in the field of atmospheric fluid dynamics. We then proceed to investigate the effects of the precession resonance in the full barotropic vorticity equation, focusing on the redistribution of energy and enstrophy in the large-scale circulation patterns.

Among the possible types of precession resonance, one of them, the phase alignment or synchronization, which consists of the case *p *=* *0 in Eq. (27), is verified in our numerical results presented in Secs. III and IV and seems to be of great relevance in several other contexts, for instance, in the dynamics of the Burgers' equation^{34,49,50} and in MHD waves in the solar tachocline.^{51} This mechanism fits in a more general context as one of the mechanisms involving the role of the dynamical phases' organization on the nonlinear energy transfer; other studies within this research line include Chian *et al.*^{31} in the context of the regularized long-wave equation and Miranda *et al.*^{32} in the context of 3D astrophysical flows. The organization of the waves' phases in those settings is often related to coherent structures in the physical space.^{35}

In the simplest setting, the precession resonance can be introduced in systems with four- or five-wave modes arranged in two triads. These reduced spectral models of (2) will be introduced in Secs. II C and II D and will be studied by numerical simulations in Sec. III with realistic atmospheric parameters.

### C. Five-wave model

As a first step in complexity to study the precession resonance, one might consider a set of two-wave triads interacting with each other through one common mode. If one truncates expansion (6) to consider only five spherical harmonics $(n1,m1),(n2,m2),(n3,m3),(n4,m4)$, and (*n*_{5}, *m*_{5}) in which modes (1, 2, 5) and (3, 4, 5) satisfy (10), it follows that the corresponding spectral amplitudes satisfy

A set of independent Manley–Rowe invariants of system (28) is given by

### D. Four-wave system

A second way to connect two interacting wave triads is by having two modes in common. In this way, if one truncates expansion (6) to consider only four spherical harmonics $(n1,m1),(n2,m2),(n3,m3)$, and (*n*_{4}, *m*_{4}) in which modes (1, 2, 3) and (1, 2, 4) satisfy conditions (10), it follows that the corresponding spectral amplitudes satisfy

where $\delta 1=\omega 1+\omega 2\u2212\omega 3$ and $\delta 2=\omega 1+\omega 4\u2212\omega 2$ are the mismatch frequencies of each triad. The invariants of this system are

As in the five-wave model, the invariants above are not enough to make system (31) integrable in general.

## III. PRECESSION RESONANCE AND LOW-FREQUENCY OSCILLATIONS IN REDUCED SPECTRAL MODELS

Here, we analyze inter-triad energy and enstrophy fluxes in the highly truncated systems (28) and (31). As discussed in Secs. I and II, resonant and nearly resonant triads play a crucial role in the initial pump wave instability, whereas the energy leakage from the primary (nearly) resonant triad associated with the spectral broadening of the solution can also be due to off-resonant triads. Thus, we have integrated these two systems for some representative examples consisting of a base triad composed of Rossby–Haurwitz waves with spherical harmonics (5, 4), (3, 1), and (7, 3). As mentioned in Sec. II A, this triad is a nearly resonant one and, as demonstrated by the numerical simulations with the barotropic vorticity equation performed by Lynch,^{24} constitutes the leading destabilizing disturbances for the unstable mode (5, 4). Other works have investigated the role of this mode,^{23,24,45} finding it to be unstable by considering a larger number of modes in the truncation of the barotropic nondivergent model (and, therefore, a larger number of destabilizing disturbances). Lynch^{24} considered in his numerical simulations a vorticity forcing resonating with the Rossby–Haurwitz mode (5, 4). After this mode is sufficiently excited, it was found that it preferentially excites modes (3, 1) and (7, 3), with a complete cycle of energy exchanges within this triad occurring until energy leaks to secondary modes. In addition, observational studies of the normal mode spectra based on reanalysis data sets have shown evidence to support the important role of mode (5, 4). Gehne and Kleeman^{53} showed that the main peak in the wavenumber-frequency spectrum of the zonal velocity field corresponds to barotropic Rossby–Haurwitz waves having zonal wavenumbers m = 3–4 and a meridional structure compatible with the spherical harmonic (5, 4), although they utilized a different basis function set given by the Hough vector harmonics.^{54–57}

In order to demonstrate the effect of the precession resonance in maximizing the inter-triad energy transfer, in the integration of the five-wave system (28), we consider a triad cluster in which the base triad is connected with modes (7, 4) and (3, 0) via the pump mode (5, 4). For the four-wave system integration, the base triad is coupled with mode (9, 2) through modes (3, 1) and (7, 3).

A schematic view of these triad systems is presented in Fig. 1. The parameter range of the experiments is chosen so that the resulting velocity fields exhibit physically realistic values (i.e., $|v\u2192|<50\u2009m/s$). While a systematic study of the efficiency landscape in a generic cluster of four or five waves is beyond the scope of this study, the reader is referenced to Murray^{50} (Chap. 6) for a more detailed study along this line.

### A. Example 1: Four-wave model

Here, we show results of numerical integrations of system (31) to illustrate the precession resonance mechanism and its ability to maximize the efficiency of energy transfer between two different wave triads and to yield low-frequency oscillations. To perform these numerical simulations, we select a representative example of four spherical harmonics obeying the four-wave system (31) given by ${(n1,m1),(n2,m2),(n3,m3),(n4,m4)}$ = ${(3,1),(7,3),(5,4),(9,2)}$ (Fig. 1). In this case, the corresponding mode eigenfrequencies (normalized by $2\Omega E$) are $\omega 1=\u22120.083,\u2009\omega 2=\u22120.053,\u2009\omega 3=\u22120.133$, and $\omega 4=\u22120.022$, corresponding to the mismatch frequencies $\delta 1=\omega 1+\omega 2\u2212\omega 3=\u22120.003$ for the first triad and $\delta 2=\omega 1+\omega 4\u2212\omega 2=\u22120.051$ for the second one. The corresponding nonlinear coupling coefficients are given by $(C2,31,C1,32,C1,23)=(\u22120.773,\u22122.497,\u22123.270)$ (triad 1) and $(C2,41,C1,24,C1,42)=(\u22120.569,\u22125.527,\u22126.096)$ (triad 2).

According to Eq. (26), we further define the precession frequencies associated with the two interacting triads of the system

where

and $\phi \u0303j=\phi j\u2212\omega jt$ are the original phases of the waves [namely, the argument of $Aj(t)e\u2212i\omega jt$]. The integrations have been performed by using an explicit eighth-order Runge–Kutta numerical scheme. The initial condition considered in our numerical simulations is given by $(A1,A2,A3,A4)=\alpha \xd7(1.0,1.0,1.0,1.0\xd710\u22123)/a$, where *α* is a free parameter to control the magnitude of the initial wave fields and, consequently, the frequency of the energy modulation of the base triad interaction. Thus, mode 4 has very small initial energy compared with the modes of triad $({1,2,3})$ (six orders of magnitude smaller), since it will be considered as the target mode to analyze the efficiency of the inter-triad energy transfer. The integrations have been performed by varying the parameter *α* from 0 to 100 in order to find the maximum efficiency of energy transfer from the primary [(1, 2, 3)] to the secondary [(1, 2, 4)] triad that characterizes the precession resonance.^{33} The initial mode amplitudes described above constitute a simple way of investigating the precession resonance mechanism in which the characteristic frequencies of the primary triad are changed without modifying the initial energy distribution of this triad. However, a similar analysis on the inter-triad energy exchange efficiency can also be made considering different regimes of energy exchange associated with the primary triad by independently varying two of its initial mode amplitudes. This analysis is presented in Appendix B.

In order to measure the efficiency of the energy transfer from the primary to the secondary triad, we define the following quantity that is bounded between 0 and 1 and is zero when all the energy of the system is in the primary triad composed of modes 1, 2, and 3:

where $\kappa j=nj(nj+1)$, with the efficiency being defined as $\epsilon =maxt(Q(t))$. Similarly, in the enstrophy case

with the efficiency being again defined as $\epsilon \u0303=maxt(Q\u0303(t))$.

Figure 2 displays the inter-triad energy [Fig. 2(a)] and enstrophy [Fig. 2(b)] transfer efficiencies, the corresponding power spectrum in the low-frequency range [Fig. 2(c)] and the precession frequencies of the two-wave triplets [Fig. 2(d)]. All these quantities are displayed as a function of the initial condition free parameter *α*. The spectral power in the low-frequency range is computed from the time evolution of a particular mode energy (mode *j*) by

where $\omega \u0303$ is defined such that $2\pi /\omega \u0303=10$ days, and the “ $\u0302$ ” indicates the time Fourier transform.

From Figs. 2(a) and 2(b), one notices that the energy/enstrophy transfer efficiencies are characterized by a broad peak centered at about $\alpha \u224815$. From this value of *α*, the efficiency decreases until $\alpha \u224826$. For larger values of *α*, the energy/enstrophy transfer efficiencies increase with the initial condition-free parameter. It is also noticeable that the corresponding power spectrum in the low-frequency range [Fig. 2(c)] exhibits a high power in low frequencies only in the first half of the peak of the energy/enstrophy transfer efficiency, followed by an abrupt decrease in the low-frequency power for increasing values of *α*. Comparing Fig. 2(d) with Figs. 2(a)–2(c) shows that the energy/enstrophy transfer efficiencies and low-frequency spectral power are determined by a combination of both precession frequencies. While the *α* range corresponding to the region with high power in the low frequencies corresponds to a case of precession resonance characterized by the phase alignment of $\Omega 1,42$, the beginning and end of the broad efficiency peak are marked by sharp changes in the regime of $\Omega 1,23$.

In order to illustrate the meaning of the precession frequencies, Fig. 3 shows the time evolution of the original individual wave phases $\phi \u0303j=\phi j\u2212\omega jt$ (left panels) and the dynamical phases of triads (1, 2, 3) and (1, 2, 4) (middle panels), together with the corresponding estimated precession frequency of triad (1, 2, 4) (right panel), referred to some of the integrations illustrated in Fig. 2. We recall that the phases of the individual wave modes are defined by $\phi j=arg(Aj),j=1,2,3,4$. Thus, the phases are unwrapped, meaning that, for each complete cycle, we sum a $2\pi $ increment to the evaluated angle. Figures 3(a), 3(b), and 3(e) refer to an integration near the precession efficiency peak ($\alpha =7.1)$, while Figs. 3(c) and 3(d) refer to an integration away from the inter-triad energy transfer efficiency peak ($\alpha =25.0$). The triad phases are given by $\Phi 123(t)=\phi 3(t)\u2212\phi 1(t)\u2212\phi 2(t)+\delta 1\u2009t$ and $\Phi 142(t)=\phi 2(t)\u2212\phi 1(t)\u2212\phi 4(t)+\delta 2\u2009t$. Graphically, the precession frequency in the time interval $[T1,T2]$ can be estimated from the mean slope of the triad phase curves

One observes that, for $\alpha =7.1$ [Figs. 3(a) and 3(b)], while the precession frequency $|\Omega 123|$ is nearly 0, the precession frequency of the second triad, $|\Omega 124|$, is much higher in comparison. With an increase in the energy level of the initial condition by setting $\alpha =25.0$ [Figs. 3(c) and 3(d)], the evolution of the individual phases changes completely [see Fig. 3(a) in comparison with Fig. 3(c)], most notably for modes 2 and 3 (mode 3 in fact reverses its propagation from westward-anticlockwise to eastward-clockwise). The resulting triad phases change their behavior; while $\Phi 123$ has a much higher slope, $\Phi 124$ decreases its slope, now growing slower than the former. This means that $|\Omega 123|>|\Omega 124|$ in this case.

To confirm whether the precession resonance is responsible for generating low-frequency modulations, Fig. 4 shows the spectral power of the solution in low frequencies (periods longer than 10 days), together with the corresponding time evolution of the mode energies, for different values of the parameter *α*. It is evident from Fig. 4 the strong correlation between the low-frequency power spectrum of the solution and the efficiency of energy transfer between the two-wave triplets. For *α* = 1, one notices the target mode being excited with very weak efficiency of about 0.3%, with a primary periodicity of about 10 days and a secondary low-frequency periodicity of about 90 days. Thus, although a low-frequency fluctuation in the four-wave system solution is generated in this case, its energy/power is very weak. For $\alpha =7.1$, which is close to the precession resonance regime, one notices a much higher inter-triad energy transfer efficiency of about 9% and a stronger spectral power in low frequencies, predominantly with an approximate 15-day periodicity. For $\alpha =25.0$, one notices a substantial decrease in the low-frequency power, with characteristic oscillations being dominated by periodicities shorter than 2 days.

### B. Five-wave model

Here, we consider a pair of wave triads consisting of modes ${(n1,m1),(n2,m2),(n3,m3),(n4,m4),(n5,m5)}={(3,0),(7,4),(3,1),(7,3),(5,4)}$ obeying system (28) (Fig. 1). In this case, the corresponding mode eigenfrequencies (normalized by $2\Omega E$) are given by $\omega 1=0.0,\u2009\omega 2=\u22120.071,\u2009\omega 3=\u22120.083,\u2009\omega 4=\u22120.053$, and $\omega 5=\u22120.133$, corresponding to the mismatch frequencies $\delta 3=\omega 1+\omega 2\u2212\omega 5=0.062$ for the first triad and $\delta 4=\omega 3+\omega 4\u2212\omega 5=\u22120.003$ for the second one. The corresponding nonlinear coupling coefficients are $(C2,51,C1,52,C1,25)=(0.655,2.116,2.771)$ [triad (1, 2, 5)] and $(C4,53,C3,54,C3,45)=(\u22120.773,\u22122.497,\u22123.270)$ [triad (3, 4, 5)]. Similarly to the previous example, the initial condition adopted in the numerical integrations is given by $(A1,A2,A3,A4,A5)=\alpha \xd7(1.48\xd7exp\u2009(1.530i)\xd710\u22128,1.71\xd7exp\u2009(0.304i)\xd710\u22128,1.03\xd7exp\u2009(1.110i)\xd710\u22125,1.22\xd7exp\u2009(5.060i)\xd710\u22125,1.17\xd7exp\u2009(3.950i)\xd710\u22125)$. Thus, modes 1 and 2 are the target modes in this case, having their energy six orders of magnitude smaller than the modes of the primary triad, (3, 4, 5). Consequently, the inter-triad energy transfer efficiency is defined by using the quantity

where $\kappa i=ni(ni+1)$. The efficiency is then defined as the maximum in time of this quantity, $\epsilon =maxt(Qt)$. For the enstrophy case, the analogue quantity is defined as

with again the efficiency being defined as the maximum in time of this quantity, $\epsilon \u0303=maxt(Q\u0303t)$. Also, similarly to the four-wave example, we define the precession frequencies by

where

Figure 5 shows the efficiency of energy [Fig. 5(a)] and enstrophy [Fig. 5(b)] transfers to the target modes, the low-frequency power spectrum computed from (37) [Fig. 5(c)], and the precession frequencies of the triads [Fig. 5(d)], with all these quantities being displayed as a function of the initial condition control parameter *α*. One can note that the maximum efficiency of energy/enstrophy transfer occurs for *α* = 31. This value corresponds to the initial condition in which the nonlinear frequency of energy modulation of triad (3, 4, 5) approximately equals the mismatch among the linear eigenfrequencies of modes 1, 2, and 5 (precession resonance regime). In addition, it is noticeable the significant positive correlation between the inter-triad energy transfer efficiency and the low-frequency power spectrum, with the low-frequency spectral power maximum coinciding with the maximum of the energy/enstrophy transfer efficiencies. From Fig. 5(d), one notices that, while the precession frequency $\Omega 3,45$ of the primary triad remains stable for all the range of *α*'s studied here, the precession frequency of the secondary triad, $\Omega 1,25$, varies as a function of the parameter *α*, with some abrupt changes that coincide with the interval with highly efficient energy and enstrophy transfers. In particular, the peak of the efficiency of energy and enstrophy transfers occurs approximately when $\Omega 1,25$ crosses 0. $\Omega 1,25\u22480$ means that there is an alignment of the phases of the respective triads. The alignment of the phases has been reported as a regime of highly efficient energy transfers in other contexts.^{34,35,58}

To further illustrate this effect, Fig. 6 shows the time evolution of the mode energies referred to numerical solutions for two values of *α* outside the precession resonance regime (*α* = 10 and $\alpha =60$) as well as the numerical solution within the precession resonance regime ($\alpha \u224827$) for comparison. One clearly notices that the solution for $\alpha \u224827$ presents a long-term modulation of about 50 days, while the others exhibit only faster oscillations. This is also confirmed by the respective power spectrum plots in the right panel of Fig. 6, where one notices a low-frequency spectral peak in the solution for the precession resonance regime $\alpha \u224827$, which is absent in the other two solutions. One also notices that the inter-triad energy transfer in this kind of triad cluster configuration is extremely inefficient outside the precession resonance regime, as can be evidenced by the low amplitude of the energy modulations associated with the solutions for *α* = 10 and *α* = 60. We remark that low frequencies are indeed allowed outside the precession resonance regime (see, for instance, the spectral peak on Fig. 6 for *α* = 10). However, the integrated power spectrum over low frequencies [Eq. (37)] tends to be weaker outside this regime. A similar analysis performed by independently varying two initial mode amplitudes of the primary triad is presented in Appendix B.

In summary, we have demonstrated in two reduced models the signatures of the precession resonance mechanism. None of these models are integrable, meaning that the nonlinear frequency Γ cannot be calculated analytically. Rather, this frequency can only be estimated numerically through the time spectrum of the solutions for the wave amplitudes. In this sense, a proxy of this spectral content is represented by the power at low frequencies given by Eq. (37). In both reduced models given by (28) and (31), the fact that the spectral power in the low-frequency range is clearly connected with the precession frequency is by itself a strong evidence of the precession resonance mechanisms. This matching between the spectral power at low frequencies and the precession frequency of the secondary triplet is also obtained by further varying the initial mode amplitudes in terms of two independent parameters, as demonstrated in Appendix B.

## IV. PRECESSION RESONANCE IN THE FULL PARTIAL DIFFERENTIAL EQUATION (PDE)

So far, we have presented the precession resonance mechanism in the context of highly truncated spectral models, consisting of four and five waves arranged in two triads. Although these highly truncated models clearly demonstrate how the precession resonance mechanism operates, the full partial differential equation (2) contains a variety of other competing interacting triads that may affect or even overshadow the solution obtained from the reduced models. Therefore, to test the robustness of the results presented in Sec. III, it is important to investigate how the precession resonance mechanism operates in the full partial differential equation (2).

### A. Numerical implementation details

Here, we use a pseudo-spectral, or collocation, method^{59} to solve the barotropic vorticity equation (1). The basic structure of the method involves the following: (i) first, to calculate the velocity $v\u2192$ in physical space from the spectral representation of the vorticity $(\zeta )$; (ii) then to calculate the product $(\zeta +f)v\u2192$ in physical space (therefore the terminology pseudo-spectral) and to convert it to spectral space; and (iii) finally, to apply the divergence ($\u2207\xb7$) in spectral space to obtain the tendency forcing $\u2212\u2207\xb7((\zeta +f)v\u2192)$. For the temporal discretization scheme, we use the classic three-stage third-order Runge–Kutta scheme.^{60} Spherical Harmonics transforms used the package provided by Schaeffer,^{61} with implementation within the shallow water modeling package of Schreiber *et al.*^{62}

In the numerical implementation, a spectral resolution of 128° (with triangular truncation), with correspondent physical resolution of 384 longitude vs 192 latitude points, ensuring no-aliasing for quadratic terms. A time step size of 2 min was used, which ensures numerical stability and almost negligible time step size errors, and simulations were performed for periods up to one year. Realistic atmospheric conditions were imposed, with an Earth's rotation rate of $7.292\xd710\u22125\u2009rad/s$ and the initial conditions defined by the imposing spherical harmonic coefficients of the vorticity expansion given by $\alpha /a$ on pre-selected modes, where *a* is the Earth radius (set to $6371.22\u2009\u2009km$) and *α* will be specified within the range of 0–30, depending on the experiment. This setup ensures that realistic atmospheric flows are generated, in terms of intensity and complexity, of what is representable with the barotropic vorticity model. No dissipation was used in these experiments, since the experiments were performed in a transient regime (i.e., before a statistically stationary state is reached). Therefore, the time elapsed was not enough for the energy to reach small scales, which guarantees the stability of the system in the timescale of our analysis. Additionally, the high-order spectral truncation, along with the small time step size, ensures an accurate representation of energy and enstrophy conservations, with variations of the order of less than $10\u22124%\u221210\u22126%$ within an integration period of a year.

### B. Numerical results

We have performed numerical experiments solving the barotropic vorticity equation initialized with equal amplitudes in the vorticity field in the modes ${(5,4),(3,1),(7,3)}$. The aim of the experiment is to understand how the energy and enstrophy redistribute themselves among neighboring modes via nonlinear interactions, and in particular, the role played by the precession resonance mechanism in the redistribution of energy/enstrophy.

Similarly to the study of reduced ordinary differential equation (ODE) models of Sec. III, we use a parameter *α* to re-scale the amplitudes (energy) of the initial condition. Figure 7 illustrates the evolution of the vorticity field, initially (t = 0) with energy only in ${(5,4),(3,1),(7,3)}$ modes, and gradually loosing coherence due to the leakage of energy to other modes for which the amplitude was initially zero. Integrations of this model are performed for a period of one year. Unlike the case of the reduced models of Sec. III, in the full spectral equations (8), part of the energy that leaves the initial triad will not return; as such, we must be aware that we might have to deal with transient behavior in the growth or decay of some of the modes. Integrations show that the mode with wavenumber ${(n,m)=(9,2)}$ is the one to receive the largest amount of energy in the range of *α*'s studied (see Fig. 9). For this reason, we choose to focus on mode $(n,m)=(9,2)$ in evaluating the efficiency of the energy/enstrophy transfer as well as the spectral characteristics of the excited oscillations.

Figure 8 shows the efficiency of the energy/enstrophy transfer from the initial triad ${(5,4),(3,1),(7,3)}$ to mode $(9,2)$. Both efficiency plots have similar shape, with a peak around *α* = 20. The values of the efficiency are small compared with the results from the reduced four-wave system in Sec. III A where the same modes were studied. The reason for this is that in the full spectral system (8), the energy has multiple paths to leak from the initial triad ${(5,4),(3,1),(7,3)}$, exciting several modes, whereas in the reduced system (31), all the energy leaving the initial triad is directed toward mode (9, 2). The power in the low frequencies is evaluated in Fig. 8(c), where we calculate the spectral power in frequencies $1/11\u2009(days)\u22121\u2264\omega \u22641/120\u2009(days)\u22121$. Unlike the reduced system analysis, we truncate the periodicities longer than 120 days to avoid the aforementioned transient behavior. The curve representing the power in the low frequencies closely follows the energy/enstrophy transfer efficiency curves, suggesting that the regime of highly efficient energy transfer also generated low-frequency oscillations in the target mode. Finally, Fig. 8(d) shows the precession frequency of the primary triad, $\Omega 1,23$, and the precession frequency $\Omega 1,42$ of the triad comprising modes ${(3,1)(7,3),(9,2)}$. First, the precession frequency of the primary triad $\Omega 1,23$ enters a phase alignment regime at around *α* = 10. The precession frequency of the secondary (target) triad ($\Omega 1,42$) decreases toward a plateau at around *α* = 13, which remains up to $\alpha \u224825$. The decrease in the precession frequency $\Omega 1,42$ therefore coincides with the amplification in the low frequencies of the system, an evidence that the dynamics of waves' phases and the precession resonance mechanism do play a role here.

Figure 9 shows examples of integrations of this system for three values of the control parameter, $\alpha =12,20,28$, displaying the time evolution of the mode energies as well as the respective spectra for the target mode (9, 2). Modes shown in the graphs are the ones that reach a critical value of the energy of at least 0.1% of the total energy of the system. As mentioned before, (9, 2) remains as the most efficiently excited mode. Other modes such as (5, 1), (3, 2), and (5, 4) also receive a substantial amount of energy through secondary interactions with the triad ${(3,1)(7,3),(9,2)}$. Paths in spectral space for the energy to reach those modes are illustrated in Fig. 1. The spectral analysis of the target mode (9, 2) shows an increase in the power in the low frequencies in the peak of the efficiency, *α* = 20. This corroborates with the finding described in Sec. III, for the reduced wave systems, in that precession resonance tends to excite lower frequency oscillations.

By the way that the aforementioned experiments were designed, energy transfers via the five-wave type of system are not allowed. From Eq. (28), if the energy is initialized only in modes (3, 4, 5), the modes (1, 2) will remain zero forever. In order to test the feasibility of the five-wave type mechanism in the full barotropic vorticity equation, we initialized the energy in the base triad ${(3,1)(7,3),(5,4)}$, plus a small perturbation in the modes (3, 0) and (7, 4). These modes are the same ones described in Sec. III B. Experiments were performed by controlling the amount of initial energy in modes ${(3,1),(7,3),(5,4)}$ and keeping the energy in modes ${(3,0)(7,4)}$ constant. The results (not displayed here) show no significant energy transfer to modes ${(3,0),(7,4)}$, while the same modes described in the previous experiment keep being excited. This suggests that the four-wave mechanism is prevalent in comparison with the five-wave mechanism in the full spectral system. A possible reason for this is that the excited modes in four-wave arrangements are coupled to the initial triads via two modes; in contrast, in five-wave arrangements they are coupled via only one mode.

In a third experiment, we initialize the energy in the modes ${(5,4),(3,1),(7,3)}$, but include perturbations in the initial condition to all modes such that $2\u2264n\u22647$. All of these modes, except those in the main interacting triad ${(5,4),(3,1),(7,3)}$, are initialized with a small initial condition defined with *α* = 1. The main purpose of this experiment is to test the robustness of the results described above. Energy time-series of modes reaching at least 0.1% of the total energy are shown in Fig. 10 for *α* = 20. Again, we see that mode (9, 2) is among the ones to receive the largest portion of the energy, although mode (5, 5) initially receives more energy in the beginning of the experiment and mode (5, 1) grows to be the one with most energy outside the initial triad toward the end of one year of experiment. Unsurprisingly, more modes receive a substantial amount of energy, 0.1%, since even low levels of initial energy in other modes will allow for more possible nonlinear interactions.

Efficiency plots for transfers toward mode (9, 2), shown in Fig. 11, qualitatively keep the pattern of the energy/enstrophy transfer efficiency without initial perturbations in Fig. 8. The same is true for the low-frequency power spectrum, which exhibits a peak coinciding with the energy/enstrophy transfer efficiency peaks. Finally, the precession frequencies $\Omega 1,23$, for modes ${(3,1)(7,3),(5,4)}$, and $\Omega 1,42$, for the triad ${(3,1)(7,3),(9,2)}$, are shown in Fig. 11(d). The overall behavior of the precession frequency $\Omega 1,23$ is very similar to the one displayed in Fig. 8, with the system entering a phase alignment regime for $\alpha \u22489$. The precession frequency $\Omega 1,42$ is initially different from the one shown in Fig. 8, starting at a smaller frequency of $\u22480.3\xd710\u22124\u2009Hz$ and subsequently increasing to a value of $\u22480.5\xd710\u22124\u2009Hz$ and plateauing at this level from $\alpha \u224810$ up to $\alpha \u224820$, after which the precession frequency increases similarly to what is observed in the unperturbed system. Since in this case we have energy in more modes, the observed characteristics of $\Omega 1,42$ is possibly due to nonlinear contributions of other neighboring triads excluding the initial triad (${(3,1)(7,3),(5,4)}$). It is important to mention that any quantitative comparison between energy transfer efficiencies in the reduced and full spectral models should be looked at with caution. In fact, unlike the reduced spectral models in which the energy of the base triad can be transferred to only one or two other modes, the full spectral model has 8384 modes, and therefore, the initial energy will have multiple ways to leak from the primary wave triad in this case. Therefore, the energy transfer efficiencies in the full barotropic vorticity equation should be expected to be much lower. However, the signatures of the precession resonance mechanism are still present in the full spectral model simulations, as can be demonstrated by the matching between the peaks of the energy/enstrophy transfer efficiency, the low-frequency power spectra, and the phase dynamics.

## V. DISCUSSION

In the present article, we studied the precession resonance mechanism in the barotropic nondivergent model. This model describes the dynamics of a two-dimensional nondivergent flow in a rotating frame and is well known to describe the essential features of the large- and planetary-scale motions of the atmosphere and ocean. The precession resonance was recently demonstrated by Bustamante *et al.*^{33} in the drift wave context as an efficient mechanism of energy transfer between different interacting triads and enhancing the cascade process responsible for the spectral energy redistribution. It consists of a balance between the nonlinear frequency associated with the energy (amplitude) modulation of an interacting wave triad and the frequency associated with the dynamical evolution of the relative phase among the components of an adjacent wave triplet. Another aspect of the precession resonance refers to its ability to yield low-frequency fluctuations in a wave system. This aspect was studied by Raphaldini *et al.*^{37} in the context of MHD Rossby waves in the solar tachocline and was proposed to be a possible mechanism for the observed long-term modulations in the main solar cycle.

To analyze the precession resonance, we first considered two reduced spectral models of the barotropic vorticity equation consisting of five- and four-wave modes describing two interacting triads coupled by one and two common modes, respectively, since these reduced ODE systems represent the possible ways to connect different nonlinearly interacting wave triplets. In the numerical integrations of these reduced spectral models, we considered a primary triad containing almost the total initial energy, and we set the initial amplitudes of the primary triad to depend on a free parameter, whose values were allowed to vary within realistic velocity magnitudes for the atmospheric flow. From these simulations, we computed the efficiency of energy and enstrophy transfers from the primary triad to the target modes of the corresponding secondary triplet. We showed that, for both five- and four-wave systems, the precession resonance is associated with the peak of both the energy and enstrophy transfer efficiencies. In addition, we computed the corresponding low-frequency range ($\omega <1/10$ $days\u22121$) power spectra associated with the time evolution of the mode energies and showed that the maxima of the low-frequency power spectrum occur for the same range of initial conditions in which the maxima of the energy/enstrophy transfer efficiency occur. This range refers to the initial amplitudes in which the energy modulation frequency of the primary triad approximately equals the mismatch among the linear eigenfrequencies of the secondary triad. By computing the fully nonlinear precession frequencies of the two-wave triplets, it was also evident that abrupt changes in these frequencies, as a function of the initial condition free parameter, often determine the beginning/end of regimes of high-energy transfer efficiency. This latter aspect of our results suggests the importance of the phase alignment for the energy transfers between two different wave triplets. The phase alignment constitutes a special case of the precession resonance referring to *p *=* *0 in Eq. (27) and is characterized by a synchronization between the wave phases within an interacting triad.

In the representative examples for the numerical integrations of the reduced spectral systems described above, we considered the modes ${(3,1),(7,3),(5,4)}$ as the primary triad that contains the initial energy of the system. As discussed before, the pump mode of this nearly resonant triad, the Rossby–Haurwitz wave mode (5, 4), was found in the literature to be unstable, and the Rossby–Haurwitz modes (3, 1) and (7, 3) account for the leading destabilizing disturbances, as demonstrated by the numerical simulations of Lynch.^{24} This base triad was connected to modes ${(3,0),(7,4)}$, in the five-wave system simulations, and to mode (9, 2), in the four-wave system numerical simulations.

To test the robustness of the precession resonance mechanism, we also performed numerical integrations with the full barotropic vorticity equation (1) by considering a regular triangular truncation with a 128° spectral resolution. As in the reduced systems, we initialized the model with the total energy on the base triad ${(5,4),(3,1),(7,3)}$, with the initial condition's energy used as a control parameter. Fluxes of energy/enstrophy to adjacent triads were evaluated as a function of the initial energy. Among the modes outside the base triad, the mode (9, 2) was found to be the most efficiently excited, as also verified in the four-wave system integrations. The efficiency of the energy/enstrophy transfer was studied and peaks of efficiency were found, an evidence that the precession resonance mechanism is in operation in the full spectral model as well. The spectral power in the low frequency was also evaluated, and similarly to the case of four-/five-wave systems, the low frequencies are amplified near the regime of efficient energy/enstrophy transfer. The robustness of the results was tested by initializing the system with a small amount of energy in the modes outside the initial triad ${(5,4),(3,1),(7,3)}$. The results previously reported were found to be robust, with mode (9, 2) remaining as the most relevant one. Quantitative differences were found between the reduced four-wave spectral model and the inter-triad interactions in the full barotropic vorticity equation. On the reduced model, the efficiency of the energy transfer to mode (9, 2) was $\u224810%$, while in the full barotropic vorticity equation the efficiency was $\u22481.5%$. This difference results from the fact that in the full PDE model, the energy has multiple ways of leaving the triad ${(5,4),(3,1),(7,3)}$, while in the reduced model it can only interact with mode (9, 2). Furthermore, the dynamics of the triad ${(9,2),(3,1),(7,3)}$ will also be affected by its interactions with other nearby interacting triads.

Therefore, the results of both the reduced and full spectral models of the barotropic vorticity equation, with parameters relevant for the Earth's atmospheric flow, demonstrate that the precession resonance mechanism enhances the low-frequency oscillation spectrum, along with promoting the maximum efficiency of energy transfer between different interacting wave triplets. Thus, based on the results, we conjecture that resonant triads might be crucial for the initial pump wave instability of a Rossby mode preferentially excited by forcing mechanisms such as baroclinic instability and thermal/orographic forcings. By contrast, the precession resonance involving two triads coupled by two common modes might play an essential role in the energy leakage from the primary triad that initiates the cascade processes associated with the spectral energy redistribution. In addition, the energy transfers between different wave triplets yield modulations in the energy (amplitude) oscillations associated with the energy exchanges within a primary triad, resulting in an amplification of the low-frequency power spectrum of the atmospheric flow. Although we also solved a full PDE model of the atmosphere, the particular setting of the simulation is still idealized, since in order to evaluate the precession resonance in triad pairs within the fully spectral model, we considered controlled initial conditions initiated in only one triad so that we could compute the efficiency of the inter-triad energy transfer. In a more realistic setting, for instance by initializing the system with a realistic spectrum, the effect of the precession resonance would be mixed within the full overall dynamical interactions and, therefore, should be evaluated using statistical measures of phase organizations, such as the order parameters introduced in Murray and Bustamante^{34} and in Arguedas-Leiva *et al.,*^{35} instead of tracking individual triad phases. However, this goes beyond the scope of this work and is intended to be addressed in a future work.

In the atmospheric flow, the so-called low-frequency variability refers to phenomena whose timescales are longer than 10 days and, therefore, is responsible for modulating the weather systems and consequently yielding climate anomalies. Although much progress has been achieved in the last decades in the understanding of the atmospheric low-frequency variability, many points behind its generation mechanisms remain as an open question in atmospheric sciences. In this direction, there has been some research effort pointing out the role of nonlinearity in the dynamics of some phenomena within this variability timescale. For example, resonant triad interactions have been evoked to explain some features of intraseasonal oscillations;^{38,39} the multiscale transfer of heat and momentum from synoptic to planetary scales has shown to explain the morphology of the Madden–Julian oscillation (MJO),^{63} the main tropical component of the intraseasonal variability; the upscale cascade related to the so-called Rossby wave breaking^{64} is believed to be important for the development of the annular modes.^{65,66} In this context, the mechanism studied here with a simplified model of the atmospheric flow was demonstrated to maximize both the cascade process that takes place after Rossby wave instability and the long-term modulation of intra-triad energy exchanges, the two types of dynamical processes that have been evoked in the literature to explain some of the low-frequency variability modes of the atmospheric circulation.

However, it is important to mention that, as we have adopted here a simplified model of the atmospheric flow, we can only explain the essence of the processes that may potentially yield the low-frequency variability modes described above. The modeling of these phenomena, such as the MJO, goes far beyond the scope of this study as it involves some complexities not considered in our model, such as other wave types (e.g., gravity and Kelvin waves), as well as other vertical modes such as the baroclinic mode associated with the circulation response to deep convection heating. In order to evaluate the relevance of the precession resonance mechanism in this context, the results presented here need to be generalized to more complex models of the atmosphere, such as the primitive equations taking into account triads involving different wave types and vertical modes. Another important issue to be investigated refers to the formation of vortices and coherent structures characteristic of geophysical flows.^{67} Preliminary work indicates that the precession resonance mechanism favors the formation of coherent vortices. This fact can be interpreted as the vortices being a result of the superposition of several modes with coherently interacting phases. This interplay between the nonlinear phase dynamics and coherent structures is also investigated in other contexts, such as astrophysical flows.^{32}

The onset of the precession resonance and the phase alignment may also be investigated in the observational reanalysis data. Indeed, reanalysis datasets may be decomposed using normal mode function projections that provide time-series for amplitudes and phases of different wave modes associated with low-frequency atmospheric oscillations, such as the Quasi-Biennial Oscillation (QBO) and the MJO.^{68–70} Raphaldini *et al.*^{71} provided evidence for the interaction of atmospheric waves that contribute to the MJO with those that contribute to the QBO. In addition, observed wavenumber-frequency spectra obtained by a normal mode projection approach^{53} have shown a spectral peak associated with a barotropic Rossby mode compatible with the Rossby–Haurwitz mode (5, 4) studied here. Consequently, this approach could, in principle, be applied for the identification of relevant interactions among barotropic Rossby modes. Events of strong energy transfers in spectral space could also be investigated from this observational normal-mode perspective, for instance, the role of breaking Rossby waves on the disruption of the QBO,^{72,73} extreme weather events,^{74,75} among other problems related to the circulation of the atmosphere.^{76} The association between precession resonance/phase organization and events of strong energy transfer could, therefore, be verified in this setting.

## ACKNOWLEDGMENTS

The work presented here was supported by FAPESP (Fundação de Amparo à Pesquisa do Estado de São Paulo), Grant Nos. 2015/50686–1, 2017/23417–5, 2016/18445–7, 2020/14162–6, and 2021/06176–0, and Coordenação de Aperfeiçoamento de Pessoal de Nível Superior-Brasil (CAPES)-Finance Code 001. We also thank the three anonymous reviewers for their valuable comments that helped to improve this manuscript.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Breno Raphaldini:** Conceptualization (lead); Formal analysis (lead); Writing – original draft (lead); Writing – review and editing (lead). **Pedro Silva Peixoto:** Formal analysis (equal); Methodology (equal); Writing – review and editing (equal). **Andre S. Teruya:** Formal analysis (equal); Methodology (equal). **Carlos Frederico Raupp:** Writing – original draft (supporting); Writing – review and editing (equal). **Miguel D. Bustamante:** Conceptualization (equal); Writing – review and editing (equal).

## DATA AVAILABILITY

Data sharing is not applicable to this article as no new data were created or analyzed in this study.

### APPENDIX A: CONSERVED QUANTITIES AND THEIR CONSTRAINTS FOR TRIAD INTERACTIONS

In this section, we follow the study of Ripa^{40} to demonstrate how the conserved quantities (4) and (5) lead to the constraints (11) for the wave triad interactions. Substituting expansion (6) into (4) and (5) and using the orthogonality of the spherical harmonics, we obtain the following Parseval's identities:

*j*indicates each spherical harmonic, viz., $j=(nj,mj)$ and $\kappa j=nj(nj+1)$. On the other hand, multiplying (8) by $Aj*$ (where $*$ denotes the complex conjugate) and the corresponding equation for $Aj*$ by

*A*, it follows

_{j}Therefore, in order for equations above to be satisfied, the following relations must be satisfied for each interacting wave triplet (*j*, *k*, *l*),

### APPENDIX B: TWO-PARAMETER EVALUATION OF THE PRECESSION RESONANCE MECHANISM

In order to further evaluate the robustness of the precession resonance mechanism in the reduced spectral models introduced in Sec. III, here we analyze the precession resonance by considering different regimes of energy exchange among the waves of the primary triad. For this purpose, instead of a single parameter *α* multiplying the initial mode amplitudes of the base triad, now two of the initial amplitudes of this triad are multiplied by different parameters *α* and *β* independently in order to yield a broader range for its characteristic frequencies.

##### 1. Four-wave model

We proceed to evaluate the efficiency of the energy and enstrophy transfers as defined in Eqs. (35) and (36), respectively. The initial condition is now defined as a function of the dimensionless parameters *α* and *β*, varying from 0 to 60. The initial condition considered in our numerical simulations is given by

Results are shown in Fig. 12, where we observe that the energy/enstrophy transfer efficiencies [Figs. 12(a) and 12(b)], as well as the low-frequency power [Fig. 12(c)], are strongly correlated with the pattern of the precession frequencies of the triads displayed in Figs. 12(e) and 12(f), demonstrating the role of the precession resonance mechanism in maximizing both the inter-triad energy transfer efficiency and the low-frequency power spectrum. Figure 12(d) illustrates the shape parameter *μ* of the Jacobi elliptic functions defined according to Eq. (24). From this figure, one notices that the variation range of the parameters *α* and *β* considered here is associated with the corresponding variation range of $0<\mu <1$. This variation range of the parameter *μ* encompasses distinct energy modulation regimes for the primary triad. This can be noticed from Fig. 13, which illustrates the behavior of the Jacobian elliptic function $sn$ for three different values of the parameter *μ*.

##### 2. Five-wave model

Similarly, we perform the 2D analysis for the energy/enstrophy transfer efficiency referred to the five-wave model (28), with the efficiencies being defined by Eqs. (38) and (39). The initial condition is defined as $(A1,A2,A3,A4,A5)=(1.48\xd7exp\u2009(1.530i)\xd710\u22128,1.71\xd7exp\u2009(0.304i)\xd710\u22128,1.03\xd7exp\u2009(1.110i)\xd710\u22125,\alpha \xd71.22\xd7exp\u2009(5.060i)\xd710\u22125,\beta \xd71.17\xd7exp\u2009(3.950i)\xd710\u22125)$.

Similarly to the four-wave example, the results displayed in Fig. 14 referred to the integrations of system (28) show that the energy/enstrophy transfer efficiencies and the low-frequency power [panels (a)–(c), respectively] are largely determined by the precession frequency of the target triad [panel (f)]. In addition, the variation range of the parameters *α* and *β* corresponds to a variation range of $0<\mu <1$, which encompasses distinct regimes of energy exchanges for the primary triad ${(3,1),(7,3),(5,4)}$.

## References

^{18,40}