# Inverse Cascades of Kinetic Energy as a Source of Intrinsic Variability: A Global OGCM Study

Guillaume Sérazin, Thierry Penduff, Bernard Barnier, Jean-Marc Molines, Brian Arbic, Malte Muller, Laurent Terray

► To cite this version:

Guillaume Sérazin, Thierry Penduff, Bernard Barnier, Jean-Marc Molines, Brian Arbic, et al. Inverse Cascades of Kinetic Energy as a Source of Intrinsic Variability: A Global OGCM Study. Journal of Physical Oceanography, 2018, 48 (6), pp.1385-1408. 10.1175/JPO-D-17-0136.1. hal-03003507

HAL Id: hal-03003507

https://hal.science/hal-03003507

Submitted on 17 Nov 2021

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_1.png" /></figure>
# Inverse Cascades of Kinetic Energy as a Source of Intrinsic Variability: A Global OGCM Study

GUILLAUME SÉRAZIN

Université Grenoble Alpes, CNRS, IRD, Grenoble INP, IGE, Grenoble, and Sciences de l'Univers au CERFACS, CERFACS/CNRS, URA1857, Toulouse, France

THIERRY PENDUFF, BERNARD BARNIER, AND JEAN-MARC MOLINES

Université Grenoble Alpes, CNRS, IRD, Grenoble INP, IGE, Grenoble, France

BRIAN K. ARBIC

Department of Earth and Environmental Sciences, University of Michigan, Ann Arbor, Michigan

MALTE MÜLLER

Research and Development Department, Norwegian Meteorological Institute, Oslo, Norway

LAURENT TERRAY

Sciences de l'Univers au CERFACS, CERFACS/CNRS, URA1857, Toulouse, France

(Manuscript received 10 July 2017, in final form 9 March 2018)

## ABSTRACT

A seasonally forced 1/12° global ocean/sea ice simulation is used to characterize the spatiotemporal inverse cascade of kinetic energy (KE). Nonlinear scale interactions associated with relative vorticity advection are evaluated using cross-spectral analysis in the frequency–wavenumber domain from sea level anomaly (SLA) time series. This analysis is applied within four eddy-active midlatitude regions having large intrinsic variability spread over a wide range of scales. Over these four regions, mesoscale surface KE is shown to spontaneously cascade toward larger spatial scales—between the deformation scale and the Rhines scale—and longer time scales (possibly exceeding 10 years). Other nonlinear processes might have to be invoked to explain the longer time scales of intrinsic variability, which have a substantial surface imprint at midlatitudes. The analysis of a fully forced 1/12° hindcast shows that low-frequency and synoptic atmospheric forcing barely affects this inverse KE cascade. The inverse cascade is also at work in a 1/4° simulation, albeit with a weaker intensity, consistent with the weaker intrinsic variability found at this coarser resolution. In the midlatitude North Pacific, the spatiotemporal cascade transfers KE from high-frequency frontal Rossby waves (FRWs), probably generated by baroclinic instability, toward the lower-frequency, westward-propagating mesoscale eddy (WME) field. The WMEs provide local gradients of potential vorticity that support these short Doppler-shifted FRWs. FRWs have periods shorter than 2 months and might be subsampled by altimetric observations, perhaps explaining why the temporal inverse cascade deduced from high-resolution models and mapped altimeter products can be quite different. The nature of the nonlinear interactions between FRWs and WMEs remains unclear but might involve wave turbulence processes.

## 1. Introduction

Baroclinic and barotropic instabilities in the ocean spontaneously generate a complex mesoscale eddy field on spatial scales of order 10–100 km. This mesoscale eddy field is now resolved (at least partly) in many ocean general circulation models (OGCMs), whose increased

Supplemental information related to this paper is available at the Journals Online website: https://doi.org/10.1175/JPO-D-17-0136.s1.

Corresponding author: Guillaume Sérazin, guillaume.serazin@legos.obs-mip.fr
resolution has led to substantially improved eddy kinetic energy (EKE) levels and mean current positions relative to observations (Penduff et al. 2010; Bryan 2013). Explicitly resolved mesoscale eddies strongly interact with the general circulation, as suggested in idealized studies (e.g., Holland 1978). These eddies not only modify the spatial structure of the large-scale flow, but they also contribute to the spontaneous generation of low-frequency oceanic intrinsic variability, such as the vacillation of eastward jets (Holland and Haidvogel 1981).

Characterizing the exact nature of eddy-driven low-frequency variability is still a topic of ongoing research. Using process-oriented quasigeostrophic (QG) models in double-gyre configurations, some authors have shown that resolving the nonlinear eddy field is a necessary condition for shaping the large-scale intrinsic variability modes of eastward jets (Dewar 2003; Hogg and Blundell 2006; Berloff et al. 2007a,b). The process involved in these four studies is the rectification of a slowly varying flow by eddy potential vorticity fluxes. Other studies, based on dynamical system theory concepts, suggest that the mesoscale activity acts as a background noise, able to trigger spontaneous transitions between basin-scale modes in the system through Hopf and homoclinic bifurcations (Simonnet and Dijkstra 2002; Dijkstra and Ghil 2005; Pierini 2006, 2011). These two paradigms provide complementary views on this eddy-driven low-frequency variability, but few studies have thus far attempted to reconcile the two of them. More recently, Arbic et al. (2012, 2014, hereinafter A12 and A14) suggested another, potentially complementary, paradigm. A12 and A14 discuss the possibility that the nonlinear advection of relative vorticity (NLA) drives intrinsic oceanic variability through a “temporal inverse cascade” of mesoscale kinetic energy (KE) alongside the more familiar spatial inverse cascade.

Turbulent cascades and inertial ranges in the spatial, or wavenumber, domain have been examined for three-dimensional flows in the pioneering work of Kolmogorov (1941) and for two-dimensional turbulence in later work (Batchelor 1953; Fjørtoft 1953; Kraichnan 1967, 1971). Adding the Coriolis force and stratification in the equations leads to the problem of geostrophic turbulence, which is still a topic of active investigation. The classic theory of geostrophic turbulence (Charney 1971; Fu and Flierl 1980; Hua and Haidvogel 1986; Salmon 1998) predicts that total (kinetic plus potential) baroclinic energy undergoes a forward cascade toward larger wave-numbers, followed by a barotropization of the flow around deformation radius-like scales. The barotropic energy, fed from the baroclinic energy, then inversely cascades toward smaller wavenumbers until it is arrested by the $\beta$ effect or bottom friction. The existence of an

inverse cascade of surface oceanic KE to larger spatial scales has been highlighted using the framework of spectral KE fluxes computed in Fourier (wavenumber) space from satellite altimeter data (Scott and Wang 2005, hereinafter SW05). Assuming that the surface signal seen by altimeters principally reflects the first baroclinic mode (Wunsch 1997), SW05 interpreted this observed inverse cascade of kinetic energy as being essentially baroclinic. Using two-layer QG simulations, Scott and Arbic (2007, hereinafter SA07) argued that baroclinic KE undergoes an inverse cascade alongside the barotropic KE inverse cascade, without violating the classical predictions of a forward cascade of total baroclinic energy.

A12 and A14 applied the framework of spectral KE fluxes utilized in SW05 and SA07 to the frequency and frequency-wavenumber domain analysis of ocean simulations. A related frequency-wavenumber approach was taken in earlier atmospheric studies, which employed cross-spectral analyses to characterize nonlinear interactions due to momentum advection (Hayashi 1980, 1982; Sheng and Hayashi 1990a,b). The frequency-wavenumber framework is useful for investigating the interactions made by NLA between high- and low-frequency motions. A12 and A14 demonstrated the existence of a spatiotemporal inverse cascade of surface KE in both QG models and in strongly eddying regions of realistic-domain high-resolution ocean models. However, A12 and A14's frequency and frequency-wavenumber analyses of gridded satellite altimetry products yielded less clear results, likely because of the relatively coarse spatial and temporal sampling of satellite altimeter data.

The nature and characteristics of the spatiotemporal cascade are still poorly known. In particular, the realistic-domain model outputs used in A12 and A14 were limited to a few years' duration and therefore were not long enough to quantify the contribution of the spatiotemporal inverse cascade beyond a few years. The role of the spatiotemporal inverse cascade in the spontaneous generation of low-frequency intrinsic variability, which imprints sea level anomalies (SLAs) on a wide range of space and time scales (Penduff et al. 2011; Sérazin et al. 2015, hereinafter S15), also needs to be diagnosed. S15 indeed suggest that the spatiotemporal inverse cascade might explain the coincidence between regions of intense mesoscale activity and regions where the turbulent ocean spontaneously generates substantial intrinsic low-frequency (>1.5 yr) small-scale (<6°) SLA variance. Here, we examine the validity of this hypothesis using much longer simulations than A12 and A14 did in order to address the contributions of NLA at longer time scales. We also examine the impact of model horizontal resolution, an important topic given that eddy-permitting ocean models
<figure><img src="image_2.png" /><figcaption>FIG. 1. Diagram describing the intercomparisons among the three simulations (dashed boxes) using time series of various lengths (e.g., “47 yr”) and temporal sampling intervals (“5d” for “5 days,” “1d” for “1 day,” and “1m” for “1 month”). See text for description of *T* and *I* experiments. *KE(ω)* and *Π<sub>KE</sub>(ω)*, respectively, denote the spectra and spectral kinetic energy fluxes in the frequency domain (see text).</figcaption></figure>

(~1/4° resolution) are being implemented in coupled climate projection models.

In this study, we argue that the spatiotemporal inverse
cascade process in the ocean involves a transfer of KE, in
both space and time, between two quasigeostrophic
dynamical regimes, characterized by two ranges of
temporal and spatial scales that coexist in realistic ed-
dying models. We thereby address the following ques-
tions: What are the space and time scales that are
involved in the transfer of KE by NLA? Does this KE
cascade spontaneously occur in the ocean, and is it
sensitive to low-frequency and synoptic atmospheric
forcing? What is the impact of the model resolution on
the KE cascade? To address these questions, we apply
the spectral framework of A12 and A14 to four mid-
latitude regions within 50- to 70-yr simulations of a
global eddying OGCM. We compare the results of two
1/12° OGCM simulations that differ only in the fre-
quency content of their atmospheric forcing in order to
investigate the impact of the low-frequency and synoptic
atmospheric variability. To assess the impact of model
resolution on the nonlinear fluxes of KE, results from
the 1/12° simulation are compared to results from a 1/4°
simulation.

Section 2 presents the numerical simulations used in
this study and describes our regions of interest, as well as
the spectral diagnostics we use. Section 3 discusses the
intrinsic nature of the temporal inverse cascade and
the associated arrest time scales. Section 4 describes the
wavenumber distribution of KE spectral transfers and
suggests a spatial arrest scale. Section 5 focuses on the
midlatitude North Pacific and distinguishes between two
dynamical regimes of the flow: westward-propagating
eddies (WMEs) and frontal Rossby waves (FRWs).
Section 6 discusses the role of those two regimes and of

the KE inverse cascades in the energetic pathways.
Section 7 investigates the impact of model resolution
and the existence of the spatiotemporal cascade in a 1/4°
simulation. Our conclusions are given in section 8.

2. Simulations and methodology

a. Numerical simulations and regions of interest

Figure 1 summarizes the model outputs, taken from three OGCM experiments (dashed boxes), that we used in this study, as well as the intercomparisons (arrows) made in this study. To isolate the intrinsic variability spontaneously generated by the eddying ocean, we performed global climatological 1/4° and 1/12° simulations driven by a repeated mean annual cycle (i.e., with no interannual or synoptic forcing), denoted here as *I* experiments. The 1/12° *I* experiment is taken here as a reference to characterize the spontaneous inverse cascades of KE, and it is compared to the 1/4° *I* experiment in section 7 to assess the impact of resolution on the temporal cascade. In addition, we performed a 1/12° hindcast forced by the full range of atmospheric time scales, thus providing a link to the observed ocean. This fully forced 1/12° simulation is denoted as the *T* experiment and is able to reproduce the low-frequency variability of sea level observed by altimetry with high accuracy (see S15’s Fig. 1). The *T* experiment is used in section 3 to investigate how the full atmospheric variability influences the transfer of KE by NLA. The interannual-to-decadal variabilities of SLA in the 1/12° simulations are quantified in terms of standard deviation in Fig. 2. As already shown in Penduff et al. (2011) and S15, the spontaneously generated SLA intrinsic variability (*I* experiment; bottom panel) is very comparable in intensity to the fully forced variability (*T* experiment;
<figure><img src="image_3.png" /><figcaption>FIG. 2. Standard deviation of low-frequency (periods >1.5 year) SLAs from the (top) 1/12° *T* experiment (fully forced variability) and the (bottom) 1/12° *I* experiment (intrinsic variability). Monthly mean SLA fields have been low passed by applying a linear convolution with a temporal Lanczos filter. The black boxes correspond to the four midlatitude regions where the spectral analysis is applied.</figcaption></figure>

top panel) in eddy-active regions, especially in the western boundary currents (WBCs) and in the Antarctic Circumpolar Current (ACC).

