# Imprint of Chaos on the Ocean Energy Cycle from an Eddying North Atlantic Ensemble

Takaya Uchida, Quentin Jamet, William K Dewar, Bruno Deremble, Andrew

Poje, Luolin Sun

► To cite this version:

Takaya Uchida, Quentin Jamet, William K Dewar, Bruno Deremble, Andrew Poje, et al. Imprint of Chaos on the Ocean Energy Cycle from an Eddying North Atlantic Ensemble. Journal of Physical Oceanography, 2024, 54 (3), pp.679-696. 10.1175/JPO-D-23-0176.1. hal-04813948

HAL Id: hal-04813948

https://hal.inrae.fr/hal-04813948v1

Submitted on 4 Dec 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Imprint of Chaos on the Ocean Energy Cycle from an
Eddying North Atlantic Ensemble

TAKAYA UCHIDA, a,b QUENTIN JAMET, b,c WILLIAM K. DEWAR, b,d BRUNO DEREMBLE, b ANDREW C. POJE, e
AND LUOLIN SUN b

(Manuscript received 29 August 2023, in final form 21 December 2023, accepted 29 December 2023)

ABSTRACT: We examine the ocean energy cycle where the eddies are defined about the ensemble mean of a partially air-sea coupled, eddy-rich ensemble simulation of the North Atlantic. The decomposition about the ensemble mean leads to a parameter-free definition of eddies, which is interpreted as the expression of oceanic chaos. Using the ensemble framework, we define the reservoirs of mean and eddy kinetic energy (MKE and EKE, respectively) and mean total dynamic enthalpy (MTDE). We opt for the usage of dynamic enthalpy (DE) as a proxy for potential energy due to its dynamically consistent relation to hydrostatic pressure in Boussinesq fluids and nonreliance on any reference stratification. The curious result that emerges is that the potential energy reservoir cannot be decomposed into its mean and eddy components, and the eddy flux of DE can be absorbed into the EKE budget as pressure work. We find from the energy cycle that while baroclinic instability, associated with a positive vertical eddy buoyancy flux, tends to peak around February, EKE takes its maximum around September in the wind-driven gyre. Interestingly, the energy input from MKE to EKE, a process sometimes associated with barotropic processes, becomes larger than the vertical eddy buoyancy flux during the summer and autumn. Our results question the common notion that the inverse energy cascade of wintertime EKE energized by baroclinic instability within the mixed layer is solely responsible for the summer-to-autumn peak in EKE and suggest that both the eddy transport of DE and transfer of energy from MKE to EKE contribute to the seasonal EKE maxima.

SIGNIFICANCE STATEMENT: The Earth system, including the ocean, is chaotic. Namely, the state to be realized is highly sensitive to minute perturbations, a phenomenon commonly known as the “butterfly effect.” Here, we run a sweep of ocean simulations that allow us to disentangle the oceanic expression of chaos from the oceanic response to the atmosphere. We investigate the energy pathways between the two in a physically consistent manner in the North Atlantic region. Our approach can be extended to robustly examine the temporal change of oceanic energy and heat distribution under a warming climate.

KEYWORDS: Eddies; Energy transport; Mesoscale processes; Ocean dynamics; Ensembles; Oceanic variability

**1. Introduction**

There has been much interest in the recent decades on
Earth's climate sensitivity, the long-term thermal response to
an increase in atmospheric carbon dioxide (Sherwood et al.
2020), motivated by the fact that our emission of anthropo-
genic carbon since the industrial revolution may be the culprit

for a warming climate (Masson-Delmotte et al. 2021). None-
theless, significant uncertainties persist in climate sensitivity
and have not been reduced since the Charney report pub-
lished in 1979 (Charney et al. 1979; Knutti et al. 2017). In un-
derstanding and quantifying the climate system, a useful
framework has been to examine the energy pathways, which
elucidates how much of the inc-oming solar radiation gets re-
tained and redistributed around the Earth system to warm or
cool the climate (Hartmann et al. 1986).

Among the components of the Earth system, the ocean is
perhaps the most significant reservoir of energy on centennial
to millennial time scales due to its large heat capacity and
density, and its ability to dissolve and store carbon and salts

Corresponding author: Takaya Uchida, tuchida@fsu.edu

<footnote>a Center for Ocean-Atmospheric Prediction Studies, Florida State University, Tallahassee, Florida</footnote>

<footnote>b Université Grenoble Alpes, CNRS, INRAE, IRD, Grenoble-INP, Institut des Géosciences de l'Environnement, Grenoble, France</footnote>

<footnote>c INRIA, ODYSSEY Group, Ifremer, Plouzané, France</footnote>

<footnote>d Department of Earth, Ocean and Atmospheric Science, Florida State University, Tallahassee, Florida</footnote>

<footnote>e Department of Mathematics, College of Staten Island, The City University of New York, New York, New York</footnote>

<footnote>Denotes content that is immediately available upon publica-
tion as open access.</footnote>
gives it a primary role in the global carbon cycle. In a seminal work, Wunsch and Ferrari (2004) attempted to provide an overview of the energy pathways for the oceans but came short in one crucial aspect: The role of mesoscale eddies as a conduit between the wind-driven general circulation and small-scale three-dimensional isotropic turbulence. Our lack of understanding on how mesoscale eddies interact with dynamics associated with other scales hinders our ability to understand the past and predict the future climate due to their disproportionately large role in globally transporting heat and carbon (Griffies et al. 2015; Gnanadesikan et al. 2015).

One of the difficulties in quantifying the impact of meso-scale eddies on global energetics lies in the identification of the eddies themselves (Wunsch 1981), which exist in a soup of anisotropic and inhomogeneous flows in space and are also nonstationary over an exceedingly large range of time scales (Uchida et al. 2022c, 2023b). Despite the lack of any homogeneous direction to average over, by necessity and practicality, eddies have often been defined via a Reynold's decomposition about a spatial and/or temporal coarse graining (e.g., Bachman et al. 2015; Aiki et al. 2016; Aoki et al. 2016; Uchida et al. 2017; Buzzicotti et al. 2023; Xie et al. 2023); this choice explicitly reduces the dimensionality of the mean field, and raises questions about the sensitivity of the results to the particularities of the averaging procedure. Aiki and Richards (2008) documented that by adjusting the temporal window over which the mean was taken, the amount of kinetic and potential energy stored in the mean and eddy reservoirs could change by up to a factor of 4. More recently, Demyshev and Dymova (2022) showed complementary results that depending on the time frame over which the averaging is taken to define the mean flow, the relative significance of energy pathways to the eddy kinetic energy (EKE) reservoir changed. As one may imagine, the amount of energy stored in each reservoir and exchanged among them is also dependent on the spatial scale taken for the decomposition (Loose et al. 2023).

Here, we take a different approach to the eddy–mean flow decomposition problem by running an ensemble of “eddy-rich” simulations of the North Atlantic Ocean using the Massachusetts Institute of Technology general circulation model (MITgcm; Marshall et al. 1997). By producing an ensemble of plausible ocean states at any given time, the flow can be averaged in the ensemble dimension and fluctuations are then defined about this ensemble mean. The ensemble mean field is itself fully spatially and temporally dependent and the only averaging parameter is the ensemble size; this is a convergent process where the “true” ensemble-mean state would be extracted upon having a sufficiently large ensemble. The mean fields defined via this decomposition can be interpreted as the oceanic response to the common atmospheric state and the eddy fields as an expression of the intrinsic variability, or chaotic component of the ocean (e.g., Chen and Flierl 2015; Sérazin et al. 2017; Leroux et al. 2018; Uchida et al. 2021a). The single-model ensemble approach has some history in the atmospheric and climate literature focusing on the dynamics and process-oriented studies (e.g., Sui et al. 1994; Lenderink et al. 2007; Nikiéma and Laprise 2013; Hersbach et al. 2015;

Xie et al. 2016; Romanou et al. 2023) but is still relatively novel in the field of oceanography (Aoki et al. 2020; Uchida et al. 2022b; Jamet et al. 2022; Zhao et al. 2023).

With this definition of mean flow and eddies, we will diagnose the ocean energy cycle (Lorenz 1955; Bleck 1985). Ensemble averaging is the only averaging operator that strictly commutes with all space–time derivatives and, therefore, allows us to close the mean–eddy energy cycle at any given time and across all available spatial scales. The goal here is to diagnose, throughout the course of a single model year, the components of the ocean energy cycle with an emphasis on the interaction between the kinetic and potential energy reservoirs. In doing so, we adopt a definition of potential energy in terms of dynamic enthalpy (DE). Formulating the problem in terms of DE has several advantages: DE provides a measure of potential energy that is both dynamically and thermodynamically consistent with the equations being solved (Eden 2015; Jamet et al. 2021), and DE does not depend on a reference state of stratification, providing a level of objectivity in defining the potential energy reservoir.

As we shall see, the energy cycle that emerges differs from the canonical Lorenz energy cycle primarily in that the energy reservoir corresponding to eddy available potential energy (APE) does not explicitly appear and instead the mean total potential energy reservoir directly interacts with both the mean and eddy kinetic energy (MKE and EKE, respectively) reservoirs. Our framework will elucidate that not only does the energy pathway from potential- to kinetic-energy reservoir consist of a local vertical buoyancy flux but also a nonlocal transport of DE, a mechanism overlooked by previous studies focusing on the energetics of North Atlantic wind-driven gyre (e.g., Kang and Curchitser 2015; Kang et al. 2016). We also find that barotropic processes contribute significantly to the summer-to-autumn EKE seasonal maxima, which have often solely been attributed to wintertime mixed-layer instability cascading upscale over time (Uchida et al. 2017; Khatri et al. 2021; Steinberg et al. 2022).

The paper is organized as follows: In section 2, we provide a brief description of the simulation and an overview on the energy cycle. Results are given in section 3 and we provide a summary in section 4.

## 2. Methods

### a. Model description

We use model outputs from a recently developed 48-member eddy-rich (1/12°) ensemble of the North Atlantic (Jamet et al. 2019a,b) partially air–sea coupled via the Cheap Atmospheric Mixed Layer model (CheapAML; Deremble et al. 2013). As our model domain is focused on the North Atlantic, the basin was configured to wrap around zonally in order to save memory allocation (e.g., Fig. 3). The dataset has been used to quantify the effect of oceanic chaos on the Atlantic meridional overturning circulation (Jamet et al. 2019b, 2020b; Dewar et al. 2022), and spatial variability of eddies, here defined about the ensemble mean and interpreted as the expression of oceanic
chaos (Uchida et al. 2022c, 2022b, 2023a,b). In this study, we shift our attention to the temporal variability by examining the energy cycle and shall refer to chaos as the physically consistent and deterministic yet inherent sensitivity of the system (here, taken as the North Atlantic) to its initial conditions (cf. Verhulst 1845; Poincaré 1890; Lorenz 1963). We will be analyzing the fifth year from ensemble initialization when the ensemble statistics have converged. The ensemble mean, being orthogonal to the spatiotemporal dimensions, commutes with the space-and-time derivatives and maintains the desirable statistical properties of nonstationarity and inhomogeneity upon a Reynolds decomposition. We use data from the year 1967 where ensemble outputs from the MITgcm diagnostics package are saved as instantaneous snapshots every five days, which allows us to close the total (mean + eddy) momentum budgets to machine precision. In other words, our analysis is somewhat restricted by the available model outputs in closing the budget. The kinetic energy (KE) budgets are subsequently constructed by taking the dot product between the horizontal momentum vector and each term in the momentum equations (see the appendix).

## b. Ocean energy cycle

The mean total kinetic energy (MTKE; $\langle |\mathbf{u}|^2 \rangle / 2$) can be decomposed into its mean and eddy kinetic energy (MKE and EKE) reservoirs as