The simulations used in this study are all part of the Drakkar<fnref n="1" /> project and all utilize the Nucleus for European Modelling of the Ocean (NEMO; Madec 2008) ocean/sea ice model. The ORCA12 configuration (Molines et al. 2014) was implemented with NEMO version 3.4 to perform the 1/12° *I* and *T* experiments, while the 1/4° *I* experiment was made using the ORCA025 configuration with NEMO version 2.3. The Drakkar forcing set DFS4.4 (Dussin and Barnier 2013) was used to force the 1/12° *T* experiment over the period 1958–2012. DFS4.4 is based on monthly precipitations and daily radiative heat fluxes from satellite observations and on atmospheric variables

(6-hourly 10-m air temperature, humidity, and winds) from ERA-40 before 31 December 2001 and from ERA-Interim afterward. The method used to derive the climatological forcing of the 1/12° *I* experiment is described in Penduff et al. (2011) and ensures that the forcing functions yield very similar mean states in the *I* and *T* experiments (also see Grégorio et al. 2015). The climatological forcing of the 1/4° *I* experiment was derived from DFS4 (Brodeau et al. 2010), which barely differs from DFS4.4. The 1/12° and 1/4° *I* experiments were integrated over 85 and 327 years, respectively.

The three simulations are discretized over 46 vertical levels with a partial cell representation of topography. The momentum advection scheme conserves energy and enstrophy (Barnier et al. 2006; Penduff et al. 2007; Le Sommer et al. 2009). A total variance diminishing (TVD) tracer advection scheme is used with an isopycnal Laplacian diffusion operator. The vertical mixing is parameterized by

<footnote n="1">¹ https://www.dakkar-ocean.eu/.</footnote>
<table><caption>TABLE 1. Characteristics of the oceanic regions where the spectral analysis is performed.</caption><thead><tr><th>Region of study</th><th>Midlatitude North Pacific</th><th>Gulf Stream</th><th>Agulhas</th><th>South Pacific ACC</th></tr></thead><tbody><tr><td>Longitudinal extension</td><td>145°E–130°W</td><td>69°–39°W</td><td>9°–67°E</td><td>165°–80°W</td></tr><tr><td>Latitudinal extension</td><td>27°–42°N</td><td>33°–44°N</td><td>45°–35°S</td><td>66°–46°S</td></tr><tr><td>Number of grid points (N<sub>x</sub> × N<sub>y</sub>)</td><td>1000 × 234</td><td>450 × 156</td><td>696 × 170</td><td>1022 × 464</td></tr><tr><td>Deformation radius (L<sub>D</sub>)</td><td>212 km</td><td>167 km</td><td>161 km</td><td>94 km</td></tr><tr><td>Rhines scale (L<sub>R</sub>)</td><td>718 km</td><td>838 km</td><td>853 km</td><td>635 km</td></tr><tr><td>Turbulent β scale (L<sub>β</sub>)</td><td>282 km</td><td>387 km</td><td>401 km</td><td>341 km</td></tr><tr><td>Zonostrophy index (R<sub>β</sub> = L<sub>R</sub>/L<sub>β</sub>)</td><td>2.5</td><td>2.2</td><td>2.1</td><td>1.9</td></tr></tbody></table>

the TKE turbulent closure model (Blanke and Delecluse 1993), and the convective adjustment is handled by enhancing the vertical mixing in case of static instability. More details concerning the 1/12° simulations may be found in Molines et al. (2014), as well as in S15 and Grégorio et al. (2015), who used the same simulations.

Most of the plots of our study are made using 5-day mean outputs from the 70-yr 1/12° *I* experiment and from the 47-yr 1/12° *T* experiment in order to describe the spatiotemporal inverse cascade of KE. As illustrated in Fig. 1, analyses applied on monthly, 5-day, and daily mean outputs of the 1/12° *I* experiment are compared in appendix A. This comparison highlights the role of high-frequency motions in the inverse spatiotemporal cascade and justifies the use of 5-day mean outputs to investigate the longest time scales fed by NLA. Monthly mean outputs from the 300-yr 1/4° *I* experiment are used to evaluate the low-frequency part of the KE spectrum.

We select four midlatitude regions, illustrated by the black boxes in Fig. 2, where the low-frequency intrinsic variability has a large imprint on SLA and contributes to most of the total variability (see also S15): the mid-latitude North Pacific including the Kuroshio extension, the Gulf Stream extension, the Agulhas region, and the South Pacific ACC. The dimensions of these regions, along with other characteristics (some of which will be described later in the text), are given in Table 1. Subsets of these four regions have been used in A14, who showed that all four regions exhibit a spatiotemporal inverse cascade (toward larger spatial and temporal scales) in hindcasts of the Hybrid Coordinate Ocean Model (HYCOM). In the present study, motivated by the likely association of long temporal scales with larger spatial scales, we extend the Kuroshio and the South Pacific ACC boxes in the zonal direction. The extension of such regions is, however, limited by the presence of coastlines and islands that must be excluded for the spectral analysis.

## b. Spectral diagnostics

The geostrophic streamfunction $\psi(x, y, t)$, where *x* and *y* are, respectively, the East–West and North–South

coordinates, and where *t* is time, is computed from the spatially and temporally detrended sea level anomalies $\eta$ via $\psi = (g/f)\eta$, where *f* is the Coriolis parameter and *g* is the gravitational acceleration. The three-dimensional Fourier transform of $\psi$ is then taken. We use the following plane wave convention with the spatial coordinates $\mathbf{r} = (x, y)$ and the spatial wavenumbers $\mathbf{\kappa} = (k, l)$:

$$ \psi(\mathbf{r}, t) = \sum_{\omega} \sum_{\mathbf{\kappa}} \hat{\psi}(\mathbf{\kappa}, \omega) e^{i2\pi(\mathbf{\kappa} \cdot \mathbf{r} - \omega t)}, \quad (1) $$

where $\omega$ denotes frequency, and $\hat{\psi}(\mathbf{\kappa}, \omega)$ denotes the three-dimensional Fourier transform of $\psi(x, y, t)$. This convention ensures that for positive frequencies, eastward (westward) motions correspond to positive (negative) zonal wavenumbers, while northward (southward) motions correspond to positive (negative) meridional wavenumbers. Following standard practice in spectral analysis, we apply windows in both space and time to ensure periodic inputs for the Fourier transforms. We utilized Tukey windows as in A14. The KE spectrum $KE(k, l, \omega)$ can be computed through multiplication of the Fourier transform of the geostrophic streamfunction by its complex conjugate $\hat{\psi}^*$:

$$ KE(k, l, \omega) = \frac{1}{2}(k^2 + l^2)\hat{\psi}\hat{n}\hat{\psi}^*. \quad (2) $$

The transfer of kinetic energy $T_{KE}(k, l, \omega)$ associated with NLA is evaluated as in A12 and A14:

$$ T_{KE}(k, l, \omega) = \mathrm{Re} \left[ J(\psi, \nabla^2 \psi) \hat{\psi}^* \right], \quad (3) $$

where the NLA term

$$ J(\psi, \nabla^2 \psi) = \frac{\partial \psi}{\partial x} \frac{\partial (\nabla^2 \psi)}{\partial y} - \frac{\partial \psi}{\partial y} \frac{\partial (\nabla^2 \psi)}{\partial x}. $$

While KE is the variance of the geostrophic velocity, $T_{KE}$ is the cospectrum of NLA and the geostrophic streamfunction. A positive (negative) value of $T_{KE}$ characterizes a local input (output) of KE in $(k, l, \omega)$ space.
Because spectral transfers $T_{KE}$ are noisy, we prefer to present the integrals of the transfers here. The fluxes of KE across frequencies $\Pi_{KE}(\omega)$ are computed by integrating the transfers $T_{KE}$ over all wavenumbers and taking the antiderivative along the $\omega$ axis:

$$ \Pi_{KE}(\omega) = \int_{\omega \le \Omega \le \omega_N} \left[ \int_{-k_N}^{k_N} \int_{-l_N}^{l_N} T_{KE}(k, l, \Omega) dk dl \right] d\Omega, \quad (4) $$

where the Nyquist wavenumbers and frequency are, respectively, $k_N$, $l_N$, and $\omega_N$. The sign of $\Pi_{KE}$ denotes the direction of the KE cascade. A negative (positive) sign is equivalent to an inverse (direct) cascade of KE toward lower (higher) frequencies. The NLA does not create KE, but redistributes it across spatiotemporal scales through triadic interactions so that the integrated effect of the NLA on the domain is zero<fnref n="2" /> [i.e., $\Pi_{KE}(0) = 0$]. The derivative of $\Pi_{KE}(\omega)$ yields the negative of the transfer $T_{KE}(\omega)$ and thus provides information on the spectral inputs (negative slope) and outputs (positive slope) of KE averaged over all spatial scales. The definition of $\Pi_{KE}(\omega)$ implies that the shape and intensity of this function may be sensitive to the magnitude of positive and negative KE transfers at high frequencies (limited by the Nyquist frequency $\omega_N$), as shown in appendix A.

We will focus on four quantities: the frequency distribution of KE $KE(\omega)$, the KE fluxes across frequencies $\Pi_{KE}(\omega)$, the wavenumber distribution of positive and negative KE transfers $T_{KE}(k, l)$, and the wavenumber distribution of KE $KE(k, l)$. Note that $KE(\omega)$, $KE(k, l)$, and $T_{KE}(k, l)$ are further normalized by sampling frequencies before plotting power spectral densities. Certain spectra are shown in variance-preserving form $\omega KE(\omega)$ with a logarithmic $x$ axis, such that the area under the curve is proportional to the variance in particular frequency ranges. These diagnostics are performed on the NEMO irregular grid, which may introduce biases in the spatial scales. However, we verified that these biases do not yield a relative error on the location of the wavenumbers larger than 1% for spatial scales larger than 100 km and for the regions that we focus on in the following.

### 3. Temporal inverse cascade of KE

The nonlinear surface KE spectral fluxes $\Pi_{KE}(\omega)$ defined in (4) are evaluated for the $I$ and $T$ experiments over the four midlatitude regions shown in Fig. 2. Figure 3 shows that all four regions exhibit an obvious cascade of KE toward low frequencies [i.e., $\Pi_{KE}(\omega) < 0$]) in both experiments. KE is transferred from frequencies where $d\Pi_{KE}(\omega)/d\omega > 0$ toward frequencies where $d\Pi_{KE}(\omega)/d\omega < 0$. The minimum of $\Pi_{KE}(\omega)$, $\Pi_{KE}^{min}$, denoted by the thin dashed lines in Fig. 3, is found at the frequency where, averaged over wavenumbers, negative KE transfers are compensated by positive KE transfers. The time scales at which $\Pi_{KE}(\omega) = \Pi_{KE}^{min}$ are rather homogeneous among the different regions, ranging from 2.5 cpy in the South Pacific ACC to 3.5 cpy for the Gulf Stream. Thus, eddies with time scales shorter than 0.3–0.4 years continually lose KE, which then feeds motions having longer time scales. Interestingly, the minimum $\Pi_{KE}^{min}$ is associated with the most energetic time scales (maxima of blue curves, Fig. 4), except in the midlatitude North Pacific region, where $\Pi_{KE}^{min}$ occurs at shorter periods than the KE maximum. The magnitude of $\Pi_{KE}^{min}$ is not directly comparable among the different regions because it depends on the underlying dynamics.

#### a. Robustness to low-frequency and synoptic atmospheric forcing

The inverse temporal cascade spontaneously occurs under purely seasonal forcing within all four regions of interest ($I$ experiment; blue curves in Fig. 3). The full forcing used in the hindcast ($T$ experiment; green curves), which includes interannual-to-decadal and synoptic time scales, barely changes the magnitude and shape of KE fluxes across frequencies (green curves in Fig. 3), that is, the inverse temporal cascade process. The largest change in the frequency at which $\Pi_{KE} = \Pi_{KE}^{min}$ is found in the Gulf Stream region and does not exceed 10%. The atmospheric variability, therefore, has no significant influence on the spontaneous extraction of high-frequency KE by NLA. This KE subsequently feeds intrinsic variability at longer time scales in midlatitude regions, regardless of the high- and low-frequency variability of the forcing.

#### b. Arrest time scale

To quantify the range of periods associated with the spontaneous KE inverse cascade, we define the time scale $\tau_{NL}$, represented by the thin solid vertical lines in Fig. 3, as the period where the flux reaches 5% of $\Pi_{KE}^{min}$. We do not claim that the time scale $\tau_{NL}$ is an exact measure of the arrest time scale, but that it gives an

<footnote n="2">² This property is true in an idealized periodic domain. It holds in our results because the periodicity of the regions of study is forced. The signals are detrended both in space and time, and the boundaries are set to zero using a window function. Note that A14 has highlighted the importance of using a temporal detrending to guarantee that the flux $\Pi_{KE}(\omega)$ converges to 0 at the lowest frequencies.</footnote>
<figure><img src="image_4.png" /><figcaption>FIG. 3. Nonlinear surface spectral KE fluxes across frequencies Π<sub>KE</sub>(ω) in the four regions of study for the *I* (blue, intrinsic variability) and *T* (green, fully forced variability) experiments. The dashed color lines indicate the value of the minimum (deepest part of the trough) and its associated time scale. The thin solid lines are similarly set at 5% of the minimum and characterize the longest time scale at which KE is injected by NLA.</figcaption></figure>

order of magnitude of the longest time scales fed by the
inverse KE cascade.<fnref n="3" />

In the *I* experiment (intrinsic variability; thin solid blue lines), τ<sub>NL</sub> depends on the region of study: it is close to 2 years in the midlatitude North Pacific, 10 years in the Gulf Stream area, 8.5 years in the Agulhas region, and 14 years in the South Pacific ACC. In other words, the inverse cascade directly feeds intrinsic variability up to interannual time scales in the Kuroshio extension, but feeds time scales that are close to (or longer than) decadal in the other boxes. Because of the sensitivity of τ<sub>NL</sub> to dataset features, the longest time scales fed by the inverse KE cascade may well be underestimated in the South Pacific ACC; a longer 1/12° experiment might be needed to better refine our estimate. The time scales τ<sub>NL</sub> do not change by more than one-tenth of a decade with the full atmospheric forcing (*T*-experiment hindcast,

thin solid green lines), except in the Gulf Stream area,
where τNL reduces to 5 years instead of 10 years. These
results thus suggest that regardless of the atmospheric
variability, mesoscale turbulence spontaneously feeds
interannual-to-decadal intrinsic variability.

It is interesting to note that in all four regions, the KE(ω) spectra in the seasonally forced simulation (blue curves, Fig. 4) exhibit substantial multidecadal intrinsic variability (with a clear plateau in classic log-log scale, not shown here), although most of the KE is concentrated between monthly and interannual time scales (as shown by the variance-preserving spectra). The longest time scales (τNL ~ 2 years) fed by the inverse cascade in the midlatitude North Pacific box are likely to be well captured by our analysis of this 70-yr dataset: the temporal inverse cascade does not feed interannual-to-decadal intrinsic variability in this region. Nonlinear processes distinct from relative vorticity advection might, therefore, have to be invoked to explain this very low-frequency intrinsic variability. Large-scale baroclinic instability, which feeds multidecadal intrinsic variability directly from available potential energy (APE; Colin de Verdière and Huck 1999; Huck et al. 2015), is a plausible example of such a process. This

<footnote n="3"><sup>3</sup>The time scale τ<sub>NL</sub> may be sensitive to the temporal resolution of the dataset because temporal averaging may disturb the shape of the spectral fluxes (see A14’s filtering experiments in their Fig. 11). Appendix A shows that τ<sub>NL</sub> does not change by more than half a decade as temporal resolution and the duration of the datasets are altered.</footnote>
<figure><img src="image_5.png" /><figcaption>FIG. 4. As in Fig. 3, but for the frequency spectra of KE. The spectra have been multiplied by $\omega$ and are therefore in variance-preserving form. The spectra have been smoothed using a 7-point Hanning window. The thin and dashed color lines have been taken from Fig. 3.</figcaption></figure>

hypothesis is also discussed in Grégorio et al. (2015), who showed that the Atlantic meridional overturning circulation is affected by multidecadal intrinsic variability. Closing the full three-dimensional energy budget at these time scales (which would require longer simulations) would help identify other generation processes, but this lies beyond the scope of this study. The multidecadal intrinsic variability might also originate from outside the boxes, possibly involving transport by currents and waves [e.g., nonlocal eddy–mean flow interactions as in Chen et al. (2014)].

## 4. Spatial inverse cascade of KE

We now integrate $T_{KE}(k, l, \omega)$ and $KE(k, l, \omega)$ spectra over positive frequencies $\omega$ in order to investigate the spectral transfer and the KE spectrum in the spatial domain. Integrating only over positive frequencies allows the distinction between westward ($k < 0$) and eastward ($k > 0$) propagations, as well as southward ($l < 0$) and northward ($l > 0$) propagations. Because the spatial transfers discussed below spontaneously occur under a seasonal forcing and are very similar with full forcing, we focus only on the 70-yr 5-day mean $l$-experiment outputs.

### a. Spatial scales of $T_{KE}(k, l)$

Figure 5 shows the KE transfers (left) and spectra (right) as a function of zonal and meridional wavenumbers ($k, l$) in the midlatitude North Pacific (top) and Gulf Stream (bottom) boxes. These regions share striking similarities: blue areas in left panels show that negative KE transfers are concentrated in eastward motions (tilted slightly poleward), with spatial scales ranging between 1 and 2–3 times the box-averaged deformation scales $L_D$ associated with the first baroclinic mode ($L_D$ is represented by dashed contours and locally computed as $L_D = (2/|f|) \int_{-H}^{0} N(z) dz$, where $H$ is the bottom depth and $N$ is the Brunt–Väisälä frequency).<fnref n="4" /> Red areas show that positive KE transfers are concentrated in westward motions (tilted slightly equatorward) with larger spatial scales, close to those of mesoscale eddies (wavelength $\lambda \sim 500$ km, twice the diameter of an eddy). Because NLA only redistributes KE across scales via triadic interactions, such distributions of negative and positive KE transfers imply a net inverse cascade of KE from small

<footnote n="4">⁴ Our definition of $L_D$ is similar to (2.4) of Chelton et al. (1998), but it has been multiplied by $2\pi$ to scale as a spatial wavelength. This convention makes the comparison possible between $L_D$ and spectra of Figs. 5 and 6, plotted in cycles per kilometer.</footnote>
<figure><img src="image_6.png" /><figcaption>FIG. 5. Nonlinear spectral transfer $T_{KE}(k, l)$ (left) and spectra $KE(k, l)$ (right) in the zonal and meridional wavenumber $(k, l)$ space from the 1/12° *I* experiment in the (top) midlatitude North Pacific region and (bottom) Gulf Stream region. Two contour values from the transfer field, with isovalues of 0.4/2 (red) and -0.4/-2 (blue) [(nW kg⁻¹)cpk⁻²] in the North Pacific/Gulf Stream, are added on the corresponding spectra of KE. The dashed and dotted curves correspond to the internal deformation scale and the Rhines scale, respectively. All results are taken from the 1/12° *I* experiment.</figcaption></figure>

to large scales, as already noticed in previous observational and OGCM-based studies (e.g., SW05; A12; A14; Tulloch et al. 2011).

A novel feature of our analysis is the east/west asymmetry of positive and negative KE transfers and the corresponding asymmetry in KE spectra. The westward-traveling features that are fed by nonlinear interactions (red contours in the right panels of Fig. 5) also have large KE levels. Therefore, the spatial inverse cascade of KE could well explain the KE spectral asymmetry, at least in the absence of other sources of KE.

The results in the Agulhas region (top of Fig. 6) are consistent with the Gulf Stream and midlatitude North Pacific cases in several aspects: positive (negative) KE

transfers are predominantly found in westward (eastward) motions slightly tilted equatorward (poleward); positive KE transfers are found at larger scales than negative KE transfers and are superimposed with a KE maximum; and eastward-propagating negative KE transfers have scales ranging between 1 and 2–3 times the local deformation scales. Eastward-propagating eddies are associated with a particularly large KE peak, which contrasts with the midlatitude North Pacific and Gulf Stream regions in which KE is essentially concentrated in westward motions. The substantial KE of eastward motions in the Agulhas region is likely due to the advection of mesoscale eddies by the strong eastward depth-mean flow found in this region, as shown by the estimation of the zonal propagation speed
<figure><img src="image_7.png" /><figcaption>FIG. 6. As in Fig. 5, but in the (top) Agulhas region and (bottom) South Pacific ACC. The contour lines show isovalues of 2/0.1 (red) and -2/-0.05 (blue) [(nW kg⁻¹)cpk²] in the Agulhas/ACC.</figcaption></figure>

of nonlinear eddies from satellite observations (Klocker and Marshall 2014). The Agulhas region indeed includes one branch of the ACC, which may act as a waveguide, pushing Rossby waves eastward with a Doppler shift induced by the depth-mean flow (Hughes 1996). Such a Doppler shift might also occur locally in the Gulf Stream and the Kuroshio (Klocker and Marshall 2014), but westward propagation seems to dominate the KE spectrum in larger regions around those currents (right panels of Fig. 5).

In the South Pacific ACC (bottom of Fig. 6), the spatial inverse cascade also transfers KE from 1 to 2–3 times the deformation scales to larger scales. However, other features of the South Pacific ACC are somewhat different from the three WBC regions. Positive and negative KE transfers are more symmetrically distributed along the east–west axis, and the South Pacific ACC spectrum exhibits a single, strong KE

maximum in mesoscale eastward motions (likely due to strong advection by the depth-mean flow, as mentioned above). The negative KE transfers (blue, $L \le 250$ km) mostly take place within these strong eastward-propagating eddies, and the positive KE transfers (red, $L \ge 250$ km) feed both strong eastward-propagating motions and weaker westward-propagating motions. The spatial inverse cascade is thus also at work in this region, but it is possible that the westward-propagating large-scale features that are fed by the process are more strongly dissipated than elsewhere. Our results in the South Pacific ACC region must be considered with caution because the injection of KE by baroclinic instability might not be well represented, as the model barely resolves the deformation scale in this region. Thus, the unstable baroclinic modes that develop in the model ACC might be somewhat larger than in the real
ocean. Unstable modes larger than the deformation scale and different nonlinear energy transfers might also arise for a physical reason: the complex bathymetry along the ACC might enhance hydrodynamic instabilities (Abernathey and Cessi 2014) and modify energy pathways (Barthel et al. 2017). Further assessing the contribution of the spatial inverse cascade in the ACC would require further diagnostics (e.g., stability analysis, spectral diagnostic of the transfer of APE to EKE, and regional study at higher resolution), which are left for future studies.

## b. Spatial arrest scale