$$K^{\#} \stackrel{\text{def}}{=} \langle |\mathbf{u}|^2 \rangle / 2, \quad (1)$$

$$\langle \mathcal{H} \rangle \stackrel{\text{def}}{=} \langle |\mathbf{u}'|^2 \rangle / 2, \quad (2)$$

where $\mathbf{u} = u\hat{\mathbf{x}} + v\hat{\mathbf{y}}$ is the horizontal momentum vector, $\hat{\mathbf{x}}$ and $\hat{\mathbf{y}}$ are the zonal and meridional unit vectors, respectively, $\langle \cdot \rangle$ is the ensemble mean operator, and $\langle \cdot \rangle^{\text{def}} = \langle \cdot \rangle - \langle \cdot \rangle$, $\langle \cdot \rangle^{\prime} = 0$. Regarding potential energy, while many possible ways to define it in oceanic primitive equations have been proposed (e.g., Oort et al. 1989; Winters et al. 1995; Aiki and Richards 2008; Molemaker and McWilliams 2010; Nycander 2010; von Storch et al. 2012; Saenz et al. 2015; Tailleux 2013, 2016; Dewar et al. 2016; Aiki et al. 2016; Kang et al. 2016; Guo et al. 2022; Yang et al. 2022; Loose et al. 2023; Demyshев and Dymova 2022; Steinberg and Eriksen 2022; Tailleux and Wolf 2023; Yang et al. 2023), defining a “reference” stratification in realistic simulations has remained subjective primarily due to the nonlinear equation of state (EOS) for seawater. The prescription of such reference state, furthermore, hinders the dynamical and thermodynamical consistency with the equations of motion being solved for a Boussinesq seawater. Namely, buoyancy (or density anomaly in defining buoyancy) must satisfy the hydrostatic pressure relation while remaining a thermodynamical function when considering the energetics. We opt for dynamic enthalpy (DE; Young 2010), which is a (if not the only) natural extension of gravitational potential energy with a nonlinear EOS, and does not depend on a reference state:

$$\tilde{h}(\Theta, S, \Phi) = \int_{\Phi_0}^{\Phi} \frac{\tilde{b}(\Theta, S, \Phi^*)}{g} d\Phi^* = \int_{z}^{0} bdz^*, \quad (3)$$

where $\Phi = \Phi_0 - gz$ is the static dynamically nonactive part of hydrostatic pressure, the superscript * indicates a dummy variable, and the tilde ($\tilde{\cdot}$) denotes a thermodynamic function. Following Young (2010), the tilde notation distinguishes thermodynamic functions from fields in space-time; that is, $\tilde{b}(\Theta, S, \Phi) = b(t, z, y, x)$. The terms $\Theta$ and $S$ are potential temperature and practical salinity, $\tilde{b} = -g[(\tilde{\rho} - \rho_0)/\rho_0]$ is buoyancy, $\tilde{\rho}$ is density based on Jackett and McDougall (1995), $\rho_0 = 999.8 \text{ kg m}^{-3}$ is the reference density prescribed in MITgcm, and $g$ is gravity. Although buoyancy in Young (2010) was defined as $\tilde{b} = -g[(\tilde{\rho} - \rho_0)/\tilde{\rho}]$, we make use of the former convention for simplicity (particularly when taking its partial derivatives) and will neglect the small differences that emerge between the two (cf. Eden 2015). We emphasize that the integration in (3) is taken by fixing $\Theta$ and $S$ in respect to $\Phi(z)$, e.g., $\Theta = \Theta(t, z, y, x)$. The term on the right-hand side of (3) in parentheses shows the integration by substituting $\Phi$ with $z$ in space-time.

## 1) FULLY NONLINEAR THERMODYNAMICS

As detailed in the appendix, the evolution equations for MKE and EKE are

$$\frac{D^{\#}}{Dt} K^{\#} = -\langle \mathbf{v}' \cdot \nabla (\langle \mathbf{u} \rangle \cdot \mathbf{u}') \rangle - \langle \mathbf{v} \rangle \cdot \nabla \langle \phi \rangle + \langle w \rangle \langle b \rangle \\ + \langle \mathbf{u}' \mathbf{v}' \rangle \cdot \nabla \langle \mathbf{u} \rangle + \langle \mathbf{u} \rangle \cdot \langle \mathcal{X} \rangle, \quad (4)$$

$$\frac{D^{\#}}{Dt} \langle \mathcal{H} \rangle = -\langle \mathbf{v}' \cdot \nabla \mathcal{H} \rangle - \langle \mathbf{v}' \cdot \nabla \phi' \rangle + \langle w' b' \rangle \\ - \langle \mathbf{u}' \mathbf{v}' \rangle \cdot \nabla \langle \mathbf{u} \rangle + \langle \mathbf{u}' \cdot \mathcal{X}' \rangle, \quad (5)$$

respectively where $D^{\#}/Dt \stackrel{\text{def}}{=} (\partial/\partial t) + \langle \mathbf{v} \rangle \cdot \nabla$ is the mean Lagrangian tendency operator, $\mathbf{v} = \mathbf{u} + w\hat{\mathbf{z}}$ the nondivergent three-dimensional momentum vector, $\phi$ the dynamically active part of hydrostatic pressure, and $\mathcal{X} (= \mathcal{F} + \boldsymbol{\epsilon})$ the net nonconservative term consisting of forcing, viscous dissipation, and contribution from the K-profile parameterization (KPP; Large et al. 1994) to the momentum equations.

On the other hand, ensemble averaging the Lagrangian tendency of total dynamic enthalpy [TDE; (A7)] under adiabatic conditions yields

$$\frac{D^{\#}}{Dt} \langle h \rangle + \langle \mathbf{v}' \cdot \nabla h' \rangle = \left\langle \frac{D}{Dt} h \right\rangle = -\langle w \rangle \langle b \rangle - \langle w' b' \rangle. \quad (6)$$

Since DE is a thermodynamic function, the mean Lagrangian tendency of mean total dynamic enthalpy (MTDE; $\langle h \rangle$) can also be expressed as

$$\frac{D^{\#}}{Dt} \langle \tilde{h} \rangle = \left\langle \frac{D^{\#}}{Dt} \tilde{h} \right\rangle = \langle \tilde{h}_{\Phi} \rangle \frac{D^{\#}}{Dt} \Phi + \left\langle \tilde{h}_{\Theta} \frac{D^{\#}}{Dt} \Theta \right\rangle + \left\langle \tilde{h}_{S} \frac{D^{\#}}{Dt} S \right\rangle \\ = -\langle w \rangle \langle \tilde{b} \rangle + \mathcal{H}, \quad (7)$$

where $\mathcal{H} \stackrel{\text{def}}{=} \langle \tilde{h}_{\Theta}(D^{\#}/Dt)\Theta \rangle + \langle \tilde{h}_{S}(D^{\#}/Dt)S \rangle$ encapsulates the chain rule in respect to potential temperature and practical salinity, and is proportional to their eddy flux and diabatic, molecular, and nonhydrostatic effects. The latter three effects
are generally ignored hereon; the adiabatic approximation is made for the thermodynamics as the tendency terms for temperature and salinity were not saved as model outputs. The subscripts (·)$_{Φ}$, (·)$_{Θ}$, and (·)$_{S}$ denote partial derivatives in thermodynamics.

Subtracting (7) from (6) leaves us with the identity

$$ \langle \mathbf{v}' \cdot \nabla h' \rangle = -\langle w' b' \rangle - \mathcal{H}. \quad (8) $$

With this, the full budgets for EKE and MTDE can be written as

$$ \begin{aligned} \frac{D^{\#}}{Dt} \langle \mathcal{H} \rangle = & -\langle \mathbf{v}' \cdot \nabla \mathcal{H} \rangle - \langle \mathbf{v}' \cdot \nabla (\phi' + h') \rangle - \mathcal{H} \\ & - \langle \mathbf{u}' \mathbf{v}' \rangle \cdot \nabla \langle \mathbf{u} \rangle + \langle \mathbf{u}' \cdot \mathcal{X}' \rangle, \end{aligned} \quad (9) $$

$$ \frac{D^{\#}}{Dt} \langle h \rangle = -\langle w \rangle \langle b \rangle + \mathcal{H}. \quad (10) $$

The curious, and interesting, aspect of (9) and (10) is that the direct transfer of energy between kinetic and potential energy reservoirs is now encapsulated in the thermodynamically consistent $\mathcal{H}$ term whose sign dictates the direction of the flux. An additional term, the divergence of the eddy flux of DE, also contributes to the EKE tendency. Written as $-\langle \nabla \cdot [\mathbf{v}'(\phi' + h')] \rangle$, this can be interpreted as contributions from the fluctuating component of the dynamically active hydrostatic pressure; such interpretation is motivated by the fact that hydrostatic pressure and DE are both vertically integrated properties of buoyancy.

A subtle difference between MTKE and MTDE is that the former is quadratic while the latter is a single-order variable, yet both have the dimension of energy; MTDE cannot be explicitly decomposed into its mean and eddies like MTKE. In the quasigeostrophic sense, MTDE is the combined reservoir of mean and eddy available potential energy (APE). Nonetheless, (4), (9) and (10) form a complete set of equations to describe the energy cycle in primitive equations.

The problem arises, however, that $\langle \tilde{b} \rangle$ and $\langle \tilde{h} \rangle$ are no longer thermodynamic functions of $\langle \Theta \rangle$ and $\langle S \rangle$ for a nonlinear EOS because the ensemble-mean operator does not commute with the nonlinearity; consequently, neither are $b'$ and $h'$ thermodynamic functions of $\Theta'$ and $S'$. In other words, the energy cycle loses its direct ties with the thermodynamics as the potential energy reservoir cannot be thermodynamically expressed.

## 2) LINEARIZATION OF THE THERMODYNAMICS

We can make further progress by appealing to the approximation that second- and higher-order terms of thermodynamics are negligible, namely $\langle \tilde{b} \rangle = \tilde{b} + \tilde{B}$:

$$ \langle \tilde{h} \rangle = \int_{\Phi_0}^{\Phi} \frac{\tilde{b} + \tilde{B}}{g} d\Phi^* \approx \tilde{H}. \quad (11) $$

The term $\tilde{H} \stackrel{\text{def}}{=} g^{-1} \int_{\Phi_0}^{\Phi} \tilde{b} d\Phi^*$ where $\tilde{b} \stackrel{\text{def}}{=} \tilde{b}(\langle \Theta \rangle, \langle S \rangle, \Phi)$ is buoyancy given by the mean potential temperature and practical salinity, and the nonlinearity in EOS, which is second order in temperature and salinity fluctuations at most, is shouldered by $\tilde{B}$ and ignored [cf. (A9) and Fig. A2]. In view of (11) and realizing

that $(D^{\#}/Dt)\langle \Theta \rangle \approx -\langle \mathbf{v}' \cdot \nabla \Theta' \rangle$ and $(D^{\#}/Dt)\langle S \rangle \approx -\langle \mathbf{v}' \cdot \nabla S' \rangle$, (8) can be simplified as

$$ \langle w' b' \rangle \approx -\langle \mathbf{v}' \cdot \nabla h' \rangle + \tilde{H}_{\Theta} |_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{S} |_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle, \quad (12) $$

where $\tilde{H}_{\Theta}|_{\langle \Theta \rangle}$ and $\tilde{H}_{S}|_{\langle S \rangle}$ are partial derivatives with respect to $\Theta$ and $S$ assessed at $\langle \Theta \rangle$ and $\langle S \rangle$. For the remainder of this article $\tilde{b}' \stackrel{\text{def}}{=} \tilde{b} - \tilde{b}' (\approx b')$ and $\tilde{h}' \stackrel{\text{def}}{=} \tilde{h} - \tilde{H} (\approx h')$ [cf. (A14) and Fig. A3; (8) and (12) are identical under a linearization of the thermodynamics].

Following through with the linear approximation and plugging (12) into (9) and (10) allows us to reformulate the evolution equations of MKE and EKE as

$$ \begin{aligned} K_t^{\#} + \langle \mathbf{v} \rangle \cdot \nabla K^{\#} \approx & -\langle \mathbf{v} \rangle \cdot \nabla \bar{\Phi} + \underbrace{\langle w \rangle \bar{b}}_{=H_{\to, \mathcal{K}^{\#}}} - \nabla \cdot \langle \mathbf{v}'(\langle \mathbf{u} \rangle \cdot \mathbf{u}') \rangle \\ & + \underbrace{\langle \mathbf{v}' \mathbf{u}' \rangle \cdot \nabla \langle \mathbf{u} \rangle}_{=K_{\to, \mathcal{H}^{\#}}} + \langle \mathbf{u} \rangle \cdot \langle \mathcal{X}' \rangle, \end{aligned} \quad (13) $$

$$ \begin{aligned} \langle \mathcal{H} \rangle_t + \langle \mathbf{v} \rangle \cdot \nabla \langle \mathcal{H} \rangle \approx & -\langle \mathbf{v}' \cdot \nabla (\mathcal{H} + \phi' + h') \rangle - \underbrace{\langle \mathbf{v}' \mathbf{u}' \rangle \cdot \nabla \langle \mathbf{u} \rangle}_{=K_{\to, \mathcal{H}^{\#}}} + \langle \mathbf{u}' \cdot \mathcal{X}' \rangle \\ & + \underbrace{\tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{S}|_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle}_{=H_{\to, \mathcal{H}^{\#}}}, \end{aligned} \quad (14) $$

and MTDE as

$$ H_t + \langle \mathbf{v} \rangle \cdot \nabla H \approx \underbrace{-\langle w \rangle \bar{b}}_{=H_{\to, \mathcal{K}^{\#}}} - \underbrace{\langle \tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{S}|_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle}_{=H_{\to, \mathcal{H}^{\#}}}, \quad (15) $$

respectively, where $\bar{b}(t, z, y, x) = \tilde{b} = \tilde{b}(\langle \Phi \rangle, \langle S \rangle, \Phi)$, $\bar{\Phi}_z \stackrel{\text{def}}{=} \bar{b}(\approx \langle \phi_z \rangle)$, and $\phi_z^{\#} \stackrel{\text{def}}{=} b^{\#}(\approx \phi_z^{\#})$, and the subscripts $(\cdot)_t$ and $(\cdot)_z$ denote partial derivatives in space-time. The energy cycle with full consideration of a nonlinear EOS and diabatic terms is derived in the appendix. A joint histogram demonstrates that (14) holds relatively well considering the simplifications we have made to the thermodynamics (Fig. 1). This implies that the eddy transport of temperature and salinity dominates over the explicit diffusive transport of the tracers.

With the energy reservoirs defined, we can identify the energy exchanges among the reservoirs. The exchange between MKE and EKE ($K_{\to, \mathcal{H}^{\#}}^{\#}$), MTDE and MKE ($H_{\to, \mathcal{K}^{\#}}$), and MTDE and EKE ($H_{\to, \mathcal{H}^{\#}}$) reservoirs, which are nonzero upon a global volume integration, are

$$ K_{\to, \mathcal{H}^{\#}}^{\#} \stackrel{\text{def}}{=} -\langle (\mathbf{u}' \mathbf{v}') \cdot \nabla \langle \mathbf{u} \rangle + (\mathbf{v}' \mathbf{v}') \cdot \nabla \langle \mathbf{v} \rangle \rangle, \quad (16) $$

$$ H_{\to, \mathcal{K}^{\#}} \stackrel{\text{def}}{=} \langle w \rangle \bar{b} \left( = \tilde{H}_{\Phi} \frac{D^{\#}}{Dt} \Phi \approx \langle w \rangle \langle b \rangle \right), \quad (17) $$

$$ H_{\to, \mathcal{H}^{\#}} \stackrel{\text{def}}{=} \tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{S}|_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle (\approx \langle w' b' \rangle + \langle \mathbf{v}' \cdot \nabla h' \rangle), \quad (18) $$

as indicated from (13) to (15). Unlike the canonical Lorenz energy cycle (Lorenz 1955; Uchida et al. 2021b), notice that there is no term corresponding to the exchange with the
<figure><img src="image_3.png" /><figcaption>FIG. 1. Joint histogram of $\langle \mathcal{K}_t \rangle$ and the sum of other terms in (14) including the advection term on the left-hand side where $\mathbf{u} \cdot \nabla_h \phi^\dagger$ is approximated by the MITgcm diagnostics output of $\mathbf{u} \cdot \nabla_h \phi'$. The histogram was computed for 1 Jan 1967 over the entire three-dimensional domain of the ensemble output. A one-to-one line is shown as the gray dashed line. The histogram was computed using the xhistogram Python package (Abernathey et al. 2021b). Spatial maps of the budget terms are given in Fig. A4.</figcaption></figure>

mean and eddy APE reservoirs but rather that MTDE directly interacts with the MKE and EKE reservoirs (Fig. 2). The term $\mathbf{v} \cdot \nabla h^\dagger$, which would seemingly identify as the exchange between mean and eddy APE reservoirs, gets directly passed onto EKE with the DE fluctuation serving as its conduit and retains no energy as $h^\dagger$ (i.e., $\langle h^\dagger \rangle \approx 0$). We emphasize that $h^\dagger$ is not equivalent to eddy APE; claiming otherwise would amount to saying that under quasi-geostrophy $[b^2/(2N^2)]' = [b^2/(2N^2)] - \langle b^2/(2N^2) \rangle$ is the eddy APE, whereas the correct eddy APE is actually $\langle b^\dagger \rangle^2/(2N^2)$, where $N^2$ is the background stratification. Furthermore, (14) demonstrates that the eddy flux divergence of DE can be consolidated as pressure work, which leaves us with the effect of eddy temperature and salinity flux divergence in (18). Figure A3 exhibits that (18) approximately holds throughout the domain. We argue that this deviation from the Lorenz cycle of a nonexplicit eddy APE reservoir results from the fact that quasi-geostrophy corresponds to the thickness-weighted primitive equations of motion in isopycnal coordinates, and not the unweighted equations in geopotential coordinates. Under quasi-geostrophy, the isopycnal layer thickness is constant, leading to quasigeostrophic (QG) variables being implicitly thickness-weighted averaged (Marshall et al. 2012; Maddison and Marshall 2013; Uchida et al. 2023a; Meunier et al. 2023). Nonetheless, an energy budget can be formulated for non-thickness-weighted primitive equations under geopotential coordinates [cf. (13)–(15) and the appendix; Eden 2015] so we proceed in examining the energy cycle under this formalism. Our notation regarding the variables is summarized in Table 1.

### 3. Results

We start by showing the winter and summertime MKE, EKE, and MTDE vertically averaged over the surface 1000 m, a depth over which the wind-driven gyre is roughly contained

<figure><img src="image_4.png" /><figcaption>FIG. 2. A schematic of the ensemble Reynold's decomposition of the energy cycle using dynamic enthalpy formulation (neglecting the external forcing and diabatic terms). The ensemble-mean total kinetic energy (MTKE) is easily split into kinetic energy of the mean (MKE) and of the eddies (EKE). The ensemble-mean total dynamic enthalpy (MTDE) cannot be split, but may be simplified if $\langle \tilde{b} \rangle \approx \tilde{b}(\langle \Theta \rangle, \langle S \rangle, \Phi)$. In this case the interaction terms, $H_{\leftrightarrow}K^\# \langle \tilde{h}(\langle \Theta \rangle + \Theta', \langle S \rangle + S', \Phi) \rangle \approx \tilde{H}$ (in blue font) and $H_{\leftrightarrow}K^\# \langle \tilde{h}(\langle \Theta \rangle + \Theta', \langle S \rangle + S', \Phi) \rangle \approx \tilde{H}$ (in red font), can be written explicitly, i.e., (17) and (18).</figcaption></figure>

(Jamet et al. 2021). The Gulf Stream, Gulf of Mexico Loop Current, Equatorial Undercurrent, and North Brazil Current are apparent in MKE (Figs. 3a,b) while EKE is more concentrated around the separated Gulf Stream and North Atlantic Current region (Figs. 3c,d). Our ensemble-based eddies typically have spatial scales of about 300 km (Uchida et al. 2022c, 2023b). One may notice that MTDE is negative, which has to do with buoyancy always taking negative values due to $\rho_0$ used in MITgcm. Conceptually, MTDE can be viewed as a well of potential energy. In the subsections below, we examine the time series of the volume-averaged energy cycle. Spatial maps of the terms in the MKE and EKE budgets are given in Figs. A1 and A4, respectively.

#### a. Entire model domain

We show in Fig. 4 the time series of MTDE, EKE, and MKE, along with the energy exchange between the reservoirs volume-averaged over the entire domain (20°S–55°N, 262°–348°E, excluding the most north and south grid points where the lateral boundary conditions were prescribed; Jamet et al. 2019b) and full water column. The magnitude of EKE is larger than MKE. EKE has two local maxima about March and November respectively while MKE seemingly lags one to two months behind EKE also with a dual peak. MTDE is orders of magnitude larger than MKE and EKE, which is consistent with our notion that most of the dynamical energy in the ocean is stored as potential energy.

Shifting our attention to Fig. 4b, the energy input from MTDE to EKE takes its maximum during March, which is similar to $\langle w' b' \rangle$, implying that baroclinic instability is active during boreal winter. The similarity implies that there is a negligible amount of net eddy influx of DE fluctuation ($\nabla \cdot \langle v' h' \rangle \le 0$) at the north and south open boundaries of our domain at 20°S and 55°N; this is by construction as the same
<table><caption>TABLE 1. Definition and description of the variables and notation.</caption><thead><tr><th>Mathematical notation</th><th>Description</th></tr></thead><tbody><tr><td>⟨K⟩<sub>def</sub> = ⟨|u|<sup>2</sup>⟩/2</td><td>Mean total kinetic energy (MTKE) reservoir</td></tr><tr><td>K<sup>#</sup><sub>def</sub> = ⟨|u⟩<sup>2</sup>⟩/2</td><td>Mean kinetic energy (MKE) reservoir</td></tr><tr><td>⟨⟨K⟩⟩<sub>def</sub> = ⟨⟨|u|⟩<sup>2</sup>⟩/2</td><td>Eddy kinetic energy (EKE) reservoir</td></tr><tr><td>˜b(Θ, S, Φ) = −g(ρ̃(Θ, S, Φ) − ρ<sub>0</sub>) / ρ<sub>0</sub> = b(t, z, y, x)</td><td>Tilde notation to distinguish a thermodynamic function (Jackett and McDougall 1995) from fields in space–time; ρ<sub>0</sub> = 999.8 kg m<sup>−3</sup></td></tr><tr><td>˜b<sup>def</sup> = ˜b(⟨Θ⟩, ⟨S⟩, Φ) ≈ ⟨˜b(Θ, S, Φ)⟩</td><td>Mean buoyancy</td></tr><tr><td>˜b<sup>†</sup> = ˜b − ˜b ≈ b'</td><td>Buoyancy perturbation</td></tr><tr><td>˜H = ∫<sub>0</sub><sup>δ</sup> g<sup>−1</sup> ˜b dΦ<sup>*</sup> ≈ ⟨˜h(Θ, S, Φ)⟩</td><td>Mean total dynamic enthalpy (MTDE) reservoir</td></tr><tr><td>˜h<sup>†</sup> = ˜h − ˜H ≈ h'</td><td>Dynamic enthalpy perturbation</td></tr><tr><td>K<sup>#</sup><sub>↔</sub>, K<sup>def</sup> = −⟨⟨u′v′⟩·∇⟨u⟩ + ⟨v′v′⟩·∇⟨v⟩⟩</td><td>MKE and EKE exchange</td></tr><tr><td>H<sub>↔</sub>, K<sup>#</sup> = ⟨w⟩b</td><td>MTDE and MKE exchange</td></tr><tr><td>H<sub>↔</sub>, K<sup>def</sup> = ˜H<sub>Θ|⟨Θ⟩</sub>⟨v′·∇Θ′⟩ + ˜H<sub>S|⟨S⟩</sub>⟨v′·∇S′⟩</td><td>MTDE and EKE exchange</td></tr></tbody></table>

lateral boundary conditions are applied across all members although surface DE conditions differ among the members (Jamet et al. 2019a,b). Energy input to EKE from MKE is also positive year around although with two local maxima around February and August–October respectively. Although noisy, energy from MTDE is fluxed to MKE for most of the year with largest values during March, and is two orders of magnitude larger than the energy fluxes to EKE (cf. Fig. A3). The EKE and MKE advective influx from the boundaries (i.e., −⟨v·∇·K⟩ and −[⟨v⟩·∇K<sup>#</sup> + ∇·⟨v′(⟨u⟩·⟨u′⟩)⟩], respectively) are negligible compared to the flux between the energy reservoirs and exhibit no seasonality (gray solid curves in Figs. 4b,f). Note that $H_{↔}K^#$ is roughly two orders of magnitude larger than the other terms shown in the right column of Fig. 4 but the mean vertical pressure work cancels it out (not shown), leaving the mean horizontal pressure work as the net contribution (−⟨u⟩·∇<sub>h</sub>Φ, where ∇<sub>h</sub> is the horizontal gradient operator; cf. Figs. A1b,c). Another interesting thing to note is the existence of two local maxima in both the EKE and $K_{↔}^#$, K<sup>def</sup> time series. The dual peak in the two time series seemingly implies that the energy input from MKE to EKE (sometimes associated with barotropic instability) dominates over baroclinic instability in modulating the temporal variability of EKE. The magnitude of $K_{↔}^#$, K<sup>def</sup> (green solid curve in Fig. 4b) being larger than $H_{↔}$ (red solid curve in Fig. 4b) later into the year also corroborates our argument.

Energy input from forcing, largely due to wind stress, to MKE is consistently positive (⟨u⟩·⟨F⟩ > 0) while it tends to damp out the eddies (⟨u'·F'⟩ < 0), although the latter is smaller by two orders of magnitude than the former (cf. Figs. A1d and A4d). The eddy killing effect by wind stress and thermal feedback is consistent with recent findings (Renault et al. 2016, 2023), albeit severely underestimated here due to us prescribing absolute wind stress instead of relative. We remind the reader that CheapAML follows the Coupled Ocean–Atmosphere Response Experiment (COARE3) prescription for wind forcing that

enters the momentum equations (Fairall et al. 2003) and hence ⟨u'·F'⟩ has a weak dependence on sea surface temperature fluctuations even when absolute wind stress is prescribed. Dissipation is consistently a sink for KE and takes the largest magnitude in winter, namely ⟨u⟩·⟨ε⟩. The EKE dissipation ⟨u'·ε'⟩ tends to mirror $H_{↔}$ and ⟨w'b'⟩. The fact that dissipation (here established by harmonic and biharmonic numerical viscosity) is a leading-order term for EKE may indicate a forward cascade of KE and/or wind-induced mixing via KPP particularly during boreal winter. It is possible that the results on dissipation are sensitive to our combination of free-slip lateral boundary conditions and no bottom drag but this does not undermine the fact that the KE budgets can be closed using DE as the potential energy reservoir.

b. Wind-driven gyre

We now focus on the subdomain of 270°–337°E, 14°–43°N and the upper 1000 m as the wind-driven gyre region, a domain somewhat similar to that in Jamet et al. (2021). The wind-driven gyre imprints itself onto MTDE as a shoaling of the potential energy well (Figs. 3e,f). The domain-averaged time series show that EKE roughly doubles MKE and takes its maximum around August/September while MKE exhibits a more pronounced dual peak. MTDE is in sync with EKE (Fig. 5a).

Unlike when averaged over the entire model domain, the energy flux from MTDE to EKE remains roughly twice as large as the eddy vertical buoyancy flux ($H_{↔}K^# > ⟨w'b'⟩$; red solid and magenta dashed curves in Fig. 5b); the eddy outflux of DE fluctuation is nonnegligible across the southern border at 14°N, namely ∇·⟨v'h'⟩ > 0 [cf. (12)], and remains relatively constant throughout the year. The relative magnitude between $H_{↔}$ and ⟨w'b'⟩ was insensitive to zonal and vertical expansion of the subdomain and extension northward (not shown). The energy flux from MTDE to MKE ($H_{↔}K^#$) tends to change sign, indicating an occasional steepening of isopycnals in the mean flow. It is again canceled out by the mean
<figure><img src="image_5.png" /><figcaption>FIG. 3. (a)–(d) MKE and EKE vertically averaged over the top 1000 m for winter (January–March) and summer (July–September) of 1967. (e), (f) MTDE, also averaged over the top 1000 m, for each season. The subdomain considered for the wind-driven gyre in section 3b is indicated with the cyan dotted lines.</figcaption></figure>

vertical pressure work in closing the MKE budget (not shown). Unlike the eddy flux of DE fluctuation, the influx of KE into the subdomain due to advection is negligible for both the eddy and mean flow (gray solid and blue dashed curves in Fig. 5b).

Interestingly, similar to Fig. 4, the energy flux from MKE to EKE becomes larger than the eddy vertical buoyancy flux over the summer and autumn (green solid and magenta

dashed curves in Fig. 5b), which implies that locally, barotropic processes are still the regulating mechanism over baroclinic. The KE maximum during boreal summer and autumn in western boundary current regions (e.g., the Kuroshio and Gulf Stream) has been observed in nature and is often explained as the time lag for the submesoscale eddies energized by wintertime baroclinic instability within the surface mixed layer to locally cascade upscale to the mesoscale (e.g., Zhai
<figure><img src="image_6.png" /><figcaption>FIG. 4. Time series of the volume-averaged energy stored in each reservoir and nondivergent terms in the KE budget: (a)–(d) EKE ($\mathcal{K}$) and (e)–(h) MKE ($K^#$). MTDE ($H$) is plotted against the left y-axis in orange in (a) with values varying about –1961 m² s⁻². The exchange between MTDE and EKE ($H \leftrightarrow \mathcal{K}$; solid red curve) and MKE and EKE ($K^# \leftrightarrow \mathcal{K}$; solid green curve) and eddy vertical buoyancy flux (magenta dashed curve) are plotted in (b) while MTDE and MKE ($H \leftrightarrow K^#$; blue solid curve) is plotted in (f). Colors representing the energy flux between reservoirs in (b) and (f) correspond to Fig. 2. The contribution from the advective terms (–$\langle \mathbf{v} \cdot \nabla \mathcal{K} \rangle$) and $–[\langle \mathbf{v} \rangle \cdot \nabla K^# + \nabla \cdot \langle \mathbf{v}'(\mathbf{u}' \cdot \mathbf{u}) \rangle]$) are in gray solid curves in (b) and (f). The diabatic terms, forcing and dissipation, are plotted in (c) and (g) and (d) and (h), respectively.</figcaption></figure>

et al. 2008; Sasaki et al. 2014; Uchida et al. 2017; Dong et al. 2020). Our results imply that seasonality in the Gulf Stream Extension is also modulated strongly by the eddy flux divergence of DE fluctuation from the region and a local KE transfer from the mean flow to eddies; the significance of the latter mechanism ($K^# \leftrightarrow \mathcal{K}$) is consistent with Uchida et al. (2021b), who diagnosed the energy cycle from a seasonally forced quasi-stigeostrophic double-gyre ensemble.

Regarding the diabatic terms, the relative contribution of viscous dissipation increases compared to when averaged over the full water column (Figs. 4d and 5d), which is attributable to KPP in the surface mixed layer. EKE dissipation tends to mirror $\langle w'b^\dagger \rangle$ with comparable magnitude (magenta dashed curve in Fig. 5b and black solid curve in Fig. 5d), implying that much of the conversion from potential to kinetic energy due to baroclinic instability is lost locally to dissipation. Dissipation, a driver for a forward cascade of EKE at our model resolution (Molemaker et al. 2010; Arbic et al. 2013), again peaks

during boreal winter, consistent with the seasonality found by Contreras et al. (2023).

We end this section by comparing the annual mean of the ensemble-based energy cycle to where the eddy–mean flow decomposition was done via a temporal filter from an arbitrary single realization, namely, about the annual mean of 1967, in order to be consistent with the atmospheric state seen by the ensemble. For the latter, the mean energy reservoirs are defined as $K^{\#t} \stackrel{\text{def}}{=} |\mathbf{u}'|^2/2$ and $\tilde{H}^t \stackrel{\text{def}}{=} \int_{\Phi_0}^{\Phi} g^{-1} \bar{b}(\bar{\Theta}', \bar{S}', \Phi') d\Phi^+$, where $\overline{(\cdot)^t}$ is the annual-mean operator, and subsequently the eddies are defined about the annual mean (e.g., $\Theta'^t \stackrel{\text{def}}{=} \Theta - \bar{\Theta}'$, $\mathcal{K}^t \stackrel{\text{def}}{=} |\mathbf{u}'|^2/2$). The magnitude of energy pathways in the ensemble and temporal framework tends to be of the same order (Fig. 6), which is comforting but also may seem to undercut the utility of ensembles. However, because the temporal framework takes the averaging over the time dimension, information on temporal variability within 1967 is lost (cf. Fig. 5). While investigation on longer time series is
<figure><img src="image_7.png" /><figcaption>FIG. 5. As in Fig. 4 but for the subdomain of wind-driven gyre. The subdomain is shown in Fig. 3 with the cyan dotted lines.</figcaption></figure>

warranted, Fig. 6 seemingly supports ergodicity for volume-integrated properties of ocean turbulence.

**4. Conclusions and discussion**

In this study, we have showcased the ocean energy cycle within the ensemble framework and geopotential coordinates (as opposed to isopycnal coordinates). To our knowledge, our study is novel in that we (i) decompose the mean and eddy energy reservoirs about the ensemble dimension for the ocean energy cycle and (ii) diagnose the potential energy in energy cycles via dynamic enthalpy (DE). The ensemble dimension being orthogonal to the space-time dimensions provides a parameter-free definition of eddy-mean flow decomposition, and preserves the nonstationary nature of the energy pathways, which we have addressed by examining the time series of the energy cycle. While the adoption of DE as potential energy is relatively scarce in the oceanographic literature (e.g., Jamet et al. 2020a), perhaps attributable to its rather computationally intensive nature, we have (i) argued that it is a natural and dynamically consistent extension of gravitational potential energy for a nonlinear equation of state (EOS) and its independence from a reference state of stratification provides a level of

objectivity in how potential energy is defined (Young 2010) and (ii) documented its utility in the energy cycle.

One bewildering aspect, which naturally results from using DE, is that the potential energy reservoir can no longer be split into its mean and eddy components and the mean total dynamic enthalpy (MTDE) reservoir directly interacts with the mean kinetic energy (MKE) and eddy kinetic energy (EKE) reservoirs. This is a stark contrast to Lorenz (1955), where the potential energy reservoir available to eddies is a quadratic term and can be explicitly identified. We argue that this discrepancy is due to the fact that quasi-geostrophy corresponds to the thickness-weighted averaged primitive equations of motion in isopycnal coordinates, and not the unweighted equations in geopotential coordinates (Marshall et al. 2012; Maddison and Marshall 2013; Meunier et al. 2023). Any dynamically consistent quantity resembling APE analogous to the quadratic form in quasi-geostrophy under a nonlinear EOS arises only upon thickness-weighted averaging the governing equations (cf. Aoki 2014; Loose et al. 2023; Uchida et al. 2022b, their appendix A).

By examining the temporal variability of the energy cycle, we have demonstrated that in addition to the well-acknowledged mechanism of baroclinic instability local in space ((w'b') in
<figure><img src="image_8.png" /><figcaption>FIG. 6. (a) Annual-mean, volume-averaged energy cycle diagram of 1967 within the wind-driven gyre based on the ensemble, and (b) temporal framework from a single realization within the ensemble. For the ensemble framework, the annual mean was taken after the energy cycle was diagnosed.</figcaption></figure>

modulating its seasonality (Uchida et al. 2017), the local transfer
from MKE to EKE ($K_{\leftarrow, \mathcal{K}}^{\#}$) and nonlocal eddy transport of
DE fluctuation (($v^\perp h^\perp$)) could also be significant factors in the
western boundary current regions, here, the separated Gulf
Stream (Figs. 4, 5, A3, and A4). In other words, the energy flux
from the potential energy reservoir to EKE is not just the verti-
cal buoyancy flux but the net residual between it and nonlocal
transport of DE. The significance of a nonlocal eddy flux of
DE fluctuation is consistent with recent studies demonstrating
that a nonlocal transport of potential vorticity is crucial for a
proper jet formation in wind-driven gyres (Uchida et al. 2022a;
Deremble et al. 2023). The term $K_{\leftarrow, \mathcal{K}}^{\#}$ being of first-order im-
portance among the energy pathways between energy reser-
voirs is complementary to Jamet et al. (2020a, their Table 2),
where it was exhibited that the energy input to the mean flow by
wind stress was lost to the eddies substantially via barotropic pro-
cesses in the subtropical North Atlantic. It is also consistent with
Kang et al. (2016) and Uchida et al. (2021b), where it was found
that barotropic pathways to the EKE reservoir can overtake bar-
oclinic pathways under increased summertime stratification. This
is, however, slightly at odds with previous studies that attribute
the summer-to-autumn maximum in KE wavenumber spectra
solely to the time lag for winter-to-summer baroclinic instabilities
within the mixed layer to energize EKE and for it to cascade up-
scale (e.g., Sasaki et al. 2014; Uchida et al. 2019; Dong et al. 2020;
Khatri et al. 2021; Steinberg et al. 2022). Subsequently, what they
also find is an increase in inverse KE cascade during late winter.
What is often not documented, however, is the spectral estimates

corresponding to $K_{\leftarrow, \mathcal{K}}^{\#}$. The difference here is that by hav-
ing closed the energy cycle at each moment in time, we are
able to document the instantaneous relative significance of
$K_{\leftarrow, \mathcal{K}}^{\#}$ to the EKE seasonality.

Future work involves (i) extending the time frame of analy-
sis beyond 1967 for a robust seasonal cycle and/or temporal
trend, (ii) investigating how the energy cycle would differ
when the equations of motion are thickness weighted (e.g.,
Bleck 1985; Aiki and Richards 2008; Loose et al. 2023; Uchida
et al. 2022b), and (iii) analyzing ocean ensembles with higher
spatial resolution (currently under production at 1/50° resolu-
tion; Uchida et al. 2023b, their supplemental material) to better
resolve the effect of eddy dynamics and topography on the en-
ergy cycle (Chassignet et al. 2023). Our decomposition about the
ensemble dimension is generic and can be directly applied to at-
mospheric ensemble simulations to diagnose the EKE, MKE,
and large-scale circulation patterns. In the context of climate, our
framework is extendable straightforwardly to fully coupled cli-
mate ensemble simulations (e.g., Maher et al. 2019; Romanou
et al. 2023), which would allow us to quantify the temporally cu-
mulative effect of anthropogenic carbon onto the ocean energy
cycle and integrate it into the climate energy cycle as a whole
(Deser et al. 2020).

*Acknowledgments.* We thank the editor Paola Cessi and two anonymous reviewers for their constructive feedback. This study is a contribution to the “Assessing the Role of forced and internal Variability for the Ocean and climate Response in a changing climate” (ARVOR) project supported by the French ‘Les Enveloppes Fluides et l’Environnement’ (LEFE) program. W. Dewar is supported through NSF Grants OCE-1829856, OCE-1941963 and OCE-2023585, and the French ‘Make Our Planet Great Again’ (MOPGA) program managed by the Agence Nationale de la Recherche under the Programme d’Investissement d’Avenir, reference ANR-18-MPGA-0002. The latter two grants served as the primary support for T. Uchida and partially for Q. Jamet. A. Poje acknowledges support from the NSF Grant OCE-2123633. High-performance computing resources on Cheyenne (https://doi.org/10.5065/D6RX99HX) used for running the Chaocean ensemble of the North Atlantic were provided by NCAR’s Computational and Information Systems Laboratory, sponsored by NSF, under the university large allocation UFSU0011. We thank Edward Peirce and Kelly Hirai for maintaining the Florida State University cluster on which the data were analyzed. Upon getting his hands wet in producing realistic ensemble simulations, T. Uchida realized the formidable dedication it takes in doing so and is indebted to Q. Jamet for generating the Chaocean dataset. The MITgcm outputs were read using the xmitgcm Python package (Abernathey et al. 2022) and postprocessed with the xgcm Python package (Abernathey et al. 2021a).

*Data availability statement.* The simulation outputs are available on the Florida State University cluster (http://ocean.fsu.edu/~qjamet/share/data/Uchida2021/). Jupyter notebooks used for diagnosing the model outputs are available via Github (Uchida 2023).
# APPENDIX

## Energy Cycle of Non-Thickness-Weighted Primitive Equations under Geopotential Coordinates

The ensemble-mean kinetic energy (MKE; $K^# = |\langle \mathbf{u} \rangle|^2/2$) equation is given as

$$ 
\begin{aligned}
\frac{D^#}{Dt} K^# &= -\langle \mathbf{u} \rangle \cdot \nabla_h \langle \phi \rangle - \langle \mathbf{u} \rangle \nabla \cdot \langle \mathbf{v}'\mathbf{u}' \rangle - \langle \mathbf{v} \rangle \nabla \cdot \langle \mathbf{v}'\mathbf{v}' \rangle + \langle \mathbf{u} \rangle \cdot \langle \mathcal{X} \rangle \\
&= -\langle \mathbf{v} \rangle \cdot \nabla \langle \phi \rangle + \langle \mathbf{w} \rangle \langle \mathbf{b} \rangle \\
&\quad - \underbrace{[\nabla \cdot \langle \mathbf{v}'(\mathbf{u} \cdot \mathbf{u}') \rangle - \langle \mathbf{u}'\mathbf{v}' \rangle \cdot \nabla \langle \mathbf{u} \rangle]}_{= \langle \mathbf{u} \rangle \nabla \cdot \langle \mathbf{v}'\mathbf{u}' \rangle + \langle \mathbf{v} \rangle \nabla \cdot \langle \mathbf{v}'\mathbf{v}' \rangle} + \langle \mathbf{u} \rangle \cdot \langle \mathcal{X} \rangle.
\end{aligned}
\quad (A1) $$

The eddy-mean flow interaction term in (A1) is rewritten in the form in square brackets because the divergence component vanishes upon a global volume integration. Figure A1 exhibits some of the terms in the MKE budget. The total kinetic energy (TKE), on the other hand, is

$$ K_t + \mathbf{v} \cdot \nabla K = -\mathbf{u} \cdot \nabla_h \phi + \mathbf{u} \cdot \mathcal{X}. \quad (A2) $$

Now, TKE can be expanded as

$$ 
\begin{aligned}
K &= \frac{1}{2} |\langle \mathbf{u} \rangle + \mathbf{u}'|^2 \\
&= K^# + \mathcal{H} + \langle \mathbf{u} \rangle \cdot \mathbf{u}',
\end{aligned}
\quad (A3) $$

where $\mathcal{H} = |\mathbf{u}'|^2/2$ is the eddy kinetic energy (EKE) so

$$ 
\begin{aligned}
\langle \mathbf{v} \cdot \nabla K \rangle &= \langle (\langle \mathbf{v} \rangle + \mathbf{v}') \cdot \nabla (K^# + \mathcal{H} + \langle \mathbf{u} \rangle \cdot \mathbf{u}') \rangle \\
&= \langle \mathbf{v} \rangle \cdot \nabla K^# + \langle \mathbf{v} \cdot \nabla \mathcal{H} \rangle + \nabla \cdot \langle \mathbf{v}'(\langle \mathbf{u} \rangle \cdot \mathbf{u}') \rangle,
\end{aligned}
\quad (A4) $$

Hence, subtracting (A1) from the ensemble mean of (A2) yields

$$ 
\begin{aligned}
\langle \mathcal{H} \rangle_t = &-\frac{\langle \mathbf{v}' \cdot \nabla \phi' \rangle + \langle \mathbf{w}' \mathbf{b}' \rangle}{\langle \mathbf{v}' \cdot \nabla \mathcal{H} \rangle + \langle \mathbf{u}' \mathbf{v}' \rangle \cdot \nabla \langle \mathbf{u}' \rangle} \\
&-\frac{\langle (\mathbf{v} \cdot \nabla \mathcal{H}) + \langle \mathbf{u}' \mathbf{v}' \rangle \cdot \nabla \langle \mathbf{u}' \rangle \rangle}{\langle \mathbf{u}' \cdot (\mathbf{v} \cdot \nabla \mathbf{u}) \rangle} \\
&+ \langle \mathbf{u}' \cdot \mathcal{X}' \rangle,
\end{aligned}
\quad (A5) $$

where we see the mean flow and eddies exchanging energy via the term $-\langle \mathbf{u}'\mathbf{v}' \rangle \cdot \nabla(\mathbf{u})$, sometimes referred to as the shear-production term in the turbulence literature.

The Lagrangian tendency of dynamic enthalpy is (Young 2010)

$$ h_t + \mathbf{v} \cdot \nabla h = \frac{Dh}{Dt} = -wb + \tilde{h}_\Theta \frac{D\Theta}{Dt} + \tilde{h}_S \frac{DS}{Dt}. \quad (A6) $$

The Lagrangian tendency in the latter two terms on the right-hand side (equivalent to $\Theta$ and $\tilde{S}$, respectively) are, in theory, proportional to molecular and/or nonhydrostatic effects and diabatic forcing. In the main text, $\Theta = S = 0$ is assumed as the necessary outputs to diagnose them were not saved. The equation for mean total dynamic enthalpy (MTDE) becomes

$$ 
\begin{aligned}
\langle h \rangle_t + \langle \mathbf{v} \rangle \cdot \nabla \langle h \rangle = \frac{D^#}{Dt} \langle h \rangle = &-\langle \mathbf{w} \rangle \langle \mathbf{b} \rangle - \langle \mathbf{w}'\mathbf{b}' \rangle - \langle \mathbf{v}' \cdot \nabla h' \rangle \\
&+ \langle \tilde{h}_\Theta \mathring{\Theta} \rangle + \langle \tilde{h}_S \mathring{S} \rangle.
\end{aligned}
\quad (A7) $$

Ensemble averaging (3) yields

$$ \langle \tilde{h} \rangle = \int_{\Phi_0}^{\Phi} \frac{\langle \tilde{b} \rangle}{g} d\Phi^*. \quad (A8) $$

Buoyancy in general is a thermodynamic variable, but the equation of state (EOS) is nonlinear. Thus, we make use of a Taylor expansion as

$$ 
\begin{aligned}
\langle \tilde{b} \rangle &= \langle \tilde{b}((\Theta') + \Theta', \langle S \rangle + S', \Phi) \rangle \\
&= \tilde{b}(\langle \Theta \rangle, \langle S \rangle, \Phi) + \tilde{b}_{\Theta}|_{(\Theta)}(\langle \Theta \rangle, \langle S \rangle, \Phi)\langle \Theta' \rangle \\
&\quad + \tilde{b}_{S}|_{(S)}(\langle \Theta \rangle, \langle S \rangle, \Phi)\langle S' \rangle + \tilde{b}_{\Theta\Theta}|_{(\Theta)}\frac{\langle \Theta'^2 \rangle}{2} \\
&\quad + \tilde{b}_{\Theta S}|_{(\Theta)}\tilde{b}(\langle \Theta' S' \rangle) + \tilde{b}_{SS}|_{(S)}\frac{\langle S'^2 \rangle}{2} + \dots,
\end{aligned}
$$

which argues for

$$ \langle \tilde{b} \rangle = \tilde{\tilde{b}} + \tilde{B}, \quad (A9) $$

where $\tilde{\tilde{b}} \stackrel{\text{def}}{=} \tilde{b}(\langle \Theta \rangle, \langle S \rangle, \Phi)$ and $\tilde{B}$ is at most second order in perturbation because the terms with single-order perturbation vanish, $\langle \Theta' \rangle = \langle S' \rangle = 0$. The term $\tilde{B}$ is only nonzero for a nonlinear EOS and generally represents a small correction. Hence, MTDE becomes

$$ \langle \tilde{h} \rangle = \int_{\Phi_0}^{\Phi} \frac{\tilde{\tilde{b}} + \tilde{B}}{g} d\Phi^* \stackrel{\text{def}}{=} \tilde{H} + \tilde{\Lambda}, \quad (A10) $$

where $\tilde{\Lambda} \stackrel{\text{def}}{=} g^{-1} \int_{\Phi_0}^{\Phi} \tilde{B} d\Phi^*$ shoulders the nonlinear effects and is ignored in (11). Buoyancy fluctuation, on the other hand, can be expanded as

$$ 
\begin{aligned}
b' &= \tilde{b} - \langle \tilde{b} \rangle \\
&= \tilde{b}_{\Theta}|_{(\Theta)}(\langle \Theta \rangle, \langle S \rangle, \Phi)\Theta' + \tilde{b}_{S}|_{(S)}(\langle \Theta \rangle, \langle S \rangle, \Phi)S' \\
&\quad + \tilde{b}_{\Theta\Theta}|_{(\Theta)}\frac{\Theta'^2 - \langle \Theta'^2 \rangle}{2} + \tilde{b}_{\Theta S}|_{(\Theta)}(\Theta' S' - \langle \Theta' S' \rangle) \\
&\quad + \tilde{b}_{SS}|_{(S)}\frac{S'^2 - \langle S'^2 \rangle}{2} + \dots,
\end{aligned}
$$

showing that it is approximated to second order by linear corrections. Thus,

$$ b' = \tilde{b} - (\tilde{\tilde{b}} + \tilde{B}), \quad (A11) $$

$$ h' = \tilde{h} - (\tilde{H} + \tilde{\Lambda}). \quad (A12) $$

In view of (A9)–(A12), (A7) can be rewritten as
<figure><img src="image_9.png" /><figcaption>FIG. A1. (a) The tendency of MKE, (b) mean vertical buoyancy flux reduced by three orders of magnitude, (c) horizontal pressure work, (d) forcing reduced by two orders of magnitude, (e) advection of MKE, and (f) net loss to EKE are shown for 1 Jan 1967. The variables are vertically averaged over the top 1000 m except for the forcing, which only takes nonzero values at the surface.</figcaption></figure>

$$ \begin{aligned} & \overline{\tilde{H}_{\Phi} \frac{D^{\#}}{Dt} \Phi + \tilde{H}_{\Theta}|_{(\Theta)} \frac{D^{\#}}{Dt} \langle \Theta \rangle + \tilde{H}_{S}|_{\langle S \rangle} \frac{D^{\#}}{Dt} \langle S \rangle + \frac{D^{\#}}{Dt} \Lambda} \\ &= \frac{D^{\#}}{Dt} H \\ &= -\langle w \rangle \langle b \rangle - \langle w' b' \rangle - \langle v' \cdot \nabla (h - H - \Lambda) \rangle + \langle \tilde{h}_{\Theta} \mathring{\Theta} \rangle + \langle \tilde{h}_{S} \mathring{S} \rangle \\ &= -\langle w \rangle (\bar{b} + \mathcal{B}) - \langle w' (b - \bar{b} - \mathcal{B}) \rangle - \langle v' \cdot \nabla (h - H - \Lambda) \rangle \\ &\quad + \langle \tilde{h}_{\Theta} \mathring{\Theta} \rangle + \langle \tilde{h}_{S} \mathring{S} \rangle, \end{aligned} $$

$$ \begin{aligned} & \therefore \tilde{H}_{\Theta}|_{(\Theta)} \frac{D^{\#}}{Dt} \langle \Theta \rangle + \tilde{H}_{S}|_{\langle S \rangle} \frac{D^{\#}}{Dt} \langle S \rangle + \frac{D^{\#}}{Dt} \Lambda \\ &= -\langle w \rangle \mathcal{B} - \langle w' (b - \bar{b} - \mathcal{B}) \rangle - \langle v' \cdot \nabla (h - H - \Lambda) \rangle \\ &\quad + \langle \tilde{h}_{\Theta} \mathring{\Theta} \rangle + \langle \tilde{h}_{S} \mathring{S} \rangle, \end{aligned} \quad (A13) $$

where we have used $H_{\Phi}(D^{\#}/Dt)\Phi = -\langle w \rangle \bar{b}$. Realizing that $D^{\#}/Dt(\Theta) = -\langle v' \cdot \nabla \Theta' \rangle + \langle \mathring{\Theta} \rangle$ and $D^{\#}/Dt(S) = -\langle v' \cdot \nabla S' \rangle + \langle \mathring{S} \rangle$ yields
$$ \begin{aligned} \langle w'(b - \bar{b} - \mathcal{B}) \rangle = & -\frac{D^{\#}}{Dt} \Lambda - \langle w \rangle \mathcal{B} - \langle \mathbf{v}' \cdot \nabla (h - H - \Lambda) \rangle \\ & + \tilde{H}_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle \\ & + (\langle \tilde{h}_{\Theta} \mathring{\Theta} \rangle - \tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathring{\Theta} \rangle) \\ & + (\langle \tilde{h}_{S} \mathring{S} \rangle - \tilde{H}_{S}|_{\langle S \rangle} \langle \mathring{S} \rangle), \end{aligned} \quad (A14) $$

Eq. (A14) allows us to unify the EKE and MTDE equations:

$$ \begin{aligned} \frac{D^{\#}}{Dt} (\langle \mathcal{H} \rangle + \Lambda) = & - \langle \mathbf{v}' \cdot \nabla [ (\phi - \bar{\phi}) - \mathcal{P} ) + \mathcal{H} + (h - H - \Lambda) ] \rangle \\ & - \langle \mathbf{v}' \mathbf{u}' \rangle \cdot \nabla \langle \mathbf{u} \rangle + \langle \mathbf{u}' \cdot \mathcal{X}' \rangle - \langle w \rangle \mathcal{B} \\ & + \tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{S}|_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle \\ & + (\langle \tilde{h}_{\Theta} \mathring{\Theta} \rangle - \tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathring{\Theta} \rangle) + (\langle \tilde{h}_{S} \mathring{S} \rangle \\ & - \tilde{H}_{S}|_{\langle S \rangle} \langle \mathring{S} \rangle), \end{aligned} \quad (A15) $$

where $\mathcal{P}_z \stackrel{\text{def}}{=} \mathcal{B}$. On the other hand, the unified MKE and MTDE equation becomes

$$ \begin{aligned} \frac{D^{\#}}{Dt} (K^{\#} + \underbrace{\langle h \rangle}_{=H} - \Lambda) = & - \langle \mathbf{v} \rangle \cdot \nabla (\bar{\phi} + \mathcal{P}) - [\nabla \cdot \langle \mathbf{v}' (\mathbf{u} \cdot \mathbf{u}') \rangle] \\ & - \langle \mathbf{v}' \mathbf{u}' \rangle \cdot \nabla \langle \mathbf{u} \rangle] + \langle \mathbf{u} \rangle \cdot \langle \mathcal{X} \rangle + \langle w \rangle \mathcal{B} \\ & - (\tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathbf{v}' \cdot \nabla \Theta' \rangle + \tilde{H}_{S}|_{\langle S \rangle} \langle \mathbf{v}' \cdot \nabla S' \rangle) \\ & + \tilde{H}_{\Theta}|_{\langle \Theta \rangle} \langle \mathring{\Theta} \rangle + \tilde{H}_{S}|_{\langle S \rangle} \langle \mathring{S} \rangle. \end{aligned} \quad (A16) $$

It appears that there is an additional energy reservoir stemming from the nonlinearities in EOS but we remind the reader that an evolution equation for $\tilde{\Lambda}$ (and $\tilde{\mathcal{B}}$) cannot be formulated based on $b'$ because $\tilde{b} - \langle \tilde{b} \rangle \neq \tilde{b} - \bar{b}$, Furthermore, the quadratic terms in $\tilde{\mathcal{B}}$ (e.g., $\tilde{b}_{\Theta\Theta}|_{\langle \Theta \rangle} (\langle \Theta^2 \rangle / 2)$) may appear analogous to APE but under a linear EOS, APE would be proportional to $\langle (\tilde{b}_{\Theta\Theta}|_{\langle \Theta \rangle})^2 \rangle$; $\tilde{b}_{\Theta\Theta}|_{\langle \Theta \rangle} (\langle \Theta^2 \rangle / 2)$ is a term tapping into the internal energy of thermodynamics and is generally much smaller than $\langle (\tilde{b}_{\Theta}|_{\langle \Theta \rangle})^2 \rangle$.

In practice, we neglect all the second- and higher-order correction terms and $b' \approx \tilde{b}^{\dagger} \stackrel{\text{def}}{=} \tilde{b} - \bar{b}$ in the main text. Subsequently, we adopt $\langle \tilde{h} \rangle \approx \tilde{H}$ and $h' \approx \tilde{h}^{\dagger} \stackrel{\text{def}}{=} \tilde{h} - \bar{H}$. The latter is used to diagnose $\nabla \cdot \langle \mathbf{v}' \mathbf{h}' \rangle$. We show in Fig. A2

<figure><img src="image_10.png" /><figcaption>FIG. A2. Comparison between buoyancy fluctuation $\tilde{b}'(=\tilde{b} - \bar{\tilde{b}})$ and $\tilde{b}_{\Theta}|_{\langle \Theta \rangle} (\langle \Theta \rangle, \langle S \rangle, \Phi)\Theta' + \tilde{b}_{S}|_{\langle S \rangle} (\langle \Theta \rangle, \langle S \rangle, \Phi)S'$. The partial derivatives respective to mean potential temperature and practical salinity were taken using the fastjmd95 Python package (Abernathey and Busecke 2020). (a) The buoyancy fluctuation $b'$, (b) $\tilde{b}_{\Theta}|_{\langle \Theta \rangle} (\langle \Theta \rangle, \langle S \rangle, \Phi)\Theta' + \tilde{b}_{S}|_{\langle S \rangle} (\langle \Theta \rangle, \langle S \rangle, \Phi)S'$, (c) their difference in percentage on 1 Jan 1967 at $z = -270$ m, and (d) a joint histogram of the former two throughout 1967 over the entire three-dimensional domain.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>FIG. A3. Comparison of approximations to $H_{\to \infty}$ on 1 Jan 1967. (a) The vertical eddy buoyancy flux, (b) convergence of eddy dynamic enthalpy flux, (c) contribution due to convergence of eddy temperature and salinity flux, (d) difference between the left- and right-hand side of (12) in percentage, (e) net eddy buoyancy flux all at $z = -270$ m, and (f) a joint histogram over the entire three-dimensional domain. Panels (a)–(c) and (e) are plotted against the same color bar. The ratio in (d) can be ill defined where the dominator is small.</figcaption></figure>

that $b'$ approximated by $\tilde{b} - \bar{\tilde{b}}$ and $\tilde{b}_{\Theta}|_{\Theta}(\langle\Theta\rangle, \langle S\rangle, \Phi)\Theta' + \tilde{b}_S|_{S}(\langle\Theta\rangle, \langle S\rangle, \Phi)S'$ are nearly identical with each other across the year throughout the basin; the differences are largely contained within the separated Gulf Stream region. Equation (A14) can be simplified by dropping the terms due to diabatic nature and nonlinearity in EOS, which leaves us with (12). Figure A3 demonstrates that (12) holds surprisingly well; we show

a histogram exhibiting that the two align mostly along a one-to-one line (Fig. A3f). The end result of neglecting the nonlinearity in the EOS and terms associated with diabatic mixing in (A15) and (A16) leaves us with (13), (14), and (15). There are two sources for eddy energy, one from the MKE and another from the MTDE field. The generalized pressure work tends to counteract the input from MTDE (Figs. A3 and A4).
<figure><img src="image_12.png" /><figcaption>FIG. A4. (a) The tendency of EKE, (b) energy influx from MTDE, (c) generalized pressure work, (d) forcing, (e) advective of EKE, and (f) shear production are shown for 1 Jan 1967. The variables are vertically averaged over the top 1000 m except for the forcing, which only takes nonzero values at the surface. The influx from MTDE and generalized pressure work tend to counteract each other. —[⟨v·∇·**F**⟩ in (e) is diagnosed as —[⟨**u**′·(**v**·∇**u**)'⟩ + (⟨**u**′**v**'⟩·∇⟨**u**⟩ + ⟨**v**′**v**'⟩·∇⟨**v**⟩)]].</figcaption></figure>

REFERENCES

Abernathey, R. P., and J. Busecke, 2020: *fastjmd95: Numba implementation of Jackett & McDougall (1995) ocean equation of state [Software]*. Zenodo, accessed 2 February 2024, https://github.com/xgcm/fastjmd95.

—, and Coauthors, 2021a: xgcm: General circulation model post-processing with xarray [Software]. Zenodo, accessed 2 February 2024, https://xgcm.readthedocs.io/en/latest/.

—, and Coauthors, 2021b: xhistogram: Fast, flexible, label-aware histograms for numpy and xarray [Software]. Zenodo, accessed 2 February 2024, https://xhistogram.readthedocs.io/en/latest/.

—, and Coauthors, 2022: xmitgcm: Read MITgcm mds binary files into xarray [Software]. Zenodo, accessed 2 February 2024, https://github.com/MITgcm/xmitgcm.

Aiki, H., and K. J. Richards, 2008: Energetics of the global ocean: The role of layer-thickness form drag. *J. Phys. Oceanogr.*, **38**, 1845–1869, https://doi.org/10.1175/2008JPO3820.1.
—, X. Zhai, and R. J. Greatbatch, 2016: Energetics of the global ocean: The role of mesoscale eddies. *Indo-Pacific Climate Variability and Predictability*, S. K. Behera and T. Yamagata, Eds., World Scientific Series on Asia-Pacific Weather and Climate, Vol. 7, 109–134, https://doi.org/10.1142/9789814696623_0004.

Aoki, K., 2014: A constraint on the thickness-weighted average equation of motion deduced from energetics. *J. Mar. Res.*, **72**, 355–382, https://doi.org/10.1357/002224014815469886.

—, A. Kubokawa, R. Furue, and H. Sasaki, 2016: Influence of eddy momentum fluxes on the mean flow of the Kuroshio Extension in a 1/10 ocean general circulation model. *J. Phys. Oceanogr.*, **46**, 2769–2784, https://doi.org/10.1175/JPO-D-16-0021.1.

—, Y. Miyazawa, T. Hihara, and T. Miyama, 2020: An objective method for probabalistic forecasting of multimodal Kuroshio states using ensemble simulation and machine learning. *J. Phys. Oceanogr.*, **50**, 3189–3204, https://doi.org/10.1175/JPO-D-19-0316.1.

Arbic, B. K., K. L. Polzin, R. B. Scott, J. G. Richman, and J. F. Shriver, 2013: On eddy viscosity, energy cascades, and the horizontal resolution of gridded satellite altimeter products. *J. Phys. Oceanogr.*, **43**, 283–300, https://doi.org/10.1175/JPO-D-11-0240.1.

Bachman, S. D., B. Fox-Kemper, and F. O. Bryan, 2015: A tracer-based inversion method for diagnosing eddy-induced diffusivity and advection. *Ocean Modell.*, **8**, 1–14, https://doi.org/10.1016/j.ocemod.2014.11.006.

Bleck, R., 1985: On the conversion between mean and eddy components of potential and kinetic energy in isentropic and isopycnic coordinates. *Dyn. Atmos. Oceans*, **9**, 17–37, https://doi.org/10.1016/0377-0265(85)90014-4.

Buzzicotti, M., B. Storer, H. Khatri, S. Griffies, and H. Aluie, 2023: Spatio-temporal coarse-graining decomposition of the global ocean geostrophic kinetic energy. *J. Adv. Model. Earth Syst.*, **15**, e2023MS003693, https://doi.org/10.1029/2023MS003693.

Charney, J. G., and Coauthors, 1979: *Carbon Dioxide and Climate: A Scientific Assessment*. National Academy of Sciences, 22 pp.

Chassignet, E. P., X. Xu, A. Bozec, and T. Uchida, 2023: Impact of the New England seamount chain on Gulf Stream pathway and variability. *J. Phys. Oceanogr.*, **53**, 1871–1886, https://doi.org/10.1175/JPO-D-23-0008.1.

Chen, R., 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.

Contreras, M., L. Renault, and P. Marchesiello, 2023: Understanding energy pathways in the Gulf Stream. *J. Phys. Oceanogr.*, **53**, 719–736, https://doi.org/10.1175/JPO-D-22-0146.1.

Demyshev, S. G., and O. A. Dymova, 2022: Analysis of the annual mean energy cycle of the Black Sea circulation for the climatic, basin-scale and eddy regimes. *Ocean Dyn.*, **72**, 259–278, https://doi.org/10.1007/s10236-022-01504-0.

Deremble, B., N. Wienders, and W. K. Dewar, 2013: CheapAML: A simple, atmospheric boundary layer model for use in ocean-only model calculations. *Mon. Wea. Rev.*, **141**, 809–821, https://doi.org/10.1175/MWR-D-11-00254.1.

—, T. Uchida, W. K. Dewar, and R. M. Samelson, 2023: Eddy-mean flow interaction with a multiple scale quasi geostrophic model. *J. Adv. Model. Earth Syst.*, **15**, e2022MS003572, https://doi.org/10.1029/2022MS003572.

Deser, C., and Coauthors, 2020: Insights from Earth system model initial-condition large ensembles and future prospects. *Nat. Climate Change*, **10**, 277–286, https://doi.org/10.1038/s41558-020-0731-2.

Dewar, W. K., J. Schoonover, T. McDougall, and R. Klein, 2016: Semicompressible ocean thermodynamics and Boussinesq energy conservation. *Fluids*, **1**, 9, https://doi.org/10.3390/fluids1020009.

—, R. Parfitt, and N. Wienders, 2022: Routine reversal of the AMOC in an ocean model ensemble. *Geophys. Res. Lett.*, **49**, e2022GL100117, https://doi.org/10.1029/2022GL100117.

Dong, J., B. Fox-Kemper, H. Zhang, and C. Dong, 2020: The seasonality of submesoscale energy production, content, and cascade. *Geophys. Res. Lett.*, **47**, e2020GL087388, https://doi.org/10.1029/2020GL087388.

Eden, C., 2015: Revisiting the energetics of the ocean in Boussinesq approximation. *J. Phys. Oceanogr.*, **45**, 630–637, https://doi.org/10.1175/JPO-D-14-0072.1.

Fairall, C. W., E. F. Bradley, J. E. Hare, A. A. Grachev, and J. B. Edson, 2003: Bulk parameterization of air-sea fluxes: Updates and verification for the COARE algorithm. *J. Climate*, **16**, 571–591, https://doi.org/10.1175/1520-0442(2003)016<0571: BPOASF>2.0.CO;2.

Gnanadesikan, A., M.-A. Pradal, and R. Abernathey, 2015: Iso-pycnal mixing by mesoscale eddies significantly impacts oceanic anthropogenic carbon uptake. *Geophys. Res. Lett.*, **42**, 4249–4255, https://doi.org/10.1002/2015GL064100.

Griffies, S. M., and Coauthors, 2015: Impacts on ocean heat from transient mesoscale eddies in a hierarchy of climate models. *J. Climate*, **28**, 952–977, https://doi.org/10.1175/JCLI-D-14-00353.1.

Guo, Y., S. Bishop, F. Bryan, and S. Bachman, 2022: A global diagnosis of eddy potential energy budget in an eddy-permitting ocean model. *J. Phys. Oceanogr.*, **52**, 1731–1748, https://doi.org/10.1175/JPO-D-22-0029.1.

Hartmann, D. L., V. Ramanathan, A. Berroir, and G. E. Hunt, 1986: Earth radiation budget data and climate research. *Rev. Geophys.*, **24**, 439–468, https://doi.org/10.1029/RG024i002p00439.

Hersbach, H., C. Peubey, A. Simmons, P. Berrisford, P. Poli, and D. Dee, 2015: ERA-20CM: A twentieth-century atmospheric model ensemble. *Quart. J. Roy. Meteor. Soc.*, **141**, 2350–2375, https://doi.org/10.1002/qj.2528.

Jackett, D. R., and T. J. McDougall, 1995: Minimal adjustment of hydrographic profiles to achieve static stability. *J. Atmos. Oceanic Technol.*, **12**, 381–389, https://doi.org/10.1175/1520-0426 (1995)012<0381:MAOHPT>2.0.CO;2.

Jamet, Q., W. K. Dewar, N. Wienders, and B. Deremble, 2019a: Fast warming of the surface ocean under a climatological scenario. *Geophys. Res. Lett.*, **46**, 3871–3879, https://doi.org/10.1029/2019GL082336.

—, —, —, and —, 2019b: Spatio-temporal patterns of chaos in the Atlantic overturning circulation. *Geophys. Res. Lett.*, **46**, 7509–7517, https://doi.org/10.1029/2019GL082552.

—, A. Ajayi, J. Le Sommer, T. Penduff, A. Hogg, and W. Dewar, 2020a: On energy cascades in general flows: A Lagrangian application. *J. Adv. Model. Earth Syst.*, **12**, e2020MS002090, https://doi.org/10.1029/2020MS002090.

—, W. K. Dewar, N. Wienders, B. Dersemble, S. Close, and T. Penduff, 2020b: Locally and remotely forced subtropical AMOC variability: A matter of time scales. *J. Climate*, **33**, 5155–5172, https://doi.org/10.1175/JCLI-D-19-0844.1.

—, B. Dersemble, N. Wienders, T. Uchida, and W. K. Dewar, 2021: On wind-driven energetics of subtropical gyres. *J. Adv.*
Model. Earth Syst., **13**, e2020MS002329, https://doi.org/10.10
29/2020MS002329.

—, S. Leroux, W. K. Dewar, T. Penduff, J. Le Sommer, J.-M.
Molines, and J. Gula, 2022: Non-local eddy-mean kinetic en-
ergy transfers in submesoscale-permitting ensemble simula-
tions. *J. Adv. Model. Earth Syst.*, **14**, e2022MS003057, https://
doi.org/10.1029/2022MS003057.

Kang, D., and E. N. Curchitser, 2015: Energetics of eddy–mean flow interactions in the Gulf Stream region. *J. Phys. Oceanogr.*, **45**, 1103–1120, https://doi.org/10.1175/JPO-D-14-0200.1.

—, —, and A. Rosati, 2016: Seasonal variability of the Gulf Stream kinetic energy. *J. Phys. Oceanogr.*, **46**, 1189–1207, https://doi.org/10.1175/JPO-D-15-0235.1.

Khatri, H., S. M. Griffies, T. Uchida, H. Wang, and D. Menemenlis, 2021: Role of mixed-layer instabilities in the seasonal evolution of eddy kinetic energy spectra in a global submesoscale permitting simulation. *Geophys. Res. Lett.*, **48**, e2021GL094777, https://doi.org/10.1029/2021GL094777.

Knutti, R., M. A. A. Rugenstein, and G. C. Hegerl, 2017: Beyond equilibrium climate sensitivity. *Nat. Geosci.*, **10**, 727–736, https://doi.org/10.1038/ngeo3017.

Large, W. G., J. C. McWilliams, and S. C. Doney, 1994: Oceanic
vertical mixing: A review and a model with a nonlocal
boundary layer parameterization. *Rev. Geophys.*, **32**, 363–403,
https://doi.org/10.1029/94RG01872.

Lenderink, G., A. van Ulden, B. van den Hurk, and E. van Meij-
gaard, 2007: Summertime inter-annual temperature variability
in an ensemble of regional model simulations: Analysis of the
surface energy budget. *Climatic Change*, **81**, 233–247, https://
doi.org/10.1007/s10584-006-9229-9.

Leroux, S., T. Penduff, L. Bessières, J.-M. Molines, J.-M. Bran-
kart, G. Sérazin, B. Barnier, and L. Terray, 2018: Intrinsic
and atmospherically forced variability of the AMOC: Insights
from a large-ensemble ocean hindcast. *J. Climate*, **31**, 1183–
1203, https://doi.org/10.1175/JCLI-D-17-0168.1.

Loose, N., S. Bachman, I. Grooms, and M. Jansen, 2023: Diagnos-
ing scale-dependent energy cycles in a high-resolution isopyc-
nal ocean model. *J. Phys. Oceanogr.*, **53**, 157–176, https://doi.
org/10.1175/JPO-D-22-0083.1.

Lorenz, E. N., 1955: Available potential energy and the mainte-
nance of the general circulation. *Tellus*, **7**, 157–167, https://
doi.org/10.3402/tellusa.v7i2.8796.

—, 1963: Deterministic nonperiodic flow. *J. Atmos. Sci.*, **20**,
130–141, https://doi.org/10.1175/1520-0469(1963)020<0130:
DNF>2.0.CO;2.

Maddison, J. R., and D. P. Marshall, 2013: The Eliassen–Palm flux tensor. *J. Fluid Mech.*, **729**, 69–102, https://doi.org/10.1017/jfm.2013.259.

Maher, N., and Coauthors, 2019: The Max Planck Institute Grand Ensemble: Enabling the exploration of climate system variability. *J. Adv. Model. Earth Syst.*, **11**, 2050–2069, https://doi.org/10.1029/2019MS001639.

Marshall, D. P., J. R. Maddison, and P. S. Berloff, 2012: A framework for parameterizing eddy potential vorticity fluxes. *J. Phys. Oceanogr.*, **42**, 539–557, https://doi.org/10.1175/JPO-D-11-048.1.

Marshall, J., C. Hill, L. Perelman, and A. Adcroft, 1997: Hydro-
static, quasi-hydrostatic and non-hydrostatic ocean modelling.
*J. Geophys. Res.*, **102**, 5733–5752, https://doi.org/10.1029/96JC
02776.

Masson-Delmotte, V., and Coauthors, 2021: *Climate Change 2021:*
*The Physical Science Basis*. Cambridge University Press, 2391
pp., https://doi.org/10.1017/9781009157896.

Meunier, J., B. Miquel, and B. Gallet, 2023: A direct derivation of the Gent-McWilliams/Redi diffusion tensor from quasi-geostrophic dynamics. *J. Fluid Mech.*, **963**, A22, https://doi.org/10.1017/jfm.2023.347.

Molemaker, M. J., and J. C. McWilliams, 2010: Local balance and cross-scale flux of available potential energy. *J. Fluid Mech.*, **645**, 295–314, https://doi.org/10.1017/S0022112009992643.

—, —, and X. Capet, 2010: Balanced and unbalanced routes to dissipation in an equilibrated Eady flow. *J. Fluid Mech.*, **654**, 35–63, https://doi.org/10.1017/S0022112009993272.

Nikiéma, O., and R. Laprise, 2013: An approximate energy cycle for inter-member variability in ensemble simulations of a regional climate model. *Climate Dyn.*, **41**, 831–852, https://doi.org/10.1007/s00382-012-1575-x.

Nycander, J., 2010: Horizontal convection with a non-linear equa-
tion of state: Generalization of a theorem of Paparella and
Young. *Tellus*, **62A**, 134–137, https://doi.org/10.1111/j.1600-08
70.2009.00429.x.

Oort, A. H., S. C. Ascher, S. Levitus, and J. P. Peixóto, 1989: New estimates of the available potential energy in the world ocean. *J. Geophys. Res.*, **94**, 3187–3200, https://doi.org/10.1029/JC094iC03p03187.

Poincaré, H., 1890: On the three-body problem and the equations of dynamics. *Acta Math.*, **13**, 1–270.

Renault, L., M. J. Molemaker, J. Gula, S. Masson, and J. C. McWilliams, 2016: Control and stabilization of the Gulf Stream by oceanic current interaction with the atmosphere. *J. Phys. Oceanogr.*, **46**, 3439–3453, https://doi.org/10.1175/JPO-D-16-0115.1.

—, S. Masson, V. Oerder, F. Colas, and J. McWilliams, 2023: Modulation of the oceanic mesoscale activity by the meso-scale thermal feedback to the atmosphere. *J. Phys. Oceanogr.*, **53**, 1651–1667, https://doi.org/10.1175/JPO-D-22-0256.1.

Romanou, A., and Coauthors, 2023: Stochastic bifurcation of the North Atlantic circulation under a mid-range future climate scenario with the NASA-GISS ModelE. *J. Climate*, **36**, 6141–6161, https://doi.org/10.1175/JCLI-D-22-0536.1.

Saenz, J. A., R. Tailleux, E. D. Butler, G. O. Hughes, and K. I. Oliver, 2015: Estimating Lorenz's reference state in an ocean with a nonlinear equation of state for seawater. *J. Phys. Oceanogr.*, **45**, 1242–1257, https://doi.org/10.1175/JPO-D-14-0105.1.

Sasaki, H., P. Klein, B. Qiu, and Y. Sasai, 2014: Impact of oce-
anic-scale interactions on the seasonal modulation of ocean
dynamics by the atmosphere. *Nat. Commun.*, **5**, 5636, https://
doi.org/10.1038/ncomms6636.

Sérazin, G., and Coauthors, 2017: A global probabilistic study of the ocean heat content low-frequency variability: Atmospheric forcing versus oceanic chaos. *Geophys. Res. Lett.*, **44**, 5580–5589, https://doi.org/10.1002/2017GL073026.

Sherwood, S. C., and Coauthors, 2020: An assessment of Earth's climate sensitivity using multiple lines of evidence. *Rev. Geophys.*, **58**, e2019RG000678, https://doi.org/10.1029/2019RG000678.

Steinberg, J. M., and C. C. Eriksen, 2022: Eddy vertical structure and variability: Deepglider observations in the North Atlantic. *J. Phys. Oceanogr.*, **52**, 1091–1110, https://doi.org/10.1175/JPO-D-21-0068.1.

—, S. T. Cole, K. Drushka, and R. P. Abernathey, 2022: Seasonality of the mesoscale inverse cascade as inferred from global scale-dependent eddy energy observations. *J. Phys. Oceanogr.*, **52**, 1677–1691, https://doi.org/10.1175/JPO-D-21-0269.1.

Sui, C. H., K. M. Lau, W. K. Tao, and J. Simpson, 1994: The tropi-
ical water and energy cycles in a cumulus ensemble model.
Part I: Equilibrium climate. *J. Atmos. Sci.*, **51**, 711–728, https://doi.org/10.1175/1520-0469(1994)051<0711:TTWAEC>2.0.CO;2.

Tailleux, R., 2013: Available potential energy and exergy in strati-fied fluids. *Annu. Rev. Fluid Mech.*, **45**, 35–58, https://doi.org/10.1146/annurev-fluid-011212-140620.

—, 2016: Generalized patched potential density and thermody-namic neutral density: Two new physically based quasi-neutral density variables for ocean water masses analyses and circula-tion studies. *J. Phys. Oceanogr.*, **46**, 3571–3584, https://doi.org/10.1175/JPO-D-16-0072.1.

—, and G. Wolf, 2023: On the links between neutral directions, buoyancy forces, energetics, potential vorticity, and lateral stirring in the ocean: A first-principles approach. arXiv, 2202.00456v4, https://doi.org/10.48550/arXiv.2202.00456.

Uchida, T., 2023: Energy_Zycle: Jupyter repository for analyzing the Chaocean energy cycle [Software]. Zenodo, accessed 2 February 2024, https://github.com/roxyboy/Energy_Zycle.

—, R. P. Abernathey, and S. Smith, 2017: Seasonality of eddy kinetic energy in an eddy permitting global climate model. *Ocean Modell.*, **118**, 41–58, https://doi.org/10.1016/j.ocemod.2017.08.006.

—, D. Balwada, R. Abernathey, G. McKinley, S. Smith, and M. Lévy, 2019: The contribution of submesoscale over mesoscale eddy iron transport in the open Southern Ocean. *J. Adv. Model. Earth Syst.*, **11**, 3934–3958, https://doi.org/10.1029/2019MS001805.

—, B. Deremble, W. K. Dewar, and T. Penduff, 2021a: Diag-nosing the Eliassen-Palm flux from a quasi-geostrophic dou-ble gyre ensemble. EarthCube Annual Meeting, Online, NSF, https://earthcube2021.github.io/ec21_book/notebooks/ec21_uchida_etal/notebooks/TU_05_Diagnosing-the-Eliassen-Palm-flux-from-a-quasi-geostrophic-double-gyre-ensemble.html.

—, —, and T. Penduff, 2021b: The seasonal variability of the ocean energy cycle from a quasi-geostrophic double gyre en-semble. *Fluids*, **6**, 206, https://doi.org/10.3390/fluids6060206.

—, —, and S. Popinet, 2022a: Deterministic model of the eddy dynamics for a midlatitude ocean model. *J. Phys. Oceanogr.*, **52**, 1133–1154, https://doi.org/10.1175/JPO-D-21-0217.1.

—, Q. Jamet, W. K. Dewar, D. Balwada, J. Le Sommer, and T. Penduff, 2022b: Diagnosing the thickness-weighted averaged eddy-mean flow interaction in an eddying North Atlantic en-semble: The Eliassen-Palm flux. *J. Adv. Model. Earth Syst.*, **14**, e2021MS002866, https://doi.org/10.1029/2021MS002866.

—, —, A. Poje, and W. K. Dewar, 2022c: An ensemble-based eddy and spectral analysis, with application to the Gulf Stream. *J. Adv. Model. Earth Syst.*, **14**, e2021MS002692, https://doi.org/10.1029/2021MS002692.

—, D. Balwada, Q. Jamet, W. K. Dewar, B. Dersemble, T. Penduff, and J. Le Sommer, 2023a: Cautionary tales from the mesoscale eddy transport tensor. *Ocean Modell.*, **182**, 102172, https://doi.org/10.1016/j.ocemod.2023.102172.

—, Q. Jamet, A. C. Poje, N. Wienders, and W. K. Dewar, 2023b: Wavelet-based wavenumber spectral estimate of eddy kinetic energy: Application to the North Atlantic. *Earth-Ativ*, https://doi.org/10.31223/X5036Q, preprint.

Verhulst, P.-F., 1845: Recherches mathématiques sur la loi d'accroissement de la population. *Mem. Acad. Roy. Belg.*, **18** (1), 1–40, https://doi.org/10.3406/marb.1845.3438.

von Storch, J.-S., C. Eden, I. Fast, H. Haak, D. Hernández-Deckers, E. Maier-Reimer, J. Marotzke, and D. Stammer, 2012: An estimate of the Lorenz energy cycle for the world ocean based on the STORM/NCEP simulation. *J. Phys. Oceanogr.*, **42**, 2185–2205, https://doi.org/10.1175/JPO-D-12-079.1.

Winters, K. B., P. N. Lombard, J. J. Riley, and E. A. D'Asaro, 1995: Available potential energy and mixing in density-stratified fluids. *J. Fluid Mech.*, **289**, 115–128, https://doi.org/10.1017/S002211209500125X.

Wunsch, C., 1981: Low frequency variability of the sea. *The Evolution of Physical Oceanography: Scientific Surveys in Honor of Henry Stommel*, C. Wunsch and B. A. Warren, Eds., MIT Press, 342–374.

—, and R. Ferrari, 2004: Vertical mixing, energy, and the general circulation of the oceans. *Annu. Rev. Fluid Mech.*, **36**, 281–314, https://doi.org/10.1146/annurev.fluid.36.050802.122121.

Xie, J., H. Liu, and P. Lin, 2023: A multifaceted isoneutral eddy transport diagnostic framework and its application in the Southern Ocean. *J. Adv. Model. Earth Syst.*, **15**, e2023MS003728, https://doi.org/10.1029/2023MS003728.

Xie, S.-P., Y. Kosaka, and Y. M. Okumura, 2016: Distinct energy budgets for anthropogenic and natural changes during global warming hiatus. *Nat. Geosci.*, **9**, 29–33, https://doi.org/10.1038/ngeo2581.

Yang, P., Z. Jing, H. Wang, L. Wu, Y. Chen, and S. Zhou, 2022: Role of frictional processes in mesoscale eddy available poten-tial energy budget in the global ocean. *Geophys. Res. Lett.*, **49**, e2021GL097557, https://doi.org/10.1029/2021GL097557.

Yang, Y., R. Chen, G. R. Flierl, and H. Zhang, 2023: A diagnostic framework linking eddy flux ellipse with eddy-mean energy ex-change. *Res. Square*, https://doi.org/10.21203/rs.3.rs-2954176/v1, preprint.

Young, W. R., 2010: Dynamic enthalpy, conservative temperature, and the seawater Boussinesq approximation. *J. Phys. Oceanogr.*, **40**, 394–400, https://doi.org/10.1175/2009JPO4294.1.

Zhai, X., R. J. Greatbatch, and J.-D. Kohlmann, 2008: On the seasonal variability of eddy kinetic energy in the Gulf Stream region. *Geophys. Res. Lett.*, **35**, L24609, https://doi.org/10.1029/2008GL036412.

Zhao, M., R. M. Ponte, and T. Penduff, 2023: Global-scale random bottom pressure fluctuations from oceanic intrinsic variability. *Sci. Adv.*, **9**,eadg0278, https://doi.org/10.1126/sciadv.adg0278.