In the three WBC regions discussed above, the spatial inverse cascade extracts KE from eastward-propagating eddies and injects it into larger-scale, westward-propagating mesoscale eddies. The KE peaks and positive KE transfers corresponding to these mesoscale eddies are found at scales that remain shorter than the Rhines scale, defined as $L_R = 2\pi\sqrt{(2U_{\text{RMS}})/\beta}$, where $U_{\text{RMS}}$ is the root-mean-square of the absolute velocity, and $\beta$ is the planetary vorticity gradient, averaged over the boxes (dotted curves). A similar result holds for the South Pacific ACC box, where no KE is transferred by NLA toward scales larger than $L_R$. These results thus suggest that $L_R$ may well describe the spatial arrest scale of the spatial inverse cascade in the four regions studied in Figs. 5 and 6. This feature is consistent with the classical theory of waves and geostrophic turbulence established in the barotropic case (Rhines 1975): NLA cascades KE toward larger scales up to $L_R$, where the linear dynamics of Rossby waves become predominant and arrest the cascade. Using a different metric to characterize the largest scales fed by NLA,<fnref n="5" /> SW05 also discussed the Rhines scale as a potential limit of the spatial inverse cascade but did not find it relevant for high latitudes.

Finally, our analysis of the 1/12° *I* experiment also confirms that the range of spatial scales fed by the spatial inverse cascade is modest (about one decade), consistent with Tulloch et al. (2011). This range is also narrower than the range of time scales fed by the temporal inverse cascade (up to one-and-a-half decades), as already noticed by A14.

# 5. Two dynamical regimes

We now focus on the midlatitude North Pacific box and split the previous analysis into two distinct frequency bands. We integrate the KE transfers over low to medium frequencies (periods >2 months) and over high frequencies (periods <2 months). The integration is done only for positive frequencies so that the information about the propagation directions is retained. This split aims at highlighting two distinct regimes of quasigeostrophic motions, between which KE is exchanged by NLA. This cutoff frequency is chosen two times higher than the frequency associated with $\Pi_{\text{KE}}^{\text{min}}$ (~4 months; see section 3 and top-left panel in Fig. 3) so that high frequencies mainly include negative KE transfers without being polluted by low frequencies. In the following, we will show that low-frequency motions, associated with KE positive transfers, contain most of the flow KE, such that the use of a very selective filter is essential to properly isolate high-frequency motions. A movie of the SLA evolution in both frequency ranges is available in the online supplemental material. Note that the spontaneously generated SLA features at time scales shorter than 2 months were not included in S15's previous study of intrinsic variability because monthly fields were used.

## a. Westward-propagating mesoscale eddies

Low-to medium-frequency motions (periods >2 months; Fig. 7's top-right panel), corresponding primarily to westward-propagating mesoscale eddies, contribute to most of the total KE, as they contain 10 times more KE than high-frequency motions (periods <2 months; Fig. 7's bottom-right panel, different color bar scale), which are primarily eastward propagating. WMEs have spatial scales ranging between 1.5 times the deformation scale and the Rhines scale (Fig. 7's top-right panel). Their westward propagation is illustrated by the black contours in the time-longitude plot of Fig. 8 and the dashed-dotted red line that depicts a phase speed of $3 \text{ cm s}^{-1}$, thereby consistent with the observed propagation speed of oceanic mesoscale eddies (Chelton et al. 2007, 2011) and slightly larger than the phase speed of long Rossby waves in this region [about $2 \text{ cm s}^{-1}$; see also Chelton and Schlax (1996)]. WMEs are fed by NLA because they are associated with the clearly predominant positive KE transfers in the same spectral regions (red values, Fig. 7's top-left panel) at time scales longer than 2 months, larger than the negative KE transfers (blue values, same plot) in eastward motions.

## b. Downstream-propagating frontal Rossby waves

At high frequencies (periods <2 months), positive KE transfers almost vanish and give way to negative transfers, which predominate in this frequency band and

<footnote n="5">⁵ The conclusions of SW05 were based on a definition of a spatial arrest scale related to the nonlinear spectral fluxes. By analogy with the temporal case in section 3, a definition of an spatial arrest scale based on the nonlinear spectral fluxes yields results that depend on the spatial and temporal resolution inherent in the construction of satellite altimeter products (see A14's filtering experiments and their Fig. 10). Our study is based on $T_{\text{KE}}(k, l)$, whose shape at large scales may not depend on small-scale KE contrary to the nonlinear fluxes and may therefore be less sensitive to the spatial and temporal resolutions of the dataset.</footnote>
<figure><img src="image_8.png" /><figcaption>FIG. 7. (left) $T_{KE}(k, l)$ and (right) $KE(k, l)$ spectra in the midlatitude North Pacific box integrated over two different frequency bands: (top) low to medium frequencies (periods >2 months) and (bottom) high frequencies (periods <2 months). All results are taken from the 1/12° *I* experiment.</figcaption></figure>

extract KE from small, eastward-propagating motions (Fig. 7's bottom-left panel). These small-scale, high-frequency motions range between 1 and 1.5 times the deformation radius (200–300 km) and propagate downstream the slowly varying westward currents (white dashed lines in Fig. 8). In Fig. 9, we focus on SLA snapshots of a particular area of the midlatitude North Pacific region (i.e., 34°–41°N, 170°E–173°W), taken at different time intervals, where two cases of propagation are highlighted. The small-scale, high-frequency SLA A1, A2, and A3 features (colors) propagate downstream along a quasi-zonal current, characterized by a strong, slowly varying gradient of SLA (black contours). The B1 and B2 features propagate along the radial SLA gradient surrounding a mesoscale eddy. The slowly varying SLA gradients are associated with strong gradients of

potential vorticity (not shown here) that seem to support the wave propagation. When averaged over the midlatitude North Pacific region, the propagation of the small-scale, high-frequency motions is dominantly eastward, as depicted by the spectral KE distribution in Fig. 7's bottom-right panel.

Far from being a model artifact, the B1 and B2 features are similar to vortex Rossby waves, which have been theorized to explain the spiral bands in atmospheric tropical cyclones (Montgomery and Kallenbach 1997; McWilliams et al. 2003) and have been recently observed on an oceanic midlatitude vortex, west of Hawaii (Chavanne et al. 2010). Such a theory is, however, not appropriate to explain the propagation of the A1, A2, and A3 features along a zonal SLA front. To provide a physical interpretation of the small-scale,
<figure><img src="image_9.png" /><figcaption>FIG. 8. Time–longitude plot of low- and medium-frequency (periods >2 months, black contours) and high-frequency (periods <2 months, colors) components of intrinsic SLA at 37.5°N in the North Pacific box from the 1/12° I experiment, shown over periods typical of the dynamics (years 16 and 17). Red dashed–dotted lines show the westward propagation of medium- and low-frequency anomalies, and white dashed lines illustrate the eastward propagation at high frequencies. Individual structures A1, A2, and A3 are marked and may be seen also in Fig. 9.</figcaption></figure>

high-frequency features displayed in Fig. 9, we have extended the vortex Rossby wave theory to any slowly varying potential vorticity front. By analogy, we have termed these features frontal Rossby waves. We assume that a slowly varying zonal front, associated with a local gradient of relative vorticity $\partial_y \bar{\zeta}$ much larger than the planetary vorticity gradient $\beta$, may support the local propagation of these FRWs. Under this assumption, the zonal phase speed of FRWs $c_p^x$ may be derived for a slowly varying zonal flow with an equivalent barotropic model:

$$c_p^x = \bar{u} - \frac{\partial_y \bar{\zeta} + \beta + \bar{u}/R_D^2}{k^2 + l^2 + 1/R_D^2}, \quad (5)$$

where $\bar{u}$ denotes the slowly varying zonal velocity, and $R_D$ denotes the internal deformation radius of the first baroclinic mode (i.e., $R_D = L_D/2\pi$). The full derivation of (5) is detailed in appendix B. For simplicity, the impact of the vertical structure in $\bar{u}$ and $\bar{\zeta}$ is ignored in this equivalent barotropic model, which takes $\bar{u}$ and $\bar{\zeta}$ from the surface fields.

We investigate the validity of (5) along the quasi-zonal front at 37.5°N, depicted as a black dashed line in Fig. 9. We discuss the different quantities of (5) over the spatiotemporal window covering the period from 1 December 2016 to 31 January 2017 of the I experiment and the longitudes from 178°E to 175°W (i.e., the area where the eastward propagation is obvious in Fig. 8). Along

this section at 37.5°N, the meridional gradient of relative vorticity $\partial_y \bar{\zeta}$ is at least 10 times larger than the planetary vorticity gradient $\beta$ and has the same order of magnitude and the same sign as the vortex stretching term $\bar{u}/R_D^2$.

Let us look more closely at the anomaly A1, which has a zonal wavelength of 230 km and will be assumed isotropic. We compute the zonal phase speed of A1 using (5), and the result is a Doppler shift of A1 that propagates eastward with a mean phase speed of $7.7 \pm 2.3 \text{ cm s}^{-1}$. The eastward zonal flow $\bar{u}$ is about $31 \pm 13 \text{ cm s}^{-1}$ and dominates the rightmost term in (5). The computed phase speed $c_p^x$ of A1 is consistent with the phase speed inferred from the time–longitude plot in Fig. 8, which is about $8 \text{ cm s}^{-1}$. Note that such phase speeds have the same order of magnitude as the phase speed of FRWs trapped on a mesoscale eddy west of Hawaii and measured with ADCP observations (i.e., $8–9 \text{ cm s}^{-1}$; Chavanne et al. 2010). The period of the anomaly A1 computed from (5) is $35 \pm 10$ days, while a direct evaluation of the slopes in Fig. 8 gives a period of about 40 days. Thus, (5) reasonably describes the features and propagation of FRWs along the slowly varying SLA front at 37.5°N.

## 6. Discussion on energetic pathways

We now discuss the energetic pathways related to our spectral analysis of NLA and the decomposition of the North Pacific midlatitude flow into two dynamical
<figure><img src="image_10.png" /><figcaption>FIG. 9. Three snapshots of intrinsic SLA at 10-day intervals in the North Pacific. Black contours correspond to low and medium frequencies (periods >2 months), and colors correspond to high frequencies (periods <2 months). The black dashed lines correspond to the section where the time-longitude plot of Fig. 8 is taken. Individual structures A1, A2, A3, B1, and B2 are marked to highlight the downstream propagation of the high-frequency anomalies; A1, A2, and A3 may be seen also in Fig. 8.</figcaption></figure>

regimes. This discussion aims at gathering recent results from the theoretical literature of rotating and geostrophic turbulence and comparing them to our spectral results from the more complex world of OGCMs.

### a. Genesis of FRWs

In section 4, we have shown that KE is extracted by NLA from spatial scales ranging between 1 and 2–3 times the local deformation scales, corresponding to the scales of FRWs termed in section 5.

While investigating the seasonal modulation of EKE in the South Pacific Subtropical Countercurrent, Qiu et al. (2008) performed a comparable (albeit limited to the wavenumber domain) spectral analysis of altimeter data and a linear baroclinic instability analysis based on a 2.5-layer reduced gravity model. They showed that the wavenumbers over which KE is extracted by NLA coincide with those where perturbations are baroclinically unstable and are likely to produce EKE from APE. Based on linearized QG equations with a continuous stratification, Tulloch et al. (2011) performed a similar
local baroclinic instability analysis on the Ocean Comprehensible Atlas (OCCA; Forget 2010). In the Gulf Stream region, they predict two unstable wavenumber ranges with scales ranging between 1 and 2 times the internal deformation scale (see their Fig. 4's top panels). One of their wavenumber ranges corresponds to eastward motions and coincides with those where KE is extracted in our *I*-experiment outputs (Fig. 5, blue area in the bottom-left panel).

In fully nonlinear QG experiments representative of midlatitude regimes, Berloff and Kamenkovich (2013a) isolate a spectral quadrant called eastward mesoscale eddies (EMEs), which has spectral and physical properties (see their Fig. 7) similar to those of the eastward FRWs described in our OGCM simulations. In an associated paper, Berloff and Kamenkovich (2013b) performed a linear stability analysis and showed that these EMEs correspond to both baroclinic and barotropic unstable normal modes of the sheared background flow composed of alternating zonal jets. In other words, those aforementioned studies strongly suggest that the eastward FRWs in our simulations are likely generated by hydrodynamic instabilities and especially by baroclinic instability.

The study of Berloff and Kamenkovich (2013b) also demonstrates that the large-scale background flow and the alternating jets can substantially modify the linear modes computed from linear stability analysis. In section 5, we showed that the fast eastward FRWs may be associated with the slower WME fronts and that their propagation could be approximated by the Rossby-like dispersion relation (5). Performing a linear instability analysis with a mean based on the slower WMEs could be valuable for investigating whether or not the slow WME field might condition the hydrodynamic instabilities and might generate those FRWs. Such an analysis lies beyond the scope of this study and is left for the future.

## b. Genesis of WMEs

In geostrophic turbulence, the presence of substantial KE at large scales is commonly explained by an inverse cascade of KE fed by smaller scales through local eddy-eddy interactions. However, our spectral results in the North Pacific suggest that this inverse cascade might arise from nonlocal, weakly nonlinear interactions occurring between a linear wavelike regime (FRW) and a nonlinear eddy-like regime (WME).

This hypothesis is consistent with recent results from numerical and experimental setups of nonstratified rotating turbulence [see Godeferd and Moisy (2015) for a review]. These experiments have pointed out the role of linear waves, which could nonlinearly interact with the

turbulent flow through weak nonlinearities inherent in wave turbulence processes (Nazarenko 2011). In particular, Bourouiba et al. (2012) show that KE transfers at intermediate Rossby number involve direct, nonlocal energy transfers from small-scale, high-frequency, three-dimensional waves to large-scale, slow, two-dimensional vortices.

The inverse cascade in geostrophic turbulence may also be compatible with weakly nonlinear interactions, as shown by Harper et al. (2014), who applied the wave turbulence framework to the classic two-layer QG model. Using a scale separation, they suggest an update of the classic picture of geostrophic turbulence (i.e., Salmon 1998) in the case of weak nonlinearity and resonant triad interactions: the inverse KE transfers are dominated by nonlocal interactions between two small-scale baroclinic components, yielding a large-scale barotropic component.

Thus, the role of weak nonlinearities and linear processes in geostrophic turbulence remains to be clarified, but should not be underestimated. To complement our study, disentangling local and/or nonlocal interactions could be achieved by decomposing the SSH field prior to computation of KE transfers and analyzing the triadic interactions between FRWs and the WMEs. The wave turbulence theory (Nazarenko 2011) could also be a relevant framework to clarify the dynamics of oceanic turbulence. Finally, the KE transfers through NLA are only one contributor to the spectral budget of the surface KE. Whether or not WMEs are mostly fed by FRWs through NLA would require diagnosis of the contribution of other terms of the KE budget, such as the conversion of APE to KE by baroclinic instability.

## c. Spatiotemporal limits of the inverse cascade

One interesting result of section 4 is that the Rhines scale $L_R$ may correctly describe the largest spatial scales at which NLA transfers KE, consistent with the theory of Rhines (1975) and suggesting that the $\beta$ effect acts as a wavenumber barrier for the inverse cascade. However, the KE does not pile up at large zonal scales in our spectral analyses as in the freely evolving experiments performed by Vallis and Maltrud (1993; see the dumbbell shape of their Fig. 5c). Instead, the distribution of KE is mostly isotropic in $k$ and $l$ on one-half of the spectrum, corresponding to dominant eddy-like motions in our regions of study.

To better understand why KE may not reach large zonal scales in our regions of study, it is interesting to compute the turbulent $\beta$ scale, introduced by Vallis and Maltrud (1993) and defined as $L_\beta = 2\pi(\epsilon/\beta^3)^{1/5}$. In the four regions of study, the turbulent energy flux $\epsilon$ may be estimated from Fig. 3 by taking the absolute value of the
minimum flux $\Pi_{KE}^{\min}$, as done by Klocker et al. (2016), using the results from the spectral analysis of Arbic et al. (2014). The length scale $L_{\beta}$ characterizes the threshold scale at which the inverse cascade in $\beta$-plane turbulence becomes anisotropic (Sukoriansky et al. 2007). This length scale is computed for each region and is given in Table 1.

The ratio of the Rhines scale $L_R$ and the turbulent $\beta$-scale $L_\beta$ yields the zonostrophy index $R_\beta$, used to characterize different regimes of the flow in $\beta$-plane turbulence (Sukoriansky et al. 2007; Galperin et al. 2010): $R_\beta \le 1.5$ characterizes a friction-dominated regime where the nonlinear scale interactions are isotropic, and $R_\beta \ge 2.5$ characterizes a zonostrophic regime with a flow organized along alternating zonal jets with a strongly anisotropic KE spectrum. The values of $R_\beta$ in our four analysis regions are given in Table 1. These values range between 1.9 and 2.5, suggesting that the midlatitude oceanic flows lie between the two regimes described above. Thus, the flows might develop a partial zonostrophic inertial range between $L_\beta$ and $L_R$, but without giving rise to well-established zonally elongated structures. The route to dissipation may indeed be efficiently provided by bottom friction before the KE reaches large zonal scales and the KE spectral distribution remains partly isotropic. Because oceanic basins have meridional boundaries that break the conservation of enstrophy, Lacasse (2002) also suggests that basin normal modes may halt the inverse cascade rather than Rossby waves. Such basin modes might thus yield an isotropic arrest and prevent the establishment of a zonostrophic regime, explaining why KE does not pile up at large zonal scales.

In section 3, we have also seen that the inverse temporal cascade may be limited in time. Is there a physical mechanism responsible for this bounding in time? Chen (2013) showed analytically that a meridional mean flow could influence the spectral characteristics of the inverse temporal cascade. The effect of the meridional mean flow on the inverse temporal cascade was confirmed in a barotropic channel [see Exp2 in Fig. 8 of Chen and Flierl (2015)], preventing KE from cascading up to the zero frequency. A consequence of this lower bound in the inverse temporal cascade is that KE does not reach large zonal scales (linked to long time scales by the Rossby wave dispersion relation). Turbulent motions do not organize around well-defined Rhines jets, but instead may exhibit weak, low-frequency latent jets at midlatitudes, as in our OGCM experiments and in altimetric observations (Maximenko et al. 2005). This consequence is consistent with a short zonostrophic inertial range discussed above. Several other parameters have been shown to inhibit the development of persistent

zonal jets, including the intensity of bottom friction and the presence of meridional boundaries (Lacasse 2002; Berloff et al. 2009, 2011).

According to those aforementioned studies, the largest spatial and temporal scales fed by the inverse cascade may be governed by several oceanic parameters, such as the injection of turbulent energy, the basin geometry, and the intensity of the large-scale drag. How and to what extent those physical parameters control the space and arrest time scales of the inverse cascade remains, however, an open question that could be investigated in idealized QG models. Note that most of the previous arguments were based on barotropic models of the flow. The study of $\beta$-plane turbulence with substantial baroclinic components is still a topic of active research [see Berloff and Kamenkovich (2018) for a review on baroclinic jets].

### d. Global picture

To sum up, the spectral analyses from the *I* experiment show that the high-frequency, small-scale FRWs are potentially generated by hydrodynamic instabilities because their corresponding KE spectrum peaks at wavenumbers that are likely to be baroclinically unstable (i.e., around the deformation scale; Fig. 7's bottom-right panel). Their KE is then spontaneously transferred by NLA and redistributed through temporal and spatial scales to feed the slower and larger WMEs (Fig. 7, top). Whether these transfers occur locally or nonlocally in the frequency–wavenumber space remains to be clarified, though the characteristics of the FRWs and WMEs suggest that nonlocal, weakly nonlinear interactions could be at work. The spatiotemporal limits of this resulting inverse cascade may be linked to the incapacity of the turbulent flow to develop a wide zonostrophic inertial range and persistent elongated zonal jets. The Rhines scales seem to describe correctly the largest spatial scales fed by NLA in our regions of study, but the longest temporal scales reached by the cascade are not well understood. In particular, these limits might be controlled by several physical parameters governing the turbulent flow (e.g., bottom friction and rate of energy input).

## 7. Model resolution and nonlinear KE transfer

Does the temporal inverse cascade spontaneously generate low-frequency intrinsic variability in the 1/4° simulation as well? The green curves in Fig. 10's top panels indeed confirm that NLA induces a cascade of KE toward longer temporal scales in the 1/4° simulation, both in the midlatitude North Pacific box and in the South Pacific ACC box. The inverse cascade is, however,
<figure><img src="image_11.png" /><figcaption>FIG. 10. (top) Parameter Π<sub>KE</sub>(ω) and (middle, bottom) KE(ω) computed from the 1/12° (blue) and 1/4° (green) 27-yr 5-day mean outputs in the (left) midlatitude North Pacific box and (right) South Pacific ACC box; the red curves in KE(ω) seen in the bottom panels are computed from the 300-yr monthly 1/4° output in order to evaluate the low-frequency tail of the intrinsic KE spectrum. (middle) The KE spectra have been multiplied by ω and are therefore in variance-preserving form. (bottom) The KE spectra are plotted in log-log scale to highlight the spectral slopes and compare them to the ω<sup>-2</sup> power law (black lines). All KE spectra were smoothed using a 7-point Hanning window.</figcaption></figure>

weaker than in the 1/12° model (blue curves in Fig. 10's top panels) by a factor of about 2 in the midlatitude North Pacific and about 4 in the South Pacific ACC. The slopes of Π<sub>KE</sub>(ω) are steeper in the 1/12° model, indicating that, on average, positive and negative KE transfers are enhanced at higher resolution.

Frequency spectra displayed in variance preserving form in Fig. 10's middle panels indicate that finer resolution increases the strength of intrinsic variability in the midlatitude North Pacific and in the South Pacific ACC

over plotted frequencies: that is, between the high-frequency mesoscale range (periods <1 year), where the increase is the most prominent, and up to interannual time scales. The increase in KE in the 1/12° simulation is consistent with the increase in intrinsic SLA variance between the 1/4° and 1/12° simulations on time scales ranging from 2 to 18 months and from 1.5 to 20 years (see S15).

In both our 1/4° and 1/12° I experiments, the spectral shape at time scales shorter than 1 year are very similar,
suggesting that the same autoregressive processes are captured by both models but depend on the regions (Fig. 10's bottom panels). Venaille et al. (2011) used a simple stochastic model to interpret the intrinsic low-frequency variability of the transport around the Zapiola anticyclone as being induced by Reynolds stresses. In their model, the spectral slope is equal to $\omega^{-2}$, characteristic of a first-order autoregressive process, and the cutoff frequency is set by bottom friction.

In the midlatitude North Pacific, the spectral slopes are close to this first-order autoregressive process ($\omega^{-2}$, black line), whereas in the South Pacific ACC, the spectral shapes are different from a $\omega^{-2}$ slope and suggest a higher-order autoregressive process. It is, moreover, plausible that the injection of high-frequency, small-scale KE by baroclinic instability is underestimated in the 1/4° model since it barely resolves the internal deformation scale at midlatitudes. Hence, a weaker injection of KE in high-frequency, small-scale motions may therefore lead to weaker Reynolds stresses, which, integrated over time, might lead to a weaker KE over the entire frequency spectrum. Such an argument provides an explanation for the weaker inverse cascades seen in the 1/4° model, although both the 1/4° and 1/12° models capture the same physical process that transfers KE through temporal scales.

The 300-yr 1/4° simulation provides an estimate of the spectrum of intrinsic variability at very low frequencies (periods >10 years; red curve in Fig. 10's bottom panels). A low-frequency plateau is reached in this long simulation for all regions of interest and seems to coincide with the asymptotic plateau of the 1/12° simulation in the South Pacific ACC and midlatitude North Pacific boxes (Fig. 10), as well as the Gulf Stream and Agulhas boxes (not shown here). This coincidence suggests that the 1/4° simulation does resolve the main nonlinear oceanic processes generating decadal-to-multidecadal intrinsic variability.

## 8. Conclusions and discussion

We have used the framework developed in A12 and A14 to diagnose from global OGCM simulations the KE transfers induced by NLA in the frequency-wavenumber domain and to investigate the propensity of NLA to generate low-frequency oceanic intrinsic variability. Compared to these previous studies, we have used longer time series (decadal to multidecadal) and performed the analysis on both a seasonally forced simulation (that isolates the spontaneous generation of intrinsic variability) and a fully forced hindcast (that aims at reproducing the observed oceanic variability). We sum up here the main results of our study.

* NLA yields a temporal inverse cascade of KE in our 1/12° eddying OGCM and generates oceanic intrinsic variability. This process spontaneously occurs in the *I* experiment, forced by a repeated mean annual cycle, and is barely affected by the interannual-to-decadal time scales, as well as the synoptic variability of the full atmospheric forcing in the *T* experiment.

* The temporal inverse cascade feeds low-frequency variability up to interannual-to-decadal time scales (depending on the region) at midlatitudes. Multi-decadal time scales might also be fed through this process in the ACC, but longer simulations (>70 years) are required to properly assess this possibility in this region.

* Intrinsic variability at longer (decadal to multidecadal) time scales is also substantial, but it is not likely produced by the transfer of KE by NLA: other KE sources such as large-scale baroclinic instability might be at work.

* The spatial inverse cascade found in our simulations is consistent with previous studies of geostrophic turbulence, but a new feature is highlighted: an east/west asymmetry between nonlinear positive and negative KE transfers is found in the midlatitude boxes including WBCs. However, this asymmetry does not show up in the ACC box. We also find evidence that the Rhines scale describes the largest spatial scales fed by the spatial inverse cascade.

* The spatial and temporal inverse cascades of KE represent two complementary views of the same nonlinear process. In the midlatitude North Pacific, the cascades involve two dynamical regimes with different temporal and spatial scales, which interact via NLA: eastward-propagating frontal Rossby waves of deformation-scale size that predominate at high frequencies (periods <2 months) and westward-propagating mesoscale structures having longer time scales and containing most of the KE. FRWs are likely to be generated by baroclinic instability and propagate downstream along potential vorticity fronts induced by these mesoscale structures; the latter are themselves fed by FRWs through wave-eddy interactions on a limited range of spatiotemporal scales (i.e., spatial scales shorter than the Rhines scale and time scales ranging up to interannual to decadal). These mesoscale features propagate westward at a phase speed comparable to, but slightly larger than, that of long Rossby waves, consistent with previous conclusions from altimetric observations.

* The temporal inverse cascade is shown to be at work in the 1/4° model, albeit with a weaker intensity, compared to the 1/12° model. This difference might be related to the underestimation of baroclinic instability
at midlatitudes in the 1/4° model, which provides weaker high-frequency, small-scale KE available for transfer by NLA.

This study suggests that high-frequency geostrophic structures (i.e., FRWs) play a key role in transferring KE toward longer time scales and larger spatial scales. In the North Pacific midlatitudes at 37.5°N and between 178°E and 175°W, the period of FRWs is about 40 days. Higher-frequency outputs may yield somewhat shorter FRW periods than those illustrated in this study: time-longitude diagrams (not shown here), plotted in the midlatitude North Pacific box using the 10-yr daily mean 1/12° *I* experiment, show some FRWs with a period of 15 days. Such high-frequency features might be undersampled in the Archiving, Validation, and Interpretation of Satellite Oceanographic Data (AVISO) gridded products (Le Traon et al. 1998) because the spatiotemporal interpolations made on the along-track signals use an *e*-folding time scale of 15 days [see the appendix of Ducet et al. (2000)] and yield an effective Nyquist period of about 30 days in the gridded products. Subsampling these FRWs certainly hampers the detection of inverse temporal cascade processes from the observations, as demonstrated by the temporal filtering experiments in A14. A simple equivalent barotropic model was used to rationalize the (eastward) phase speed and period of simulated FRWs in the midlatitude North Pacific. FRWs and their contribution in KE pathways need, however, to be better characterized in QG and realistic models, which can provide dynamical links among baroclinic instability, eddy-mean flow interaction, and vertical structure.

An east/west asymmetry of the inverse cascades was highlighted through the use of three-dimensional Fourier transforms that allow the separation of eastward versus westward (or equatorward vs poleward) motions. We believe this spectral framework may help answer some open questions on the energy pathways in geostrophic turbulence. The study of A14 could be extended to investigate the role of the *β* effect on the inverse cascade characteristics. A separation between the barotropic and baroclinic streamfunctions, as done in Straub and Nadiga (2014) for the wavenumber spectra, could also be performed using this spectral framework. This three-dimensional Fourier transform could also serve to apply the wave turbulence framework (Nazarenko 2011) and study more precisely the nature of the interactions between FRWs and WMEs. The KE budget in open domains would also require a diagnosis of cross-boundary fluxes, as these may contribute nonlocally to the mean KE in key regions (Chen et al. 2014).

The temporal inverse cascade is not the only mechanism that generates low-frequency intrinsic variability in the ocean. More studies using a hierarchy of models, including OGCMs, are needed to better understand the other

mechanisms at work in the real ocean. In particular, the role of NLA is probably not limited to the transfer of KE across scales and may involve the nonlinear rectification of a slowly varying flow by mesoscale eddies that might modulate the position and intensity of WBCs (Hogg and Blundell 2006; Berloff et al. 2007a,b). The integration of stochastic Reynolds stresses is also a possible mechanism that may spontaneously modulate the transport (Venaille et al. 2011) and the position (Chapman and Morrow 2014) of oceanic currents close to topographic anomalies. The role of large-scale baroclinic instability in the spontaneous generation of decadal-to-multidecadal intrinsic variability in the eddying regime is also a topic of active research (Colin de Verdière and Huck 1999; Huck et al. 2015).

As already discussed in S15, the 1/4° global model generates a substantial amount of intrinsic variability and is consistent with 1/12° results regarding the basic features of the temporal inverse cascade. This motivates the use of the 1/4° model for studying oceanic intrinsic variability through fully forced ensemble simulations (Penduff et al. 2014; Bessières et al. 2017), as well as for climate studies through fully coupled global models.

*Acknowledgments.* We thank the two anonymous reviewers for their careful reading, constructive comments, and suggestions that lead to significant improvements in the manuscript. The dataset produced for this study is available upon request (contact thierry.penduff@cnrs.fr). This work is a contribution to the OCCIPUT, CHAOCEAN, and PIRATE projects. OCCIPUT (http://meom-group.github.io/projects/occiput/) is funded by the French Agence Nationale de la Recherche (ANR) through Contract ANR-13-BS06-0007-01. The CHAOCEAN (https://sealevel.jpl.nasa.gov/science/ostscienceteam/scientistlinks/scientificinvestigations2013/penduff/) and PIRATE (https://sealevel.jpl.nasa.gov/science/ostscienceteam/scientistlinks/scientificinvestigations2017/penduff/) projects are funded by the Centre National d'Etudes Spatiales (CNES) through the Ocean Surface Topography Science Team (OSTST). This work benefited from the Drakkar international coordination network (GDRI) established among the Centre National de la Recherche Scientifique (CNRS), the National Oceanography Centre in Southampton (NOCS), GEOMAR in Kiel, and IFREMER. For this work, Drakkar also benefited from a grant from the Groupe Mission Mercator Coriolis (GMMC) through the LEFE program of the Institut National des Sciences de l'Univers (INSU). The computations presented in this study were performed at the Centre Informatique National de l'Enseignement Supérieur (CINES) under the allocation made by GENCI x2013010727. The ORCA12-GJM02 simulation was performed as part of the Grands
Challenges GENCI/CINES 2013. GS is supported by CNES and Région Midi-Pyrénées; TP, BB, JMM by CNRS; and LT by CERFACS. BKA's contributions to this study were supported by U.S. National Science Foundation Grants OCE-0960820 and OCE-1351837. MM was supported by a subcontract from Office of Naval Research (ONR) Grant N00014-11-0487 to the University of Victoria.

# APPENDIX A

## Sensitivity of the Results to the Duration and Resolution of SLA Fields

This appendix compares $\Pi_{KE}(\omega)$ computed from the outputs of the 1/12° I experiment (Fig. A1) considered over different durations and temporal sampling. It discusses the limits of each dataset and issues a caveat on the use of monthly means for the spectral analysis of the KE cascade.

In the midlatitude North Pacific box, the minimum $\Pi_{KE}^{\min}$ (dashed vertical lines) is displaced toward lower frequencies (from 0.4 to 0.3 years), and its intensity is reduced by 30% when the temporal resolution is degraded from daily (blue curve) to 5-day (red and green curves). The change in intensity is smaller in the South Pacific ACC box, where $\Pi_{KE}^{\min}$ is reduced by 15% and is slightly displaced from 0.35 to 0.5 years. In both regions, the positive slopes of $\Pi_{KE}(\omega)$ that characterize frequencies at which KE is extracted by NLA are well captured at 5-day resolution. As the daily dataset includes shorter time scales, the range of frequencies at which KE is extracted becomes broader at this temporal resolution, compared to the 5-day resolution; spectral positive KE transfers also increase at monthly to interannual time scales [steeper negative slopes of $\Pi_{KE}(\omega)$], but this change may be artificially due to the temporal detrending imposing the constraint $\Pi_{KE}(0) = 0$. The time scale $\tau_{NL}$ evaluated in the 10-yr daily dataset differs by only one-tenth of a decade from the 5-day datasets.

Usage of a 27-yr record instead of a 70-yr record barely changes the intensity and the shape of $\Pi_{KE}(\omega)$ in both regions over monthly and annual time scales, probably due to the predominance of positive KE transfers at time scales shorter than 1 year. There are, however, slight differences at interannual time scales (green and red lines). These slight differences at low frequencies have a substantial effect on the time scales $\tau_{NL}$, which differ by two-tenths of a decade in the mid-latitude North Pacific and four-tenths of a decade in the South Pacific ACC. Longer time series provide more realizations of low-frequency motions, hence a better estimate of the left tail of KE spectra and fluxes. In the

ACC box in particular, 10-yr records (blue) are certainly not long enough to capture all the spectral positive KE transfers at interannual time scales.

Using monthly (cyan) instead of 5-day (green) time series strongly distorts the shape of $\Pi_{KE}(\omega)$ and eventually leads to a spurious forward cascade at low frequencies in the North Pacific box. This may be explained by the omission of high-frequency motions in the dataset from which KE is extracted. In the South Pacific ACC box, considering monthly outputs (cyan) strongly underestimates the strength of the inverse cascade (60% smaller than the 5-day mean dataset).

In conclusion, these results highlight the important role of high-frequency dynamics for diagnosing KE pathways. A long dataset is also preferable to include potential low-frequency sources of KE. A long dataset also yields a better estimate of the low-frequency tail of the spectrum. Excluding low frequencies may also change the estimate of $\tau_{NL}$ by a few tenths of a decade. Ideally, one would prefer to use a 70-yr daily mean dataset to include all time scales in the analysis. Such a global dataset is, unfortunately, unavailable to us and would require huge amounts of computer resources. The 70-yr 5-day mean record, which is motivated by the question of the longest time scales fed by the temporal inverse KE cascade, provides a good compromise between the inclusion of high-frequency dynamics and the low-frequency tail of the spectrum. We do not recommend the use of monthly mean datasets to study the temporal inverse cascade, as it strongly alters the shape of nonlinear fluxes.

# APPENDIX B

## Derivation of Frontal Rossby Wave Dispersion Relation

The conservation of potential vorticity $q$ in a flat-bottom 1.5-layer model on a $\beta$ plane is expressed as

$$ \frac{Dq}{Dt} = \frac{\partial q}{\partial t} + \mathbf{u} \cdot \nabla q = 0. \qquad (B1) $$

The potential vorticity $q$ is decomposed as

$$ q = \zeta + \beta y - \psi/R_D^2, \qquad (B2) $$

where $\zeta$ is the relative vorticity, $\psi$ is the streamfunction, $\beta$ is the planetary vorticity gradient, and $R_D$ is the internal deformation radius. We use a low-pass Lanczos filter defined as

$$ \bar{X}(t) = \frac{1}{2T} \int_{t-T}^{t+T} X(\tau)w(\tau) d\tau, \qquad (B3) $$
<figure><img src="image_12.png" /><figcaption>FIG. A1. Nonlinear surface spectral kinetic energy fluxes Π<sub>KE</sub>(ω) from the 1/12° I experiment in the (top) midlatitude North Pacific and (bottom) South Pacific ACC boxes for different temporal resolutions (daily, 5-day, monthly); 5-day outputs are also compared for time series of different durations (27 and 70 years). The dashed color lines indicate the time scale at which Π<sub>KE</sub>(ω) reaches its minimum (deepest part of the trough). The thin color lines are similarly set at 5% of the minimum and characterize the longest time scale at which KE is injected by NLA.</figcaption></figure>

where *T* is the length of the window, and *w* is the low-
pass Lanczos weighting function. Applying this low-
pass filtering on (B1) yields the slowly varying flow
equation

$$
\frac{\partial \bar{q}}{\partial t} + \bar{u} \cdot \nabla \bar{q} = - \overline{u'} \cdot \nabla q', \quad (B4)
$$

where overbars and prime, respectively, denote slowly
varying and perturbation quantities. The perturbation
equation is derived by subtracting (B4) from (B1) and
may be written as

$$
\frac{\partial q'}{\partial t} + \mathbf{u}' \cdot \nabla \bar{q} + \bar{u} \cdot \nabla q' + \mathbf{u}' \cdot \nabla q' = \overline{\mathbf{u}'} \cdot \overline{\nabla q'}. \quad (B5)
$$
Neglecting the second-order term $\mathbf{u}' \cdot \nabla q'$, and replacing $\bar{q}$ and $q'$ by their corresponding values derived from (B2), leads to

$$ \frac{\partial}{\partial t}(\xi' - \psi'/R_D^2) + \bar{\mathbf{u}} \cdot \nabla(\xi' - \psi'/R_D^2) + \mathbf{u}' \cdot (\nabla\bar{\xi} + \beta\mathbf{j} - \nabla\bar{\psi}/R_D^2) = \mathcal{F}, \quad (B6) $$

where $\mathcal{F} = \bar{\mathbf{u}}' \cdot \nabla q'$ is a slowly varying forcing term, which denotes the coupling between the fast and the slow fields, and $\mathbf{j}$ is the zonal unit vector. Injecting harmonic solutions of the form $\psi' = \bar{\psi}e^{i(kx+ly-\omega t)}$ at particular frequency $\omega$ and wavenumbers $\kappa = (k, l)$, and noting that $\xi' = -\kappa^2\psi'$ and $\mathbf{u}' = (-il\psi', ik\psi')$, finally yields the dispersion relation of FRWs:

$$ \omega = \bar{\mathbf{u}} \cdot \kappa - \frac{(k \partial_y \bar{\xi} - l \partial_x \bar{\xi}) + k\beta + \bar{\mathbf{u}} \cdot \kappa / R_D^2}{\kappa^2 + 1/R_D^2}. \quad (B7) $$

For a purely zonal flow with no meridional dependency, the propagation dispersion (B7) yields the zonal phase speed of perturbations:

$$ c_p^x = \bar{u} - \frac{\partial_y \bar{\xi} + \beta + \bar{u}/R_D^2}{\kappa^2 + 1/R_D^2}, \quad (B8) $$

where $\bar{u}$ is the slowly varying zonal velocity. If one removes the meridional gradient of relative vorticity $\partial_y \bar{\xi}$, (B8) becomes the same as the classical Rossby wave phase speed equation for a flow with a zonal mean current and a finite deformation radius [see (5.188a) in Vallis (2006)]. Equation (B8) shows that the phase speed $c_p^x$ may be modified by the local slowly varying meridional gradient of relative vorticity $\partial_y \bar{\xi}$ when it is taken into account.

## REFERENCES

* Abernathey, R., and P. Cessi, 2014: Topographic enhancement of eddy efficiency in baroclinic equilibration. *J. Phys. Oceanogr.*, **44**, 2107–2126, https://doi.org/10.1175/JPO-D-14-0014.1.

* Arbic, B. K., R. B. Scott, G. R. Flierl, A. J. Morten, J. G. Richman, and J. F. Shriver, 2012: Nonlinear cascades of surface oceanic geostrophic kinetic energy in the frequency domain. *J. Phys. Oceanogr.*, **42**, 1577–1600, https://doi.org/10.1175/JPO-D-11-0151.1.

* —, M. Müller, J. G. Richman, J. F. Shriver, A. J. Morten, R. B. Scott, G. Sérazin, and T. Penduff, 2014: Geostrophic turbulence in the frequency–wavenumber domain: Eddy-driven low-frequency variability. *J. Phys. Oceanogr.*, **44**, 2050–2069, https://doi.org/10.1175/JPO-D-13-054.1.

* Barnier, B., and Coauthors, 2006: Impact of partial steps and momentum advection schemes in a global ocean circulation model at eddy-permitting resolution. *Ocean Dyn.*, **56**, 543–567, https://doi.org/10.1007/s10236-009-0180-y.

* Barthel, A., A. M. Hogg, S. Waterman, and S. Keating, 2017: Jet–topography interactions affect energy pathways to the

deep Southern Ocean. *J. Phys. Oceanogr.*, **47**, 1799–1816, https://doi.org/10.1175/JPO-D-16-0220.1.

* Batchelor, G. K., 1953: *The Theory of Homogenous Turbulence*. Cambridge University Press, 197 pp.

* Berloff, P., and I. Kamenkovich, 2013a: On spectral analysis of mesoscale eddies. Part I: Linear analysis. *J. Phys. Oceanogr.*, **43**, 2505–2527, https://doi.org/10.1175/JPO-D-12-0232.1.

* —, and —, 2013b: On spectral analysis of mesoscale eddies. Part II: Nonlinear analysis. *J. Phys. Oceanogr.*, **43**, 2528–2544, https://doi.org/10.1175/JPO-D-12-0233.1.

* —, and —, 2018: Dynamics of baroclinic multiple zonal jets. *Zonal Jets: Phenomenology, Genesis, Physics*, B. Galperin and P. Read, Eds., Cambridge University Press, in press.

* —, W. Dewar, S. Kravtsov, and J. McWilliams, 2007a: Ocean eddy dynamics in a coupled ocean–atmosphere model. *J. Phys. Oceanogr.*, **37**, 1103–1121, https://doi.org/10.1175/JPO3041.1.

* —, A. M. Hogg, and W. Dewar, 2007b: The turbulent oscillator: A mechanism of low-frequency variability of the wind-driven ocean gyres. *J. Phys. Oceanogr.*, **37**, 2363–2386, https://doi.org/10.1175/JPO3118.1.

* —, I. Kamenkovich, and J. Pedlosky, 2009: A model of multiple zonal jets in the oceans: Dynamical and kinematical analysis. *J. Phys. Oceanogr.*, **39**, 2711–2734, https://doi.org/10.1175/2009JPO4093.1.

* —, S. Karabasov, J. T. Farrar, and I. Kamenkovich, 2011: On latency of multiple zonal jets in the oceans. *J. Fluid Mech.*, **686**, 534–567, https://doi.org/10.1017/jfm.2011.345.

* Bessières, L., and Coauthors, 2017: Development of a probabilistic ocean modelling system based on NEMO 3.5: Application at eddying resolution. *Geosci. Model Dev.*, **10**, 1091–1106, https://doi.org/10.5194/gmd-10-1091-2017.

* Blanke, B., and P. Delecluse, 1993: Variability of the tropical Atlantic Ocean simulated by a general circulation model with two different mixed-layer physics. *J. Phys. Oceanogr.*, **23**, 1363–1388, https://doi.org/10.1175/1520-0485(1993)023<1363:VOTTAO>2.0.CO;2.

* Bourouiba, L., D. N. Straub, and M. L. Waite, 2012: Non-local energy transfers in rotating turbulence at intermediate Rossby number. *J. Fluid Mech.*, **690**, 129–147, https://doi.org/10.1017/jfm.2011.387.

* Brodeau, L., B. Barnier, A.-M. Treguier, T. Penduff, and S. Gulev, 2010: An ERA40-based atmospheric forcing for global ocean circulation models. *Ocean Modell.*, **31**, 88–104, https://doi.org/10.1016/j.ocemod.2009.10.005.

* Bryan, F. O., 2013: Introduction: Ocean modeling—Eddy or not. *Ocean Modeling in an Eddying Regime*, Geophys. Monogr., Vol. 177, Amer. Geophys. Union, 1–3.

* Chapman, C. C., and R. Morrow, 2014: Variability of Southern Ocean jets near topography. *J. Phys. Oceanogr.*, **44**, 676–693, https://doi.org/10.1175/JPO-D-13-0182.1.

* Charney, J. G., 1971: Geostrophic turbulence. *J. Atmos. Sci.*, **28**, 1087–1095, https://doi.org/10.1175/1520-0469(1971)028<1087:GT>2.0.CO;2.

* Chavanne, C. P., P. Flament, D. S. Luther, and K.-W. Gurgel, 2010: Observations of vortex Rossby waves associated with a mesoscale cyclone. *J. Phys. Oceanogr.*, **40**, 2333–2340, https://doi.org/10.1175/2010JPO4495.1.

* Chelton, D. B., and M. G. Schlax, 1996: Global observations of oceanic Rossby waves. *Science*, **272**, 234–238, https://doi.org/10.1126/science.272.5259.234.

* —, R. A. de Szoeke, M. G. Schlax, K. El Naggar, and N. Siwertz, 1998: Geographical variability of the first baroclinic Rossby radius of deformation. *J. Phys. Oceanogr.*, **28**,
433–460, https://doi.org/10.1175/1520-0485(1998)028<0433: GVOTFB>2.0.CO;2.

—, M. G. Schlax, R. M. Samelson, and R. A. de Szoeke, 2007: Global observations of large oceanic eddies. *Geophys. Res. Lett.*, **34**, L15606, https://doi.org/10.1029/2007GL030812.

—, —, and —, 2011: Global observations of nonlinear mesoscale eddies. *Prog. Oceanogr.*, **91**, 167–216, https://doi.org/10.1016/j.poccean.2011.01.002.

Chen, R., 2013: Energy pathways and structures of oceanic eddies from the ECCO2 state estimate and simplified models. Ph.D. thesis, Massachusetts Institute of Technology, 206 pp., http://dspace.mit.edu/handle/1721.1/79154.

—, and G. R. Flierl, 2015: The contribution of striations to the eddy energy budget and mixing: Diagnostic frameworks and results in a quasigeostrophic barotropic system with mean flow. *J. Phys. Oceanogr.*, **45**, 2095–2113, https://doi.org/10.1175/JPO-D-14-0199.1.

—, —, and C. Wunsch, 2014: A description of local and nonlocal eddy–mean flow interaction in a global eddy-permitting state estimate. *J. Phys. Oceanogr.*, **44**, 2336–2352, https://doi.org/10.1175/JPO-D-14-0009.1.

Colin de Verdière, A., and T. Huck, 1999: Baroclinic instability: An oceanic wavemaker for interdecadal variability. *J. Phys. Oceanogr.*, **29**, 893–910, https://doi.org/10.1175/1520-0485 (1999)029<0893:BIAOWF>2.0.CO;2.

Dewar, W. K., 2003: Nonlinear midlatitude ocean adjustment. *J. Phys. Oceanogr.*, **33**, 1057–1082, https://doi.org/10.1175/1520-0485(2003)033<1057:NMOA>2.0.CO;2.

Dijkstra, H. A., and M. Ghil, 2005: Low-frequency variability of the large-scale ocean circulation: A dynamical systems approach. *Rev. Geophys.*, **43**, RG3002, https://doi.org/10.1029/2002RG000122.

Ducet, N., P. Y. Le Traon, and G. Reverdin, 2000: Global high-resolution mapping of ocean circulation from TOPEX/Poseidon and ERS-1 and -2. *J. Geophys. Res.*, **105**, 19477–19498, https://doi.org/10.1029/2000JC900063.

Dussin, R., B. Barnier, L. Brodeau, and J. M. Molines, 2016: The making of Drakkar forcing set DFS5. Drakkar–MyOcean Rep., 34 pp., https://www.drakkar-ocean.eu/forcing-the-ocean.

Fjørtoft, R., 1953: On the changes in the spectral distribution of kinetic energy for twodimensional, nondivergent flow. *Tellus*, **5**, 225–230, https://doi.org/10.3402/tellusa.v5i3.8647.

Forget, G., 2010: Mapping ocean observations in a dynamical framework: A 2004–06 ocean atlas. *J. Phys. Oceanogr.*, **40**, 1201–1221, https://doi.org/10.1175/2009JPO4043.1.

Fu, L.-L., and G. R. Flierl, 1980: Nonlinear energy and entrophy transfers in a realistically stratified ocean. *Dyn. Atmos. Oceans*, **4**, 219–246, https://doi.org/10.1016/0377-0265(80)90029-9.

Galperin, B., S. Sukoriansky, and N. Dikovskaya, 2010: Geophysical flows with anisotropic turbulence and dispersive waves: Flows with a β-effect. *Ocean Dyn.*, **60**, 427–441, https://doi.org/10.1007/s10236-010-0278-2.

Godeferd, F. S., and F. Moisy, 2015: Structure and dynamics of rotating turbulence: A review of recent experimental and numerical results. *Appl. Mech. Rev.*, **67**, 030802, https://doi.org/10.1115/1.4029006.

Grégorio, S., T. Penduff, G. Sérazin, J.-M. Molines, B. Barnier, and J. Hirschi, 2015: Intrinsic variability of the Atlantic meridional overturning circulation at interannual-to-multidecadal time scales. *J. Phys. Oceanogr.*, **45**, 1929–1946, https://doi.org/10.1175/JPO-D-14-0163.1.

Harper, K. L., S. V. Nazarenko, S. B. Medvedev, and C. Connaughton, 2014: Wave turbulence in the two-layer

ocean model. *J. Fluid Mech.*, **756**, 309–327, https://doi.org/10.1017/jfm.2014.465.

Hayashi, Y., 1980: Estimation of nonlinear energy transfer spectra by the cross-spectral method. *J. Atmos. Sci.*, **37**, 299–307, https://doi.org/10.1175/1520-0469(1980)037<0299: EONETS>2.0.CO;2.

—, 1982: Space-time spectral analysis and its applications to atmospheric waves. *J. Meteor. Soc. Japan*, **60**, 156–171, https://doi.org/10.2151/jmsj1965.60.1_156.

Hogg, A. M., and J. R. Blundell, 2006: Interdecadal variability of the Southern Ocean. *J. Phys. Oceanogr.*, **36**, 1626–1645, https://doi.org/10.1175/JPO2934.1.

Holland, W. R., 1978: The role of mesoscale eddies in the general circulation of the ocean—Numerical experiments using a wind-driven quasi-geostrophic model. *J. Phys. Oceanogr.*, **8**, 363–392, https://doi.org/10.1175/1520-0485(1978)008<0363: TROMEI>2.0.CO;2.

—, and D. B. Haidvogel, 1981: On the vacillation of an unstable baroclinic wave field in an eddy-resolving model of the oceanic general circulation. *J. Phys. Oceanogr.*, **11**, 557–568, https://doi.org/10.1175/1520-0485(1981)011<0557:OTVOAU>2.0.CO;2.

Hua, B. L., and D. B. Haidvogel, 1986: Numerical simulations of the vertical structure of quasi-geostrophic turbulence. *J. Atmos. Sci.*, **43**, 2923–2936, https://doi.org/10.1175/1520-0469(1986)043<2923:NSOTVS>2.0.CO;2.

Huck, T., O. Arzel, and F. Sévellec, 2015: Multidecadal variability of the overturning circulation in presence of eddy turbulence. *J. Phys. Oceanogr.*, **45**, 157–173, https://doi.org/10.1175/JPO-D-14-0114.1.

Hughes, C. W., 1996: The Antarctic Circumpolar Current as a waveguide for Rossby waves. *J. Phys. Oceanogr.*, **26**, 1375–1387, https://doi.org/10.1175/1520-0485(1996)026<1375: TACCAA>2.0.CO;2.

Klocker, A., and D. P. Marshall, 2014: Advection of baroclinic eddies by depth mean flow. *Geophys. Res. Lett.*, **41**, 3517–3521, https://doi.org/10.1002/2014GL060001.

—, —, S. R. Keating, and P. L. Read, 2016: A regime diagram for ocean geostrophic turbulence. *Quart. J. Roy. Meteor. Soc.*, **142**, 2411–2417, https://doi.org/10.1002/qj.2833.

Kolmogorov, A. N., 1941: Dissipation of energy in locally isotropic turbulence. *Dokl. Akad. Nauk SSSR*, **32**, 16–18.

Kraichnan, R. H., 1967: Inertial ranges in two-dimensional turbulence. *Phys. Fluids*, **10**, 1417, https://doi.org/10.1063/1.1762301.

—, 1971: Inertial-range transfer in two- and three-dimensional turbulence. *J. Fluid Mech.*, **47**, 525–535, https://doi.org/10.1017/S0022112071001216.

Lacascé, J. H., 2002: On turbulence and normal modes in a basin. *J. Mar. Res.*, **60**, 431–460, https://doi.org/10.1357/00224002762231160.

Le Sommer, J., T. Penduff, S. Theetten, G. Madec, and B. Barnier, 2009: How momentum advection schemes influence current-topography interactions at eddy permitting resolution. *Ocean Modell.*, **29**, 1–14, https://doi.org/10.1016/j.jocemod.2008.11.007.

Le Traon, P. Y., F. Nadal, and N. Ducet, 1998: An improved mapping method of multisatellite altimeter data. *J. Atmos. Oceanic Technol.*, **15**, 522–534, https://doi.org/10.1175/1520-0426(1998)015<0522:AIMMOM>2.0.CO;2.

Madec, G., 2008: NEMO ocean engine. IPSL Note du Pole de Modélisation 27, 217 pp.

Maximenko, N. A., B. Bang, and H. Sasaki, 2005: Observational evidence of alternating zonal jets in the world ocean. *Geophys. Res. Lett.*, **32**, L12607, https://doi.org/10.1029/2005GL022728.
McWilliams, J. C., L. P. Graves, and M. T. Montgomery, 2003: A formal theory for vortex Rossby waves and vortex evolution. *Geophys. Astrophys. Fluid Dyn.*, **97**, 275–309, https://doi.org/10.1080/0309192031000108698.

Molines, J.-M., B. Barnier, T. Penduff, A. M. Treguier, and J. Le Sommer, 2014: ORCA12.L46 climatological and interannual simulations forced with DFS4.4: GJM02 and MJM88. Drakkar Group Experiment Rep. GDRI-DRAKKAR-2014-03-19, 50 pp., https://www.drackar-ocean.eu/publications/reports/orca12_reference_experiments_2014.

Montgomery, M. T., and R. J. Kallenbach, 1997: A theory for vortex Rossby-waves and its application to spiral bands and intensity changes in hurricanes. *Quart. J. Roy. Meteor. Soc.*, **123**, 435–465, https://doi.org/10.1002/qj.49712353810.

Nazarenko, S., 2011: *Wave Turbulence*. Lecture Notes in Physics Series, Vol. 825, Springer, 279 pp.

Penduff, T., J. Le Sommer, B. Barnier, A.-M. Treguier, J. Molines, and G. Madec, 2007: Influence of numerical schemes on current-topography interactions in 1/4° global ocean simulations. *Ocean Sci.*, **3**, 509–524, https://doi.org/10.5194/os-3-509-2007.

—, M. Juza, L. Brodeau, G. C. Smith, B. Barnier, J.-M. Molines, A.-M. Treguier, and G. Madec, 2010: Impact of global ocean model resolution on sea-level variability with emphasis on interannual time scales. *Ocean Sci.*, **6**, 269–284, https://doi.org/10.5194/os-6-269-2010.

—, —, B. Barnier, J. Zika, W. K. Dewar, A.-M. Treguier, J.-M. Molines, and N. Audiffren, 2011: Sea level expression of intrinsic and forced ocean variabilities at interannual time scales. *J. Climate*, **24**, 5652–5670, https://doi.org/10.1175/JCLI-D-11-00077.1.

—, and Coauthors, 2014: Ensembles of eddying ocean simulations for climate. *CLIVAR Exchanges*, No. 19, International CLIVAR Project Office, Southampton, United Kingdom, 26–29.

Pierini, S., 2006: A Kuroshio extension system model study: Decadal chaotic self-sustained oscillations. *J. Phys. Oceanogr.*, **36**, 1605–1625, https://doi.org/10.1175/JPO2931.1.

—, 2011: Low-frequency variability, coherence resonance, and phase selection in a low-order model of the wind-driven ocean circulation. *J. Phys. Oceanogr.*, **41**, 1585–1604, https://doi.org/10.1175/JPO-D-10-05018.1.

Qiu, B., R. B. Scott, and S. Chen, 2008: Length scales of eddy generation and nonlinear evolution of the seasonally modulated South Pacific Subtropical Countercurrent. *J. Phys. Oceanogr.*, **38**, 1515–1528, https://doi.org/10.1175/2007JPO3856.1.

Rhines, P. B., 1975: Waves and turbulence on a beta-plane. *J. Fluid Mech.*, **69**, 417–443, https://doi.org/10.1017/S0022112075001504.

Salmon, R., 1998: *Lectures on Geophysical Fluid Dynamics*. Oxford University Press, 400 pp.

Scott, R. B., and F. Wang, 2005: Direct evidence of an oceanic inverse kinetic energy cascade from satellite altimetry. *J. Phys. Oceanogr.*, **35**, 1650–1666, https://doi.org/10.1175/JPO2771.1.

—, and B. K. Arbic, 2007: Spectral energy fluxes in geostrophic turbulence: Implications for ocean energetics. *J. Phys. Oceanogr.*, **37**, 673–688, https://doi.org/10.1175/JPO3027.1.

Sérazin, G., T. Penduff, S. Grégorio, B. Barnier, J.-M. Molines, and L. Terray, 2015: Intrinsic variability of sea level from global ocean simulations: Spatiotemporal scales. *J. Climate*, **28**, 4279–4292, https://doi.org/10.1175/JCLI-D-14-00554.1.

Sheng, J., and Y. Hayashi, 1990a: Estimation of atmospheric energetics in the frequency domain during the FGGE year. *J. Atmos. Sci.*, **47**, 1255–1268, https://doi.org/10.1175/1520-0469(1990)047<1255:EOAEIT>2.0.CO;2.

—, and —, 1990b: Observed and simulated energy cycles in the frequency domain. *J. Atmos. Sci.*, **47**, 1243–1254, https://doi.org/10.1175/1520-0469(1990)047<1243:OASECI>2.0.CO;2.

Simonnet, E., and H. A. Dijkstra, 2002: Spontaneous generation of low-frequency modes of variability in the wind-driven ocean circulation. *J. Phys. Oceanogr.*, **32**, 1747–1762, https://doi.org/10.1175/1520-0485(2002)032<1747:SGOLFM>2.0.CO;2.

Straub, D. N., and B. Nadiga, 2014: Energy fluxes in the quasi-geostrophic double gyre problem. *J. Phys. Oceanogr.*, **44**, 1505–1522, https://doi.org/10.1175/JPO-D-13-0216.1.

Sukoriansky, S., N. Dikovskaya, and B. Galperin, 2007: On the arrest of inverse energy cascade and the Rhines scale. *J. Atmos. Sci.*, **64**, 3312–3327, https://doi.org/10.1175/JAS4013.1.

Tulloch, R., J. Marshall, C. Hill, and K. S. Smith, 2011: Scales, growth rates, and spectral fluxes of baroclinic instability in the ocean. *J. Phys. Oceanogr.*, **41**, 1057–1076, https://doi.org/10.1175/2011JPO4404.1.

Vallis, G. K., 2006: *Atmospheric and Oceanic Fluid Dynamics*. Cambridge University Press, 964 pp.

—, and M. E. Maltrud, 1993: Generation of mean flows and jets on a beta plane and over topography. *J. Phys. Oceanogr.*, **23**, 1346–1362, https://doi.org/10.1175/1520-0485(1993)023<1346:GOMFAJ>2.0.CO;2.

Venaille, A., J. Le Sommer, J.-M. Molines, and B. Barnier, 2011: Stochastic variability of oceanic flows above topography anomalies. *Geophys. Res. Lett.*, **38**, L16611, https://doi.org/10.1029/2011GL048401.

Wunsch, C., 1997: The vertical partition of oceanic horizontal kinetic energy. *J. Phys. Oceanogr.*, **27**, 1770–1794, https://doi.org/10.1175/1520-0485(1997)027<1770:TVPOOH>2.0.CO;2.