# Understanding Energy Pathways in the Gulf Stream

MARCELA CONTRERAS<sup>○</sup>,<sup>a</sup> LIONEL RENAULT,<sup>a</sup> AND PATRICK MARCHESIELLO<sup>a</sup>

<sup>a</sup> LEGOS, Université de Toulouse, CNES-CNRS-IRD-UPS, Toulouse, France

(Manuscript received 13 July 2022, in final form 7 November 2022)

**ABSTRACT:** The Gulf Stream (GS) is one of the strongest ocean currents on the planet. Eddy-rich resolution models are needed to properly represent the dynamics of the GS; however, kinetic energy (KE) can be in excess in these models if not dissipated efficiently. The question of how and how much energy is dissipated and in particular how it flows through ocean scales thus remains an important and largely unanswered question. Using a high-resolution (~2 km) ocean model [Coastal and Regional Ocean Community (CROCO)], we characterize the spatial and temporal distribution of turbulent cascades in the GS based on a coarse-grained method. We show that the balanced flow is associated with an inverse cascade while the forward cascade is explained by ageostrophic advection associated with frontogenesis. Downscale fluxes are dominant at scales smaller than about 20 km near the surface and most intense at the GS North Wall. There is also strong seasonal variability in KE flux, with the forward cascade intensifying in winter and the inverse cascade later in spring. The forward cascade, which represents an interior route to dissipation, is compared with both numerical and boundary dissipation processes. The contribution of interior dissipation is an order of magnitude smaller than that of the other energy sinks. We thus evaluate the sensitivity of horizontal momentum advection schemes on energy dissipation and show that the decrease in numerical dissipation in a high-order scheme leads to an increase in dissipation at the boundaries, not in the downscale flux.

**KEYWORDS:** Nonlinear dynamics; Small-scale processes; Ocean models

## 1. Introduction

The Gulf Stream is a strong, deep, and persistent western boundary current. It is characterized by a separation near Cape Hatteras and penetration into the North Atlantic basin, where it acts as the main conduit for the upper branch of the meridional overturning circulation. According to Sverdrup's linear theory, the volume transport of western boundary currents should vary with the intensity of the basinwide wind curvature. However, many observations and high-resolution simulations contradict this simple picture, suggesting a more complex one with a large influence of topography and eddies (Spall 1996; Chassignet and Marshall 2008; Özgökmen and Chassignet 2002; Bryan et al. 1995, 2007; Couvelard et al. 2008; Gula et al. 2015, 2016; Renault et al. 2016b, 2019; Debreu et al. 2022). The limitations of low-resolution (>100 km) ocean models in reproducing Gulf Stream dynamics have been clearly demonstrated and significant biases have been observed near the separation of the Gulf Stream as well as in its northwestward penetration (Chassignet and Xu 2017; Bryan et al. 2007). On the other hand, if the refinement of the computational grid to reach an eddy-rich resolution largely improves the Gulf Stream representation and in particular its separation, a truly realistic representation seems only possible with adequate energy sinks.

At large scale, balanced motions dominate the ocean. Mesoscale eddies, generated by baroclinic and barotropic instabilities of the mean Gulf Stream (drawing energy from its potential energy), are also dominated by geostrophic balance. The quasigeostrophic theory then predicts that their kinetic

energy is essentially transferred to larger scales (i.e., an inverse cascade; Charney 1971; Arbic et al. 2013; Renault et al. 2019). To maintain the energy balance, the energy input in the Gulf Stream system must be dissipated. Some of the energy is dissipated at the surface (Renault et al. 2016b, 2019) and bottom boundaries (top and bottom drag, respectively), but another, more uncertain, portion may follow an interior route to dissipation when the momentum balance is broken—i.e., a direct transfer to smaller scales or forward cascade. For the models to produce a forward cascade, quasi-geostrophy must be relaxed (Molemaker et al. 2010), which is most likely to occur in the detached Gulf Stream jet (Jamet et al. 2021), allowing ageostrophic motions to reverse the spectral energy fluxes (e.g., Brüggemann and Eden 2015; Klein et al. 2008; Capet et al. 2008a; Molemaker et al. 2010). However, the persistence (in space and time) and effectiveness of this dissipation route is unclear.

Submesoscale oceanic fronts and eddies are mostly unbalanced motions and carry a substantial ageostrophic component. They are a relatively recent discovery, and their interest as a research topic has grown rapidly (McWilliams 2016). They occur on an intermediate scale of the order of 0.1–10 km horizontally, 10–100 m vertically, and from hours to days temporally, that is, smaller and more rapidly evolving than the mesoscale eddies and with vertical velocities that are one or two orders of magnitude larger (Capet et al. 2008b; Su et al. 2020; Siegelman 2020; McWilliams 2021). These can affect momentum, buoyancy, nutrient transport and biogeochemistry (Lévy et al. 2018; Uchida et al. 2019; Kessouri et al. 2020), and gas exchange between the ocean and atmosphere (Su et al. 2018). Important here, they may also provide an oceanic interior energy route to dissipation (Gula et al. 2016; Schubert et al. 2020), which may partially justify the use of a high turbulent viscosity to control Gulf Stream dynamics in eddy-rich models

<footnote>*Corresponding author:* Marcela Contreras, marcela.contreras@legos.obs-mip.fr</footnote>
(Chassignet et al. 2003; Chassignet and Marshall 2008; Chassignet and Xu 2017; Schoonover et al. 2016). The principle invoked is that, at submesoscales, the effect of rotation and stratification being weaker, the geostrophic or gradient-wind balance is disrupted, allowing the divergent flow component to transfer energy to smaller scales (Capet et al. 2008b; Klein et al. 2008). The proposed ageostrophic mechanisms are inertia-gravity wave emission by eddies (Barkan et al. 2015); frontogenesis (Capet et al. 2008c); submesoscale quasigeostrophic instabilities at finite Rossby number (Capet et al. 2016); and at even smaller scales, ageostrophic instabilities such as centrifugal, symmetric, gravitational, or Kelvin-Helmholtz instabilities (McWilliams 2016; Thomas et al. 2013; Haine and Marshall 1998).

However, submesoscale processes can also reinforce mesoscale eddies, extending the inverse kinetic energy (KE) cascade to a scale of a few tens of kilometers (Klein et al. 2019). This extension is particularly evidenced for frontogenesis (Klein et al. 2019) and for submesoscale quasigeostrophic instabilities that draw energy from the mesoscale potential energy at spatial scales of $O(10)$ km and time scales of $O(1)$ day$^{-1}$, for example, mixed layer instabilities (MLIs; Schubert et al. 2020; Boccaletti et al. 2007) and Charney baroclinic instabilities (Capet et al. 2016). The broadening of the inverse cascade suggests that the KE of the ocean is less controlled by viscosity than one might think by considering only the effect of submesoscale processes on the forward cascade. In any case, the overall effect of submesoscales on the oceanic turbulent cascade is complex and may vary in space and time.

This study presents a submesoscale simulation at 2-km resolution, where strain-induced frontogenesis and mixed layer instabilities are resolved, but ageostrophic instabilities and associated secondary frontogenesis are not. The simulation is carried out over the Gulf Stream for a period of 5 years with three objectives. First, we assess the temporal and spatial distribution of cross-scale KE fluxes in the Gulf Stream. The evaluation is enabled by the use of a coarse-grained approach (Alue et al. 2018). Second, we explore the submesoscale mechanisms and the contributions of rotational and divergent velocity components (associated with the balanced and unbalanced motions, respectively), in particular the extent to which they affect energy fluxes over the Gulf Stream. Last, we quantitatively evaluate the different energy dissipation processes. We compare the magnitude of the forward cascade with the dissipation by the top and bottom drag and by numerical discretization effects. Next, we analyze the sensitivity of horizontal momentum advection schemes on energy dissipation. The paper is organized as follows: section 2 describes the model configuration and method. In section 3, the spatial and season variability of cross-scale KE fluxes is assessed. Section 4 examines the possible mechanisms involved in the energy transfer. The quantification of energy pathways is discussed in section 5, followed by the conclusions.

## 2. Method

### a. Model configuration

The oceanic simulation is performed with the Coastal and Regional Ocean Community model (CROCO; Shchepetkin

and McWilliams 2005; Debru et al. 2012). CROCO is a free-surface, terrain-following coordinate model with split-explicit time stepping. In the present study, the version with Boussinesq and hydrostatic approximations is used, and the equations are discretized with high-order computational methods. A third-order predictor-corrector time step algorithm, and fifth-order upstream biased momentum advection (instead of the more standard third-order scheme) allows reducing numerical dispersion and diffusion to achieve an effective resolution of about 5 times the horizontal resolution (Soufflet et al. 2016; Ménesguen et al. 2018). The advection of horizontal tracers is discretized using the rotated split third-order upstream scheme (Marchesiello et al. 2009; Lemarié et al. 2012). A nonlocal planetary K-profile boundary layer (KPP) scheme (Large et al. 1994) parameterizes the vertical subgrid-scale eddy effects at the surface, bottom, and interior of the ocean.

The simulation domain extends from 22.5° to 48.84°N and from 36° to 82°W (Fig. 1) with a horizontal resolution of 1/42° (~2.2 km). The topography is obtained from SRTM15-plus (Shuttle Radar Topography Mission), which is a global bathymetric dataset with a nominal resolution of 15 arc s (~0.5 km) (http://topex.ucsd.edu/WWW_html/srtm30_plus.html). To reduce errors associated with the pressure gradient in a terrain-following coordinate model—in addition to the high-order correction scheme of Shchepetkin and McWilliams (2003)—the bottom topography is interpolated onto the model grid and its slope is selectively smoothed when the r-factor $\delta h/h$ exceeds 0.2, where $h$ is the bottom depth and $\delta h$ is its horizontal variation between two grid points (Penven et al. 2005).

The simulation is performed over the period 2005–09 after a 5-yr spinup. The initial field and lateral boundary forcing are derived from the 1/12° daily mean Mercator Glorys12V1 product (Lellouche et al. 2018). The boundary forcing is introduced through open boundary conditions that consist of an active-passive 2D radiation scheme for the baroclinic mode (including temperature $T$ and salinity $S$) and a modified Flather-type scheme for the barotropic mode (Marchesiello et al. 2001). We use 50 $\sigma$ levels in the vertical direction (Shchepetkin and McWilliams 2009), with stretching parameters $h_{\text{cline}} = 200$ m, $\theta_b = 2$, and $\theta_s = 7$.

Surface momentum, heat, and freshwater fluxes are estimated using the COARE bulk formula (Fairall et al. 2003). Surface atmospheric fields are derived from the hourly Climate Forecast System Reanalysis (CFSR; Saha et al. 2010), which has a spatial resolution of ≈35 km. The surface ocean current feedback to the atmosphere (that controls the top drag) (Renault et al. 2016b) is parameterized using the stress-correction approach described in Renault et al. (2020), allowing realistic representation of the (sub)mesoscale momentum exchange between the ocean and the atmosphere. In our analysis, 3-h averaged output fields are used.

The simulation is generally in good agreement with the observations. As an example, Fig. 2 compares the simulated mean and mesoscale circulations with those observed by the AVISO product (Ducet et al. 2000) and drifters (Laurindo et al. 2017). To ensure a fair model-data comparison when comparing with AVISO, we first apply a 50-km Gaussian spatial filter and a 7-day average on the simulated Absolute
<figure><img src="image_1.png" /><figcaption>FIG. 1. Snapshots of (a) surface relative vorticity (normalized by f) and (b) sea surface temperature in March 2007 from CROCO.</figcaption></figure>

Dynamic Topography (ADT; Arbic et al. 2013). Figures 2a and 2b show the mean ADT as estimated from CROCO and AVISO as well as the mean Gulf Stream path highlighted using the contour of 0.5 m s⁻¹ of the mean geostrophic currents. The mean circulation is well reproduced by the model and in particular the Gulf Stream path shows very good agreement with the observations.

The surface geostrophic eddy kinetic energy [EKE($u_g$)] is computed as a proxy of the mesoscale activity from daily geostrophic surface current perturbations in CROCO and AVISO (Figs. 2b,c). The perturbations are estimated as the deviation of the 3-month average. The Gulf Stream path stability is associated with the EKE: the larger the EKE, the more unstable the Gulf Stream trajectory. In agreement

<figure><img src="image_2.png" /><figcaption>FIG. 2. (a),(b) Mean dynamic topography and (c),(d) EKE estimated using geostrophic currents from (top) CROCO and (bottom) AVISO. The magenta and green contours represent the mean Gulf Stream path (0.5 m s⁻¹ contour of the mean geostrophic current) from AVISO and CROCO, respectively. Also shown is the time average of EKE estimated using the surface total currents from (e) CROCO and (f) drifter data. The magenta and green contours again represent the mean Gulf Stream path but from the mean total currents (0.5 m s⁻¹ contour) from drifters and CROCO, respectively.</figcaption></figure>
with previous studies, CROCO and observations show a stable trajectory separating straight from Cape Hatteras, with a larger EKE magnitude over the Gulf Stream core after separation. The simulated EKE($u_g$) shows very close spatial patterns and amplitudes to those observed by AVISO (Figs. 2c,d). Similar results are found when comparing the model EKE [EKE($u$)] estimated from the total surface currents with drifter observations from the Global Drifter Program (Laurindo et al. 2017) (Figs. 2e,f)—note that the drifters are also included in the AVISO estimate of the mean dynamic topography.

## b. Coarse-grained method

In the literature, the cross-scale KE flux in the ocean is generally evaluated using a spectral analysis (e.g., Scott and Wang 2005; Arbic et al. 2013; Marchesiello et al. 2011; Renault et al. 2019). Based on this approach, several studies demonstrate the existence of direct or inverse energy cascades. However, it has several drawbacks. First, it assumes a statistically homogeneous field and a value for the small-scale flux (Schubert et al. 2020). Second, the spectral approach underestimates the energy flux due to the windowing required, and the result is very sensitive to the resolved scales of the model or altimeter data and the filtering methods used on these data (Arbic et al. 2013).

The coarse-grained method (Leonard 1975; Germano 1992) is an interesting alternative to assess the cross-scale KE fluxes in the ocean. Aluie et al. (2018) are the first to apply this method to ocean dynamics, followed by Schubert et al. (2020) and Srinivasan et al. (2019) to specifically address the interior route to dissipation. The coarse-grained approach has several advantages over spectral analysis, including relaxing the homogeneous and isotropic field assumption, avoiding windowing procedures, and, most importantly, having access to spatial patterns of energy fluxes.

The coarse-grained approach is based on a “filtering” framework, which separates the signal into large and small scales around a given scale $L$. Following Aluie et al. (2018) and Schubert et al. (2020), a low-pass spatial filter is defined using the convolution of a horizontal field $\bar{F}(x, y)$ as $F(x, y) = C \times F(x, y)$, where $C$ is a top-hat kernel defined as

$$C(r) = \begin{cases} 1/A & \text{if } |r| < L/2 \\ 0 & \text{otherwise} \end{cases},$$

$A = \pi L^2/4$ is the circular normalization area of diameter $L$, and $r$ is the radial vector.

By applying this convolution to the equation of motion, we can estimate the term that represents the scale transfer of kinetic energy:

$$\Pi = - \rho_0 [(\bar{u}^2 - \bar{u}^2)\bar{u}_x + (\bar{u}\bar{v} - \bar{u}\bar{v})\bar{u}_y + \bar{v}_x] + (\bar{v}^2 - \bar{v}^2)\bar{v}_y, \quad (1)$$

where $\rho_0$ is the reference density and $\mathbf{u} = (u, v, w)$ are the oceanic components of the surface current in the zonal, meridional, and vertical direction [see Aluie et al. (2018) for the full derivation]. This term quantifies the energy transferred from

the $L$ scale to larger or smaller scales (i.e., a cross-scale flux at $L$), due to nonlinear interactions.

In this study, the coarse-grained approach is applied to the total surface current as well as its rotational ($\mathbf{u}_r$) and divergent ($\mathbf{u}_d$) components using the Helmholtz decomposition. The latter consists of a separation of the total currents in terms of the streamfunction $\psi$ and potential velocity $\phi$:

$$u = u_r + u_d = \psi_y + \phi_x \quad \text{and} \quad v = v_r + v_d = -\psi_x + \phi_y.$$

We obtain $\mathbf{u}_d$ by solving $\phi$ of a Poisson equation. We then estimate $\mathbf{u}_r$ as the difference between $\mathbf{u}$ and $\mathbf{u}_d$. Components $\mathbf{u}_r$ and $\mathbf{u}_d$ are associated with the balanced and unbalanced motions, which allows us to separate the contributions of these components to the energy cascade.

# 3. Description of the turbulent cascade

## a. Cross-scale KE fluxes

The kinetic energy flux $\Pi$ is estimated over the Gulf Stream region using the coarse-grained approach at scales of 9, 22, 61, and 105 km (Fig. 3). The 9-km scale is the closest wavelength to the effective resolution of the simulation, 22 km corresponds approximately to a transition scale between mesoscale and submesoscale currents in this region, and 61 and 105 km represent scales of high mesoscale activity. Positive (negative) values indicate a forward (inverse) cascade, i.e., a transfer of kinetic energy to smaller (larger) scales.

At the 9-km scale (Fig. 3a), $\Pi$ is mostly positive, revealing the presence of a systematic forward cascade at this scale. In contrast, at the 22-km scale, $\Pi$ is characterized by a dipole located right on the Gulf Stream path (Fig. 3b), with a stronger downscale flux on its northern edge and a weaker upscale flux on its southern edge. At larger scales (61 and 105 km), the inverse cascade becomes dominant (Figs. 3c,d), which is consistent with an intense mesoscale eddy—mean flow interaction (Arbic et al. 2013; Renault et al. 2019). However, at these scales, an intense narrow band of positive values remains near Cape Hatteras, where the Gulf Stream separates from the coast.

On all scales analyzed, $\Pi$ is most intense over the Gulf Stream and its surroundings, i.e., where the eddy activity is important (Fig. 1a). However, the various mechanisms that drive the kinetic energy flux—barotropic and baroclinic instabilities, filamentation and frontogenesis, mixed layer and ageostrophic instabilities—may bring a spatial dependency. Of particular importance is the Gulf Stream North Wall where significant density and temperature gradients (Fig. 1b) would allow for an intense frontogenesis activity (McWilliams et al. 2019). In the South Atlantic Bight, the topography plays an additional role, favoring barotropic destabilization of the flow (Gula et al. 2015) and there is a prominent topographic feature called the Charleston Bump where mesoscale frontal eddies are generated (Olson et al. 1983; Gula et al. 2015). Over the Charleston Bump, Fig. 3 reveals a large forward cascade, not only at small scales (<22 km) but also at large scales (>61 km), even though the inverse cascade intensifies in other regions. The next section will focus on the different components of the
<figure><img src="image_3.png" /><figcaption>FIG. 3. Time-averaged (2005–09) cross-scale surface kinetic energy flux estimated from (left) total currents (Π), (center) rotational (balanced) currents [Π(**u**<sub>r</sub>)], and (right) the interaction of balanced and unbalanced flow [Π<sub>CT</sub> = Π − Π(**u**<sub>r</sub>) − Π(**u**<sub>d</sub>)]. The cross-scale fluxes are estimated at 9, 22, 61, and 105 km.</figcaption></figure>

flow, particularly the role of unbalanced motions in the down-
scale KE flux.

b. Rotational and divergent components

To understand the role of balanced and unbalanced mo-
tions in energy transfer at different scales, we decompose
the surface currents into rotational and divergent components
(i.e., Helmholtz decomposition, see section 2), where the rota-
tional (divergent) velocity component is associated with the
balanced (unbalanced) motion. Note that geostrophic currents

are nondivergent and included in the rotational component.
Using these components, we estimated Π from (1) using rota-
tional and divergent components.

Figures 3e-h depict the kinetic energy flux estimated using the rotational velocity component [Π(**u**<sub>r</sub>)]. It reveals that, at all scales (stronger with increasing scale) and almost everywhere, the balanced kinetic energy flux is dominated by up-scale fluxes (Figs. 3e-h)—there are some exceptions where Π(**u**<sub>r</sub>) is positive for a few very localized regions such as near Cape Hatteras and over the Charleston Bump. These results
confirm the important role of balanced motions in the inverse cascade of kinetic energy. It is known that when balanced motions dominate, particularly the geostrophy that prevails at scales larger than 10 km, the effect of rotation and stratification is to inhibit vertical velocities and promote quasi-two-dimensional turbulence. However, even at small scales, the nondivergent part of the flow can hardly drive a forward cascade.

The KE flux $\Pi(\mathbf{u}_d)$ estimated using the divergent velocity component alone (physically representing the advection of unbalanced motion by itself) is an order of magnitude smaller than both $\Pi$ and $\Pi(\mathbf{u}_r)$ (not shown). Therefore, the total cross-scale KE flux is not explained by the simple addition of rotational $\Pi(\mathbf{u}_r)$ and divergent $\Pi(\mathbf{u}_d)$ components, but by an additional cross term representing the interaction of balanced and unbalanced motions. The cross-term contribution to the kinetic energy flux of (1) is thus estimated as $\Pi_{CT} = \Pi - \Pi(\mathbf{u}_r) - \Pi(\mathbf{u}_d)$. In Figs. 3i-l, $\Pi_{CT}$ is generally positive at all scales and the most intense values are located over the Gulf Stream. Note that patches of negative $\Pi_{CT}$ are apparent, particularly along the southern edge of the Gulf Stream at large scale (>61 km; Figs. 3k,l).

To further quantify the contribution of $\Pi(\mathbf{u}_r)$, $\Pi(\mathbf{u}_d)$, and $\Pi_{CT}$ to the total KE flux, we calculate the percent contribution of each component as $\%\Pi_X = |\Pi_X|/|[\Pi(\mathbf{u}_r)] + [\Pi(\mathbf{u}_d)] + |\Pi_{CT}||\Pi(\mathbf{u}_d)||\Pi(\mathbf{u}_r)|$, where $\Pi_X$ can be $\Pi(\mathbf{u}_r)$, $\Pi(\mathbf{u}_d)$, and $\Pi_{CT}$. Averaging $\%\Pi_X$ over the entire domain shows that at scales greater than or equal to 22 km $\Pi(\mathbf{u}_r)$ explains more than 50% of $\Pi$ and $\Pi_{CT}$ explains about 39%. At 9 km, the contribution of $\Pi_{CT}$ increases to 68% and $\Pi(\mathbf{u}_r)$ reduces to 26%. The contribution of $\Pi(\mathbf{u}_d)$ to the total cross-scale KE flux is never larger than 6% at 9 km and falls below 4% at larger scales. Note that over areas of strong eddy activity in the Gulf Stream, the contribution of $\Pi(\mathbf{u}_r)$ can be much higher than in the domain average (close to 90% at 105-km scales), whereas $\Pi_{CT}$ has a stronger contribution over regions of weak eddy activity.

Previous studies have examined the origin of the forward cascade (Brüggemann and Eden 2015; Capet et al. 2008a; Molemaker et al. 2010). In particular, Capet et al. (2008a) note that in surface quasigeostrophic models that do not support momentum advection by ageostrophic currents, the forward cascade is weaker than in primitive equation models. They suggest that the forward cascade is associated with this missing part of the advection. Here, to confirm this hypothesis and extend it to balanced and unbalanced motions, $\Pi_{CT}$ is further decomposed into the advection of rotational by divergent currents ($\Pi_{CT_d}$) and the advection of divergent by rotational currents ($\Pi_{CT_r}$), as follows:

$$ \Pi_{CT_d} = X[(\mathbf{u}_d), (\mathbf{u}_d, \mathbf{u}_r)] \cdot (\bar{\mathbf{u}}_{FX}, \bar{\mathbf{u}}_{RY}) \quad \text{and} \quad (2) $$

$$ \Pi_{CT_r} = Y[(\mathbf{u}_r), (\mathbf{u}_d, \mathbf{u}_r)] \cdot (\bar{\mathbf{u}}_{dx}, \bar{\mathbf{u}}_{dy}), \quad (3) $$

where X and Y are nonlinear operators, functions of $\mathbf{u}_d$ or $\mathbf{u}_r$ and their combinations. Note that $\Pi_{CT}$ can be defined explicitly as $\Pi_{CT} = \Pi_{CT_r} + \Pi_{CT_d}$ or implicitly as $\Pi_{CT} = \Pi - \Pi(\mathbf{u}_d) + \Pi(\mathbf{u}_r)$, both yielding the same result. The decomposition is applied at scales of 9, 22, and 61 km (Fig. 4; similar results are found at

105 km). Figure 4 reveals that $\Pi_{CT}$ at 9, 22, and 61 km is primarily explained by $\Pi_{CT_d}$ (Figs. 4a-c), i.e., by the advection of balanced momentum by unbalanced currents (unbalanced advection). The $\Pi_{CT_d}$ shows dominant downscale fluxes at all scales, especially over the Gulf Stream North Wall, a region of strong density gradients. At large scale, upscale fluxes appear on the southern edge, while weaker positive and negative values are also present in $\Pi_{CT_r}$ (Figs. 4d-f). Averaging over the entire region, $\Pi_{CT_d}$ explains a large majority of $\Pi_{CT}$ at all scales. Similar results were found with geostrophic and ageostrophic decompositions (not shown). Note that the advection of unbalanced by balanced currents may also play some role in specific regions such as the Charleston Bump. In summary, unbalanced advection has almost exclusive control over the forward cascade, while balanced motion drives the inverse cascade.

### c. Seasonal variability

A significant seasonal variability of the cross-scale KE flux has been estimated, based on high-resolution simulations, over the Agulhas Current (Schubert et al. 2020) and the subtropical western North Pacific (Qiu et al. 2014). Here, we estimate $\Pi$ for the Gulf Stream system (Fig. 5) in winter (January–March) and summer (July–September). The fluxes are generally more intense in winter, and the patterns are similar to the annual mean during this season. This is particularly true at the 9- and 22-km scale (Figs. 5a-d), where the summer fluxes are much reduced. At 61 and 105 km, in winter, the narrow band of downscale fluxes stretching from Cape Hatteras is more evident than in the annual mean and actually disappears in summer.

The inverse cascade given by $\Pi$ has a lesser seasonal cycle than appears from $\Pi(\mathbf{u}_r)$ (Fig. 6). This difference is explained by a compensation between $\Pi(\mathbf{u}_r)$ and $\Pi_{CT}$, which are respectively responsible for upscale and downscale fluxes. In winter, the inverse cascade driven by $\Pi(\mathbf{u}_r)$ is intensified (Figs. 6a-d) but the forward cascade due to $\Pi_{CT}$ is also intensified (Figs. 6i-l), balancing part of the upscale fluxes. The opposite is true in summer [weakening of both $\Pi(\mathbf{u}_r)$ in Figs. 6e-h and $\Pi_{CT}$ in Figs. 6m-p] with a similar moderate effect on upscale fluxes.

Previous studies suggest that seasonal variability not only alters the intensity of KE fluxes, but may also affect other indicators of seasonality. For example, Schubert et al. (2020) observe that the scale at which the peak upscale fluxes occur varies with the season, as does the scale at which the change of sign of the KE flux occurs (cross-over scale between downscale and upscale fluxes). To better analyze the temporal and spatial variability of $\Pi$, we focus on two regions, north and south of the Gulf Stream (Fig. 7; the areas of integration are drawn in Fig. 5). These two regions are representative of differences observed on the maps, focusing on them allows assessing the evolution of $\Pi$ as a function of wavelength (as in a wavelength spectrum). Figure 7 reveals that $\Pi$ has lower values in the southern region than in the northern region. In the northern region, the cross-over scale is also higher than in the southern region (~30–50 and ~10–30 km, respectively).
<figure><img src="image_4.png" /><figcaption>FIG. 4. Time-averaged (2005–09) cross-term contribution to the cross-scale KE flux, decomposed into terms associated with the (a)–(c) advection by unbalanced flow ($\Pi_{CT_d}$) and (d)–(f) advection by balanced flow ($\Pi_{CT_r}$). The cross-scale flux is estimated at (top) 9, (middle) 22, and (bottom) 61 km.</figcaption></figure>

The cross-over scale shows seasonal variability: in the north- ern region, it is shorter in spring and longer in winter; in contrast, in the southern region, it is shorter in fall–winter and longer in summer.

In both regions, the forward cascade peaks in winter, while the inverse cascade has its larger values in spring. However, in the northern region, the inverse cascade is similar in winter, summer, and fall. In contrast, in the southern region, at scales
<figure><img src="image_5.png" /><figcaption>FIG. 5. (top) Winter [January–March (JFM)] and (bottom) summer [July–September (JAS)] averaged cross-scale kinetic energy fluxes estimated from total surface currents at (a), (b) 9; (c), (d) 22; (e), (f) 61; and (g), (h) 105 km. Blue outlines indicate the northern and southern regions analyzed in Fig. 7, below.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>FIG. 6. Winter (JFM) mean of (a)–(d) $\Pi(\mathbf{u}_r)$ and (i)–(l) $\Pi_{CT}$, and summer (JAS) mean of (e)–(h) $\Pi(\mathbf{u}_r)$ and (m)–(p) $\Pi_{CT}$. Cross-scale energy fluxes are estimated at (left) 9, (left center) 22, (right center) 61, and (right) 105 km.</figcaption></figure>

larger than 62 km, the upscale flux is weaker in autumn, and in summer it shows a more marked decrease toward smaller scales. Interestingly, in the northern region, the inverse cascade peaks at ≈150 km and does not have significant seasonal variation. In the southern region, the peak scale is larger in summer (101-km scale) and smaller in winter and autumn (62-km scale). These results suggest that the mechanisms involved in the kinetic energy flux are different in the two regions (this will be further discussed in section 4).

## 4. Mechanisms

The objective of this section is to evaluate the mechanisms that determine the temporal and spatial variability of $\Pi$, as a measure of the KE cascade. We will be particularly interested in the mechanisms of forward cascade.

### a. Frontogenesis

As shown in section 3, the forward cascade is produced primarily by the interaction between balanced and unbalanced motions ($\Pi_{CT}$), and in particular by the advection of momentum by unbalanced flows ($\Pi_{CT_d}$).

Frontogenesis occurs in regions where there is a horizontal buoyancy gradient in a background horizontal deformation flow. The deformation intensifies the horizontal buoyancy

gradient and disrupts the geostrophic balance along the front. To restore balance, a secondary ageostrophic circulation appears, which acts by restratifying the subsurface. The secondary circulation mechanism involves the generation of cyclonic vorticity and downwelling in the cold region and anticyclonic vorticity and upwelling of weaker amplitude in the warm region (McWilliams 2016). Based on submesoscale permitting simulations, Capet et al. (2008a,d); Klein et al. (2008) suggest that the forward cascade is essentially associated with this frontogenetic process.

The Gulf Stream presents favorable conditions for frontogenesis due to intense temperature gradients, especially at the North Wall and strong deformation flow (McWilliams et al. 2019). Figures 8a-c represent the annual, winter, and summer mean relative vorticity (i.e., normalized by the Coriolis parameter $f$). From the Charleston Bump to the Gulf Stream postseparation, the vorticity is characterized by the presence of a dipole. On the northern (southern) edge of the Gulf Stream, the water masses are denser (lighter) and a cyclonic (anticyclonic) band is present, which is consistent with the observed frontal structure. Interestingly, the cyclonic side is spatially collocated with the 60 mW m⁻² km⁻¹ contour of $\Pi_{CT_d}$ at the 22-km scale (the contour shows the maximum values of $\Pi_{CT_d}$). Note that similar consistency is found using contours of the same magnitude for $\Pi_{CT_d}$ at the scales of 9 and 61 km.
<figure><img src="image_7.png" /><figcaption>FIG. 7. Spatially averaged surface kinetic energy fluxes in the (a) northern and (b) southern regions, defined by the blue outlines in Fig. 5.</figcaption></figure>

Following Hoskins and Bretherton (1972), the frontogenesis function is defined as $F_s = Q \cdot \nabla_h \rho$, where Q is the Q vector. This function quantifies the rate of increase in the horizontal buoyancy gradient at the surface. Positive values indicate the development of frontogenesis and negative values indicate the development of frontolysis, i.e., a decrease of the buoyancy gradient. Figures 8d–f show that frontogenesis dominates in the region, especially at the North Wall of the Gulf Stream. $F_s$ also presents a strong spatial and temporal coherence with the 60 mW m$^{-2}$ km$^{-1}$ contour of $\Pi_{CT_d}$ at the 22-km scale.

These results are in agreement with previous studies. In particular, Capet et al. (2008a) associate the ageostrophic advection responsible for the forward cascade with the secondary circulation generated by frontogenesis. Here, we confirm and extend their proposition to unbalanced (nondivergent) rather than ageostrophic motions. Recently, Srinivasan et al. (2022) use submesoscale-resolving simulations of a North Atlantic region between Greenland and Iceland and compute as we do the cross-scale energy flux using a Helmholtz decomposition. They find similar results to ours and conclude, with the help of an asymptotic theory, that the primary mechanism for the forward energy flux at fronts is frontogenesis.

While frontogenesis is an effective means of transferring energy to smaller scales, other processes can be considered, such as ageostrophic frontal instabilities (Molemaker et al. 2010;

Gula et al. 2016; Thomas et al. 2013). Frontal instabilities occur when chaotic advection by mesoscale eddies, intensified by submesoscale strain (due, for example, to MLI), drives the process to frontal collapse (Callies et al. 2016; Hoskins and Bretherton 1972). Because our simulations lack the spatial resolution to adequately reproduce these instabilities, we cannot analyze their role in the forward cascade. We expect their contribution to be significant, but perhaps not by increasing the energy flux as much as by taking over, at finer scales, the processes at work in finite Rossby number quasigeostrophic dynamics (McWilliams 2016). This hypothesis is also suggested by Marchesiello et al. (2011), whose model solutions from submesoscale-permitting to submesoscale-resolving show a continuation of the forward cascade but no increase in flux intensity.

North of the Gulf Stream, the relation between $F_s$ and $\Pi_{CT_d}$ appears to be lost: $F_s$ has large values (especially in summer; Fig. 8f) while $\Pi_{CT_d}$ is weak (Fig. 4b). One possible explanation is that in this region, the mesoscale EKE is weak relative to the Gulf Stream region and the frontogenesis function may reflect other processes than mesoscale confluent flow deformation.

b. Ekman buoyancy flux

Downstream of the Charleston Bump along the Gulf Stream, both the forward cascade ($\Pi$ positive) and frontogenesis ($F_s$) show strong seasonal variability. However, it is unclear why frontogenesis intensifies so strongly in winter, which raises the question: what processes explain the seasonal variability of frontogenesis? Atmospheric forcing and, in particular, wind may be good candidates to explain such a variability. Thomas and Lee (2005) suggest that, depending on wind direction, wind-driven Ekman currents can alter the ageostrophic secondary circulation. In particular, when the winds blow in the direction of the frontal jet (downfront winds), they can lead to an intensification of the front. To assess this hypothesis, following Thomas and Taylor (2010), the Ekman buoyancy flux (EBF) is defined as $\text{EBF} = M_e \cdot \nabla_h b_{z=0}$, where $M_e = \tau \times \mathbf{z}/(\rho_0 f)$ is the Ekman transport and $\tau$ is the wind stress vector. A positive EBF indicates an intensification of frontogenesis by Ekman transport.

Figures 8g-i reveal the presence of a dipolar structure in EBF downstream of the Charleston Bump, where positive and negative values are observed on the northern and southern edges, respectively, of the Gulf Stream. This structure correlates well with that of $\Pi_{CT_d}$ and $F_s$, not only in space but also in time. In particular, EBF is positive over the Gulf Stream North Wall in the area covered by the 60 mW m$^{-2}$ km$^{-1}$ contour of $\Pi_{CT_d}$ at the 22-km scale (Fig. 8g). Therefore, wind-driven mixing contributes here to the transfer of energy to smaller scales. On the other hand, EBF is negative on the southern edge of the Gulf Stream, where $\Pi_{CT_d}$ is small or even weakly negative (Figs. 4a-c). These trends are stronger in winter, when the wind and associated EBF are more intense than in summer, when the dipole disappears (Figs. 8h-i).

The Ekman buoyancy flux may also have a role in driving instabilities that promote or continue the forward cascade initiated by frontogenesis. Thomas et al. (2013) show for the
<figure><img src="image_8.png" /><figcaption>FIG. 8. Time-averaged (2005–09) surface (a)–(c) relative vorticity (normalized by f), (d)–(f) frontogenesis function ($F_s$: kg$^2$m$^{-8}$s$^{-1}$), and (g)–(i) Ekman buoyancy flux (EBF). The black contour represents the 60 mW m$^{-2}$km$^{-1}$ value of $\Pi_{CTd}$ at 22 km.</figcaption></figure>

Gulf Stream that a cross-front advection by the Ekman flow of negative potential vorticity favors the development of symmetric instabilities that would drive energy farther toward the dissipation scale. D'Asaro et al. (2011) had drawn similar conclusions for the Kuroshio current. Again, due to the resolution of our model, we cannot verify these results.

### c. Upscaling effect of submesoscale eddies

In geostrophic turbulence theory, an inverse cascade is produced by the merging of mesoscale coherent eddies starting from the scales of eddy sources (Vallis 2017). However, sub-mesoscale-permitting models reveal that the inverse cascade is intensified by submesoscale processes and extended to scales down to tenths of a kilometer (Klein et al. 2019). MLIs are a particularly effective example, and the intensification of mixed layer eddies in winter due to mixed layer deepening can be absorbed by mesoscale eddies, strengthening the latter and producing larger upscale fluxes (Khatri et al. 2021; Qiu et al. 2014; Schubert et al. 2020).

To investigate the submesoscale effect on the mesoscale activity, we use the seasonal analyses presented in section 3c. In the southern region, the seasonality of $\Pi$ in Fig. 7b suggests two possible mechanisms related to MLIs and their role in the inverse cascade. First, consistent with Schubert et al. (2020), the maximum upscale flux occurs at smaller scales in winter as

compared with summer (~50 vs ~100 km, respectively). Second, the cross-over scale between downscale and upscale fluxes is ~13 km, which is similar to the most unstable MLI wavelength (15 km; Fox-Kemper et al. 2008). The cross-over scale of KE fluxes increases in summer to about 30 km, a seasonality perhaps weaker than reported by Schubert et al. (2020) for the Agulhas system. In the northern region, the upscale flux is more intense and its seasonality shows a significantly stronger flux in spring and cross-over scale varying from about 30 km in spring-summer to about 50 km in fall-winter.

These results agree with Schubert et al. (2020) and their interpretation that submesoscale flows affect the mesoscale seasonal cycle. The peak of MLI activity occurs in winter when the mixed layer is deeper, but the effect on mesoscale eddies (an absorption effect according to Schubert et al. 2020) develops during this season to peak in spring. We thus confirm for the Gulf Stream a phase shift of a few months in the effect of submesoscales on mesoscale activity.

### 5. Quantification of energy pathways

Mesoscale-resolving simulations were associated in the past with excessive mesoscale activity, reflecting a missing energy sink (e.g., Renault et al. 2019). Two mechanisms that may
have been underestimated can explain this: (i) interior dissipation driven by a forward energy cascade and (ii) boundary layer dissipation produced by top and bottom drag. In addition to these mechanisms, numerical models impose another energy sink, which is the dissipation required to dampen numerical dispersion errors. Since no explicit momentum dissipation operator is used in our simulation, numerical dissipation is essentially handled by the diffusive part of an upstream momentum advection scheme, whose viscosity varies with velocity powers (Soufflet et al. 2016). In what follows, the respective roles of interior, boundary, and numerical dissipation are evaluated, followed by a sensitivity analysis to advection schemes.

### a. Interior, boundary, and numerical dissipation

Dissipation of mesoscale currents by top drag can be estimated using the eddy geostrophic wind work $F_e K_e = \langle \mathbf{u}'_g \tau' \rangle - \Pi_{\tau_{22km}}$ where angle-bracket and prime operators indicate the average over 3 months and its fluctuation, respectively, and $\Pi_{\tau_{22km}} = \bar{\tau} \cdot \bar{u}_g - \bar{\tau} \cdot \bar{u}_g$ is the coarse-grained wind work at 22-km scale. The $F_e K_e$ is calculated with geostrophic currents, excluding the Ekman or any ageostrophic contribution. Subtracting $\Pi_{22km}$ from $F_e K_e$ allows us to remove also the contribution of submesoscale currents and thus to further isolate the energy transfer from the mesoscale field to the atmosphere. Note, however, that $\Pi_{\tau_{22km}}$ is only about 10% of $F_e K_e$. Dissipation by bottom drag is computed from the covariance of the bottom currents and bottom stress anomalies $F_b K_b = \langle \mathbf{u}'_b \tau'_b \rangle$. Last, the numerical dissipation is $H_{\text{Diff}} = \int_{\text{bottom}}^{\text{surface}} \mathbf{u} \cdot D \, dz$, where the horizontal momentum diffusion $D$ is diagnosed as either $A_{\text{UP5}} - A_{\text{C6}}$, that is, the difference between the UP5 advection term and an analogous advection term computed with a centered sixth-order nondissipative scheme, or $A_{\text{UP3}} - A_{\text{C4}}$, that is, the difference between the UP3 advection term and an analogous advection term computed with a centered fourth-order nondissipative scheme. Negative values of $F_e K_e$, $F_b K_b$, and $H_{\text{Diff}}$ indicate a dissipation of energy. These terms are directly comparable to the downscale flux at 9 km integrated between 100-m depth and the surface: $I_{\text{Diss}} = -\int_{-100\text{m}}^{\text{surface}} \Pi_{9\text{km}} dz$, that is, where the forward cascade dominates. Note that $I_{\text{Diss}}$ is defined with a negative sign for downscale fluxes to represent a dissipation of energy. We consider 9 km to calculate the energy flux to better compare it with the numerical dissipation, which occurs mostly below this scale [the effective resolution of the model defined in Soufflet et al. (2016)].

Figures 9a–d depict the four energy sinks described above. The largest energy dissipation is from the top drag (negative $F_e K_e$), i.e., the transfer from oceanic mesoscale currents to the atmosphere, over the Gulf Stream (Fig. 9a), consistent with Renault et al. (2016a). The positive values of $F_e K_e$ over the shelf areas represent the generation of wind-driven shallow water currents (e.g., Renault et al. 2009). Not surprisingly, $F_b K_b$ is negative everywhere, with higher magnitudes in regions where strong currents interact with the topography, e.g., over the Charleston Bump and New England Seamounts

(Fig. 9b). The numerical dissipation $H_{\text{Diff}}$ is generally smaller than $F_e K_e$ and $F_b K_b$ and, like $F_b K_b$, tends to be larger over the rugged topography along the Gulf Stream path. Finally, the interior dissipation $I_{\text{Diss}}$ is one order of magnitude smaller than the other energy sinks (note the scale factor applied in Fig. 9d) and, consistent with our previous results, has its largest values in the Gulf Stream North Wall (see Fig. 3a).

For a more quantitative comparison, Fig. 10 presents the spatial average of $F_e K_e$, $F_b K_b$, $H_{\text{Diff}}$, and $I_{\text{Diss}}$ over the entire domain, the Gulf Stream region after separation (black contour in Fig. 9 using the $500 \text{ cm}^2 \text{ s}^{-2}$ EKE contour estimated from AVISO) and the Gulf Stream before separation (yellow contour in Fig. 9 delineating the area between isobaths 200 and 780 m and between 28.6° and 36.5°N). Averaging is performed over regions deeper than 200 m (thick dark green line in Fig. 9) to represent only mesoscale eddy dissipation. In the full domain, $F_b K_b$ is the main dissipation process, while $H_{\text{Diff}}$ and $F_e K_e$ represent 74% and 43% of $F_b K_b$, respectively, and $I_{\text{Diss}}$ is one order of magnitude smaller. On the Gulf Stream after separation (where the largest energy sinks are located), $F_e K_e$ has a similar magnitude to $F_b K_b$ due to strong mesoscale activity at the surface, while $H_{\text{Diff}}$ is smaller with 53% of $F_b K_b$. Again, the forward cascade is an order of magnitude smaller than the other terms. West of Cape Hatteras, in the region around the Charleston Bump, dissipation by bottom drag is the main energy sink as expected. $H_{\text{Diff}}$ is still high with about 30% of $F_b K_b$, due to the presence of strong currents, but both $F_e K_e$ and $I_{\text{Diss}}$ are one order of magnitude smaller than these two terms. The low values of $F_e K_e$ can be explained by a re-energization of the ocean by the wind over the Charleston Bump (see Fig. 9).

### b. Sensitivity to momentum advection schemes

In our simulation (as in Ménesguen et al. 2018), horizontal momentum advection is discretized with a fifth-order (UP5) rather than the more usual third-order (UP3) upstream-biased advection (Shchepetkin and McWilliams 1998). UP3 has the advantage of a lower computational cost due to a decrease in the required computations and exchanges between parallel subdomains due to a relatively compact stencil. However, this is achieved at the cost of a higher truncation error, and thus higher numerical dissipation and lower effective resolution relative to UP5 (Soufflet et al. 2016; Ménesguen et al. 2018).

In a little more detail, the dominant truncation error term of UP3 for a one-dimensional problem in the $x$ direction is

$$ \frac{|u| \Delta x^3 \partial^4 u}{12 \partial x^4}, $$

which appears as a hyperdiffusion operator of order 4 and hyperviscosity $|u|\Delta x^3/12$ ($\Delta x$ is the mesh size). Soufflet et al. (2016) demonstrate that this numerical diffusion has the ability to precisely dampen the second-largest error term, which is dispersive and would produce unwanted and potentially explosive noise in the solution. The same qualities are found for UP5, but in this case the dominant error is
<figure><img src="image_9.png" /><figcaption>FIG. 9. Time-averaged (2005–09) kinetic energy dissipation: (top) $F_eK_e$, (top middle) $F_bK_b$, (bottom middle) $H_{Diff}$, and (bottom) $10I_{Diss}$, computed from simulations with (a)–(d) UP5 and (e)–(h) UP3 numerical advection schemes. The eastern Gulf Stream region is bounded by the 500 cm$^2$ s$^{-2}$ EKE contour from AVISO (black contour). The western region is the area between the 200- and 780-m isobaths and between 28.6° and 36.5°N (yellow contour). The green contours indicate the 1000-, 3000-, and 200-m (thicker line) isobaths. Note that the energy dissipation color bar is nonlinear.</figcaption></figure>

$$ \frac{|u|\Delta x^5 \partial^6 u}{60 \partial x^6}, $$

that is, a hyperdiffusion of order 6 and hyperviscosity $|u|\Delta x^5/60$. In this case, numerical diffusion decreases faster with resolution and concentrates at scales closer to the grid scale, that is, about $5\Delta x$ instead of $10\Delta x$, according to the dispersion relation of these schemes and the KE spectra performed with the idealized simulations of Ménesguen et al. (2018).

In the following, to evaluate the impact of the momentum horizontal advection schemes on the different energy sinks, the very same simulation (initially with UP5) is performed for a period of 5 years but using UP3. $F_eK_e$, $F_bK_b$, $H_{Diff}$, and $I_{Diss}$ are then estimated from that simulation and compared with the UP5 simulation. Figures 9e–h show that the spatial patterns of the different energy sinks are relatively similar in UP5 and UP3 simulations. However, as expected, the average

over the entire domain (excluding the shelf seas) indicates that the numerical dissipation in the UP3 simulation increases by 17.5% relative to UP5 (Fig. 10). The effect on the mean EKE is milder, with only about 3% decrease (Figs. 10a,b). The reason is that the numerical schemes also indirectly modify the other energy sinks, with compensation effects. Going from UP5 to UP3, the domain-average $F_eK_e$ and $F_bK_b$ are reduced by about 5% and 3%, respectively (Fig. 10c). This decrease is related to the greater damping of the eddies, which results in a lower loss by top drag and bottom drag. Similar results are found over the Gulf Stream boxes, with $F_eK_e$ particularly affected after GS separation.

Quite surprisingly, the dissipation associated with the forward cascade ($I_{Diss}$) increases by 24% with the more diffusive UP3 scheme. This is counterintuitive, as we would expect less resolved energy to produce less energy flux. A possible explanation is that the reduction of mean currents with UP3 (by
<figure><img src="image_10.png" /><figcaption>FIG. 10. Time-averaged (2005–09) (left), (left center) EKE and (right center), (right) KE dissipation computed from simulations with UP5 (blue bar) and UP3 (orange bar), and spatially averaged over (a)–(d) the entire domain, (e)–(h) the eastern GS region of the Gulf Stream defined by the black outline in Fig. 9, and (i)–(l) the western GS region defined by the yellow outline in Fig. 9; (b), (f), and (j) are zoomed-in panels of (a), (e), and (i) respectively, and (d), (h), and (l) are zoomed-in panels of $I_{\text{Diss}}$ in (c), (g), and (k), respectively. In each panel, the error bar shows the standard error estimated by the bootstrap method. Note that the y axis is different in each region.</figcaption></figure>

about 1.5%) also reduces the current feedback effect, resulting in a stronger mean surface stress. As mentioned in section 4b, a larger surface stress would favor frontogenesis, and thus enhance the forward cascade. EBF was estimated and found to be 5% higher with UP3 than UP5 over the entire domain and more than 6% in the Gulf Stream region. However, other explanations for a more intense downscale flux when numerical dissipation is increased are possible. In Marchesiello et al. (2011), for example, an increase of grid resolution extends the spectrum range of forward cascade to smaller scales, but its intensity tends to decrease slightly. In this case, there is no current feedback to the atmosphere, and the system behaves as if an extradissipation near the grid scale creates a vacuum for cross-scale fluxes. We leave this question for further studies in more idealized framework.

The general conclusion of this section is that numerical dissipation, even at a relatively high order of discretization, is still much larger than the resolved downscale flux. We would expect them to be similar for numerical dissipation to act as turbulent dissipation at the grid scale. Considering depths below 100 m, where the kinetic energy more generally follows an inverse cascade (e.g., Klein et al. 2008; Schubert et al. 2020), the problem is even more challenging and would require not only further reduction of numerical errors in advection schemes, but also consideration of subgrid-scale backscatter models to inject the missing energy.

## 6. Conclusions

In this study, we analyze the cross-scale kinetic energy flux in the Gulf Stream, using a coarse-grained approach on a sub-mesoscale permitting (2 km) ocean model. This method allows us to describe the spatial and temporal variations of energy fluxes. From our results, a forward cascade tends to dominate everywhere near the surface at scales below 22 km, while the inverse cascade, generally associated with geostrophic turbulence, extends down to scales of at least 30 km. At larger scales, the inverse cascade dominates everywhere except on a narrow band extending from the Charleston Bump to the separation region along the Gulf Stream North Wall. To better understand what drives the kinetic energy flux, the total surface current was decomposed into its rotational (balanced) and divergent (unbalanced) components using a Helmholtz decomposition. The results show that the inverse and forward cascades are driven by the balanced and unbalanced flow, respectively. The forward cascade is temporally and spatially consistent with strong frontogenesis on the GS North Wall and associated ageostrophic secondary circulation. Next, we show that the kinetic energy flux is marked by seasonal variability. In winter, the forward cascade at 22 km is intensified everywhere, possibly due to the development of mixed layer instabilities, whereas the inverse cascade is intensified later in spring as submesoscales reinforce the mesoscale activity, consistent with previous studies.
Eddy-rich resolution models can exhibit excessive kinetic energy if not dissipated efficiently. How much energy is dissipated and how it flows through ocean scales is still an open question. In this study, we compare various energy sinks associated with boundary processes (top and bottom drag) and the interior route with dissipation through a downscale flux of resolved energy (forward cascade). Another process, numerical dissipation, is spurious in the sense that its role is to provide numerical integrity to the model, possibly at the expense of physical accuracy, but we generally assume that it is not a dominant process. It turns out that, in our model, the interior route to dissipation is one order of magnitude smaller than the dissipation produced by boundary processes and numerical discretization. Two important conclusions for numerical modeling can be drawn from these results. (i) The main energy dissipation processes in the real ocean are at the boundaries, and both top and bottom drags must be carefully accounted for in the models. Neglecting, in particular, the top drag process leads to an excess of eddy kinetic energy by ≈30% and overestimation of the inverse cascade (Renault et al. 2019; Aluie et al. 2018). (ii) The turbulent viscosity used in the explicit or implicit diffusion operators of current models largely overestimates the dissipation required by the physical mechanisms.

We tested the sensitivity of numerical momentum advection schemes using fifth- and third-order upstream-biased schemes. As expected, UP3 is more dissipative than UP5, but the difference is somewhat compensated by a reduced dissipation in the boundary layers. Weaker surface dissipation with UP3 can even enforce the mean surface stress, increase frontogenesis and the forward cascade. Therefore, efforts to reduce numerical dissipation to more physically meaningful values (with higher-order schemes and higher resolution) may give more weight to both the top drag and the bottom drag, perhaps simultaneously increasing the inverse cascade, but not necessarily the forward cascade. However, another dissipation process that UP5/UP3 may account for (besides eddy damping) is the breaking of internal waves generated by storms or tides. This process is not included in our estimate of downscale fluxes but may be important, although it is also probably greatly overestimated by current advection schemes (Jouanno et al. 2016).

Previous studies have estimated the rate of energy dissipation in symmetrically unstable fronts through observations. Using wintertime surveys of the Gulf Stream, Thomas et al. (2013) estimate that the dissipation associated with symmetric instabilities—following larger-scale energy-draining processes—is between $O(1)$ and $O(10)$ mW m$^{-2}$ depending on the intensity of the front. In the Kuroshio region, D'Asaro et al. (2011) found dissipation rates to be $O(100)$ mW m$^{-2}$ for a symmetrically unstable front, the high intensity explained by conditions favorable to a large EBF. Using the winter snapshot of our simulation, we observe that the downscale energy flux in the front, which would be dissipated by symmetric instabilities if the resolution allowed, is consistent with the magnitudes estimated by Thomas et al. (2013) and can sometimes be as high as $O(100)$ mW m$^{-2}$ (not shown). These magnitudes are comparable to the dissipation produced by bottom and top drag and numerical dissipation, so the interior dissipation may be locally substantial. However, due to the large spatial and temporal variability of the

fronts, the interior dissipation decreases significantly when estimating its long-term average contribution (i.e., greater than one month), making it significantly lower than other energy sinks.

Our study has several caveats. The spatial resolution of our simulations is too coarse to allow for ageostrophic submesoscale instabilities. Previous studies have suggested that symmetric (Dong et al. 2021; Thomas et al. 2013) or centrifugal (Gula et al. 2016) instabilities may contribute to the forward cascade. However, this would occur below the 20-km scale that appears to mark the beginning of a pseudodissipation range associated with frontogenesis. It is therefore unclear whether frontal instabilities would simply take over the direct cascade initiated by frontogenesis without much change in the dissipation rate, or whether they would restructure the turbulent forward cascade and the associated dissipation rate. Second, we do not explicitly consider the effects of air-sea interactions, but rather rely for our simulation on a parameterization of the current feedback. Considering that top and bottom drags are key processes for mesoscale energy dissipation, a coupled model could improve the representation of air-sea exchanges and their effect on mesoscale and submesoscale activity. However, our validation of mean and eddy flows against satellite and drifter data seems already beyond the accuracy of these observations, and it is unclear how much improvement can be gained or how to evaluate it with current observations.

Current satellite observations do not allow for an accurate characterization of the turbulent cascade (Klein et al. 2019). Satellite products such as AVISO are highly smoothed and can only resolve eddies with a radius greater than about 40 km and a lifetime greater than one week (Chelton et al. 2011). Arbic et al. (2013) and Renault et al. (2019) show that this limitation strongly affects the representation of the turbulent cascade, with a large underestimation of the cross-scale KE fluxes and a shift to lower wavenumbers of the forward cascade. Progress should be made as mesoscale and submesoscale motions with scales down to 30 km should be observable by the upcoming SWOT (Surface Water and Ocean Topography) altimetry mission (Morrow et al. 2019). However, the ageostrophic flow will be incompletely monitored and entangled with internal gravity waves (Klein et al. 2019). Future satellite missions—such as Odysea (Bourassa et al. 2016; Rodríguez et al. 2018)—will aim to consistently measure total surface currents (geostrophic and ageostrophic currents) and surface stress. This may contribute to our understanding of the energy pathways in the ocean.

*Acknowledgments.* Author Contreras was supported by Becas Chile (CONICYT-PFCHA/Doctorado Becas Chile/2020-72210196). We appreciate support from the Centre National d'Étude Spatial through the TOSCA projects I-CASCADE and CARAMBA, the LEFE VENUS project, and from the GdR and GdRI CROCO. The model and calculations were executed at GENCI (project 13051).

*Data availability statement.* The data to reproduce the figures are available online (https://doi.org/10.6084/m9.fishshare.20193092.v2).
REFERENCES

Aluie, H., M. Hecht, and G. K. Vallis, 2018: Mapping the energy cascade in the North Atlantic Ocean: The coarse-graining approach. *J. Phys. Oceanogr.*, **48**, 225–244, https://doi.org/10.1175/JPO-D-17-0100.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.

Barkan, R., K. B. Winters, and S. G. L. Smith, 2015: Energy cascades and loss of balance in a reentrant channel forced by wind stress and buoyancy fluxes. *J. Phys. Oceanogr.*, **45**, 272–293, https://doi.org/10.1175/JPO-D-14-0068.1.

Boccaletti, G., R. Ferrari, and B. Fox-Kemper, 2007: Mixed layer instabilities and restratification. *J. Phys. Oceanogr.*, **37**, 2228–2250, https://doi.org/10.1175/JPO3101.1.

Bourassa, M. A., E. Rodriguez, and D. Chelton, 2016: Winds and currents mission: Ability to observe mesoscale air/sea coupling. *2016 IEEE Int. Geoscience and Remote Sensing Symp.*, Beijing, China, Institute of Electrical and Electronics Engineers, 7392–7395, https://doi.org/10.1109/IGARSS.2016.7730928.

Brüggemann, N., and C. Eden, 2015: Routes to dissipation under different dynamical conditions. *J. Phys. Oceanogr.*, **45**, 2149–2168, https://doi.org/10.1175/JPO-D-14-0205.1.

Bryan, F. O., C. W. Böning, and W. R. Holland, 1995: On the midlatitude circulation in a high-resolution model of the North Atlantic. *J. Phys. Oceanogr.*, **25**, 289–305, https://doi.org/10.1175/1520-0485(1995)025<0289:OTMCIA>2.0.CO;2.

—, M. W. Hecht, and R. D. Smith, 2007: Resolution convergence and sensitivity studies with North Atlantic circulation models. Part I: The western boundary current system. *Ocean Modell.*, **16**, 141–159, https://doi.org/10.1016/j.ocemod.2006.08.005.

Callies, J., G. Flierl, R. Ferrari, and B. Fox-Kemper, 2016: The role of mixed-layer instabilities in submesoscale turbulence. *J. Fluid Mech.*, **788**, 5–41, https://doi.org/10.1017/jfm.2015.700.

Capet, X., P. Klein, B. L. Hua, G. Lapeyre, and J. C. Mcwilliams, 2008a: Surface kinetic energy transfer in surface quasi-geostrophic flows. *J. Fluid Mech.*, **604**, 165–174, https://doi.org/10.1017/S0022112008001110.

—, J. C. McWilliams, M. J. Molemaker, and A. F. Shchepetkin, 2008b: Mesoscale to submesoscale transition in the California Current system. Part I: Flow structure, eddy flux, and observational tests. *J. Phys. Oceanogr.*, **38**, 29–43, https://doi.org/10.1175/2007JPO3671.1.

—, —, —, and —, 2008c: Mesoscale to submesoscale transition in the California Current system. Part II: Frontal processes. *J. Phys. Oceanogr.*, **38**, 44–64, https://doi.org/10.1175/2007JPO3672.1.

—, —, —, and —, 2008d: Mesoscale to submesoscale transition in the California Current system. Part III: Energy balance and flux. *J. Phys. Oceanogr.*, **38**, 2256–2269, https://doi.org/10.1175/2008JPO3810.1.

—, G. Roullet, P. Klein, and G. Maze, 2016: Intensification of upper-ocean submesoscale turbulence through Charney baroclinic instability. *J. Phys. Oceanogr.*, **46**, 3365–3384, https://doi.org/10.1175/JPO-D-16-0050.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.

Chassignet, E. P., and D. P. Marshall, 2008: Gulf Stream separation in numerical ocean models. *Ocean Modeling in an*

Eddying Regime, Geophys. Monogr., Vol. 177, Amer. Geophys. Union, 39–62, https://doi.org/10.1029/177GM05.

—, and X. Xu, 2017: Impact of horizontal resolution (1/12° to 1/50°) on Gulf Stream separation, penetration, and variability. *J. Phys. Oceanogr.*, **47**, 1999–2021, https://doi.org/10.1175/JPO-D-17-0031.1.

—, L. T. Smith, G. R. Halliwell, and R. Bleck, 2003: North Atlantic simulations with the Hybrid Coordinate Ocean Model (HYCOM): Impact of the vertical coordinate choice, reference pressure, and thermobaricity. *J. Phys. Oceanogr.*, **33**, 2504–2526, https://doi.org/10.1175/1520-0485(2003)033<2504:NASWTH>2.0.CO;2.

Chelton, D. B., M. G. Schlax, and R. M. Samelson, 2011: Global observations of nonlinear mesoscale eddies. *Prog. Oceanogr.*, **91**, 167–216, https://doi.org/10.1016/j.pocean.2011.01.002.

Couvelard, X., P. Marchesiello, L. Gourdeau, and J. Lefèvre, 2008: Barotropic zonal jets induced by islands in the south-west Pacific. *J. Phys. Oceanogr.*, **38**, 2185–2204, https://doi.org/10.1175/2008JPO3903.1.

D'Asaro, E., C. Lee, L. Rainville, R. Harcourt, and L. Thomas, 2011: Enhanced turbulence and energy dissipation at ocean fronts. *Science*, **332**, 318–322, https://doi.org/10.1126/science.1201515.

Debreu, L., P. Marchesiello, P. Penven, and G. Cambon, 2012: Two-way nesting in split-explicit ocean models: Algorithms, implementation and validation. *Ocean Modell.*, **49–50**, 1–21, https://doi.org/10.1016/j.ocemod.2012.03.003.

—, N. K.-R. Kevlahan, and P. Marchesiello, 2022: Improved Gulf Stream separation through brinkman penalization. *Ocean Modell.*, **179**, 102121, https://doi.org/10.1016/j.ocemod.2022.102121.

Dong, J., B. Fox-Kemper, H. Zhang, and C. Dong, 2021: The scale and activity of symmetric instability estimated from a global submesoscale-permitting ocean model. *J. Phys. Oceanogr.*, **51**, 1655–1670, https://doi.org/10.1175/JPO-D-20-0159.1.

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.

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.

Fox-Kemper, B., R. Ferrari, and R. Hallberg, 2008: Parameterization of mixed layer eddies. Part I: Theory and diagnosis. *J. Phys. Oceanogr.*, **38**, 1145–1165, https://doi.org/10.1175/2007JPO3792.1.

Germano, M., 1992: Turbulence: The filtering approach. *J. Fluid Mech.*, **238**, 325–336, https://doi.org/10.1017/S0022112092001733.

Gula, J., M. J. Molemaker, and J. C. Mcwilliams, 2015: Gulf Stream dynamics along the southeastern U.S. seaboard. *J. Phys. Oceanogr.*, **45**, 690–715, https://doi.org/10.1175/JPO-D-14-0154.1.

—, —, and —, 2016: Submesoscale dynamics of a Gulf Stream frontal eddy in the South Atlantic bight. *J. Phys. Oceanogr.*, **46**, 305–325, https://doi.org/10.1175/JPO-D-14-0258.1.

Haine, T. W. N., and J. Marshall, 1998: Gravitational, symmetric, and baroclinic instability of the ocean mixed layer. *J. Phys. Oceanogr.*, **28**, 634–658, https://doi.org/10.1175/1520-0485(1998)028<0634:GSABIO>2.0.CO;2.

Hoskins, B. J., and F. P. Bretherton, 1972: Atmospheric frontogenesis models: Mathematical formulation and solution. *J. Atmos.*
Sci., **29**, 11–37, https://doi.org/10.1175/1520-0469(1972)029<0011:AFMMFA>2.0.CO;2

Jamet, Q., B. Deremble, N. Wienders, T. Uchida, and W. Dewar, 2021: On wind-driven energetics of subtropical gyres. *J. Adv. Model. Earth Syst.*, **13**, e2020MS002329, https://doi.org/10.1029/2020MS002329.

Jouanno, J., X. Capet, G. Madec, G. Roullet, and P. Klein, 2016: Dissipation of the energy imparted by mid-latitude storms in the Southern Ocean. *Ocean Sci.*, **12**, 743–769, https://doi.org/10.5194/os-12-743-2016.

Kessouri, F., D. Bianchi, L. Renault, J. C. McWilliams, H. Frenzel, and C. A. Deutsch, 2020: Submesoscale currents modulate the seasonal cycle of nutrients and productivity in the California Current system. *Global Biogeochem. Cycles*, **34**, e2020GB006578, https://doi.org/10.1029/2020GB006578.

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.

Klein, P., B. L. Hua, G. Lapeyre, X. Capet, S. Le Gentil, and H. Sasaki, 2008: Upper ocean turbulence from high-resolution 3D simulations. *J. Phys. Oceanogr.*, **38**, 1748–1763, https://doi.org/10.1175/2007JPO3773.1.

—, and Coauthors, 2019: Ocean-scale interactions from space. *Earth Space Sci.*, **6**, 795–817, https://doi.org/10.1029/2018EA000492.

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.

Laurindo, L. C., A. J. Mariano, and R. Lumpkin, 2017: An improved near-surface velocity climatology for the global ocean from drifter observations. *Deep-Sea Res. I*, **124**, 73–92, https://doi.org/10.1016/j.dsr.2017.04.009.

Lellouche, J.-M., and Coauthors, 2018: Recent updates to the Coernicus Marine Service global ocean monitoring and forecasting real-time 1/12° high-resolution system. *Ocean Sci.*, **14**, 1093–1126, https://doi.org/10.5194/os-14-1093-2018.

Lemarié, F., L. Debreu, A. Shchepetkin, and J. C. McWilliams, 2012: On the stability and accuracy of the harmonic and biharmonic isoneutral mixing operators in ocean models. *Ocean Modell.*, **52–53**, 9–35, https://doi.org/10.1016/j.ocemod.2012.04.007.

Leonard, A., 1975: Energy cascade in large-eddy simulations of turbulent fluid flows. *Adv. Geophys.*, **18**, 237–248, https://doi.org/10.1016/S0065-2687(08)60464-1.

Lévy, M., P. J. S. Franks, and K. S. Smith, 2018: The role of submesoscale currents in structuring marine ecosystems. *Nat. Commun.*, **9**, 4758, https://doi.org/10.1038/s41467-018-07059-3.

Marchesiello, P., J. C. McWilliams, and A. Shchepetkin, 2001: Open boundary conditions for long-term integration of regional oceanic models. *Ocean Modell.*, **3**, 1–20, https://doi.org/10.1016/S1463-5003(00)00013-5.

—, L. Debreu, and X. Couvelard, 2009: Spurious diapycnal mixing in terrain-following coordinate models: The problem and a solution. *Ocean Modell.*, **26**, 156–169, https://doi.org/10.1016/j.ocemod.2008.09.004.

—, X. Capet, C. Menkes, and S. C. Kennan, 2011: Submesoscale dynamics in tropical instability waves. *Ocean Modell.*, **39**, 31–46, https://doi.org/10.1016/j.ocemod.2011.04.011.

McWilliams, J. C., 2016: Submesoscale currents in the ocean. *Proc. Roy. Soc., 472A*, 20160117, https://doi.org/10.1098/rspa.2016.0117.

—, 2021: Oceanic frontogenesis. *Annu. Rev. Mar. Sci.*, **13**, 227–253, https://doi.org/10.1146/annurev-marine-032320-120725.

—, J. Gula, and M. J. Molemaker, 2019: The Gulf Stream north wall: Ageostrophic circulation and frontogenesis. *J. Phys. Oceanogr.*, **49**, 893–916, https://doi.org/10.1175/JPO-D-18-0203.1.

Ménesguen, C., S. Le Gentil, P. Marchesiello, and N. Ducousso, 2018: Destabilization of an oceanic meddy-like vortex: Energy transfers and significance of numerical settings. *J. Phys. Oceanogr.*, **48**, 1151–1168, https://doi.org/10.1175/jpo-d-17-0126.1.

Molemaker, M. J., J. C. McWilliams, 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/S00221102009993272.

Morrow, R., and Coauthors, 2019: Global observations of fine-scale ocean surface topography with the Surface Water and Ocean Topography (SWOT) mission. *Front. Mar. Sci.*, **6**, 232, https://doi.org/10.3389/fmars.2019.00232.

Olson, D. B., O. B. Brown, and S. R. Emmerson, 1983: Gulf Stream frontal statistics from Florida straits to Cape Hatteras derived from satellite and historical data. *J. Geophys. Res.*, **88**, 4569–4577, https://doi.org/10.1029/JC088iC08p04569.

Özgökmen, T. M., and E. P. Chassignet, 2002: Dynamics of two-dimensional turbulent bottom gravity currents. *J. Phys. Oceanogr.*, **32**, 1460–1478, https://doi.org/10.1175/1520-0485(2002)032<1460:DOTDTB>2.0.CO;2.

Penven, P., V. Echevin, J. Pasapera, F. Colas, and J. Tam, 2005: Average circulation, seasonal cycle, and mesoscale dynamics of the Peru Current system: A modeling approach. *J. Geophys. Res.*, **110**, C10021, https://doi.org/10.1029/2005JC002945.

Qiu, B., S. Chen, P. Klein, H. Sasaki, and Y. Sasai, 2014: Seasonal mesoscale and submesoscale eddy variability along the North Pacific subtropical countercurrent. *J. Phys. Oceanogr.*, **44**, 3079–3098, https://doi.org/10.1175/JPO-D-14-0071.1.

Renault, L., B. Dewitte, M. Falvey, R. Garreaud, V. Echevin, and F. Bonjean, 2009: Impact of atmospheric coastal jet off central Chile on sea surface temperature from satellite observations (2000–2007). *J. Geophys. Res.*, **114**, C08006, https://doi.org/10.1029/2008JC005083.

—, M. J. Molemaker, J. Gula, S. Masson, and J. C. McWilliams, 2016a: 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.

—, —, J. C. McWilliams, A. F. Shchepetkin, F. Lemarié, D. Chelton, S. Illig, and A. Hall, 2016b: Modulation of wind work by oceanic current interaction with the atmosphere. *J. Phys. Oceanogr.*, **46**, 1685–1704, https://doi.org/10.1175/JPO-D-15-0232.1.

—, P. Marchesiello, S. Masson, and J. C. McWilliams, 2019: Remarkable control of western boundary currents by eddy killing, a mechanical air-sea coupling process. *Geophys. Res. Lett.*, **46**, 2743–2751, https://doi.org/10.1029/2018GL081211.

—, S. Masson, T. Arsouze, G. Madec, and J. C. McWilliams, 2020: Recipes for how to force oceanic model dynamics. *J. Adv. Model. Earth Syst.*, **12**, e2019MS001715, https://doi.org/10.1029/2019MS001715.

Rodríguez, E., A. Wineteer, D. Perkovic-Martin, T. Gál, B. W. Stiles, N. Niamsuwan, and R. R. Monje, 2018: Estimating ocean vector winds and currents using a Ka-band pencil-
beam Doppler scatterometer. *Remote Sens.*, **10**, 576, https://doi.org/10.3390/rs10040576.

Saha, S., and Coauthors, 2010: The NCEP Climate Forecast System Reanalysis. *Bull. Amer. Meteor. Soc.*, **91**, 1015–1058, https://doi.org/10.1175/2010BAMS3001.1.

Schoonover, J., and Coauthors, 2016: North Atlantic barotropic vorticity balances in numerical models. *J. Phys. Oceanogr.*, **46**, 289–303, https://doi.org/10.1175/JPO-D-15-0133.1.

Schubert, R., J. Gula, R. J. Greatbatch, B. Baschek, and A. Biastoch, 2020: The submesoscale kinetic energy cascade: Mesoscale absorption of submesoscale mixed layer eddies and frontal downscale fluxes. *J. Phys. Oceanogr.*, **50**, 2573–2589, https://doi.org/10.1175/JPO-D-19-0311.1.

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.

Shchepetkin, A. F., and J. C. McWilliams, 1998: Quasi-monotone advection schemes based on explicit locally adaptive dissipation. *Mon. Wea. Rev.*, **126**, 1541–1580, https://doi.org/10.1175/1520-0493(1998)126<1541:OMASBO>2.0.CO;2.

—, and J. McWilliams, 2003: A method for computing horizontal pressure-gradient force in an oceanic model with a non-aligned vertical coordinate. *J. Geophys. Res.*, **108**, 3090, https://doi.org/10.1029/2001JC001047.

—, and J. C. McWilliams, 2005: The Regional Oceanic Modeling System (ROMS): A split-explicit, free-surface, topography-following-coordinate oceanic model. *Ocean Modell.*, **9**, 347–404, https://doi.org/10.1016/j.ocemod.2004.08.002.

—, and —, 2009: Correction and commentary for “ocean forecasting in terrain-following coordinates: Formulation and skill assessment of the regional ocean modeling system” by Haidvogel et al., *J. Comp. Phys.* **228**, 8985–9000, https://doi.org/10.1016/j.jcp.2009.09.002.

Siegelman, L., 2020: Energetic submesoscale dynamics in the ocean interior. *J. Phys. Oceanogr.*, **50**, 727–749, https://doi.org/10.1175/JPO-D-19-0253.1.

Soufflet, Y., P. Marchesiello, F. Lemarié, J. Jouanno, X. Capet, L. Debreu, and R. Benshila, 2016: On effective resolution in ocean models. *Ocean Modell.*, **98**, 36–50, https://doi.org/10.1016/j.ocemod.2015.12.004.

Spall, M. A., 1996: Dynamics of the Gulf Stream/deep western boundary current crossover. Part I: Entrainment and recirculation. *J. Phys. Oceanogr.*, **26**, 2152–2168, https://doi.org/10.1175/1520-0485(1996)026<2152:DOTGSW>2.0.CO;2.

Srinivasan, K., J. C. McWilliams, M. J. Molemaker, and R. Barkan, 2019: Submesoscale vortical wakes in the lee of topography. *J. Phys. Oceanogr.*, **49**, 1949–1971, https://doi.org/10.1175/JPO-D-18-0042.1.

—, R. Barkan, and J. C. McWilliams, 2022: A forward energy flux at submesoscales driven by frontogenesis. *J. Phys. Oceanogr.*, **53**, 287–305, https://doi.org/10.1175/JPO-D-22-0001.1.

Su, Z., J. Wang, P. Klein, A. F. Thompson, and D. Menemenlis, 2018: Ocean submesoscales as a key component of the global heat budget. *Nat. Commun.*, **9**, 775, https://doi.org/10.1038/s41467-018-02983-w.

—, H. Torres, P. Klein, A. F. Thompson, L. Siegelman, J. Wang, D. Menemenlis, and C. Hill, 2020: High-frequency submesoscale motions enhance the upward vertical heat transport in the global ocean. *J. Geophys. Res. Oceans*, **125**, e2020JC016544, https://doi.org/10.1029/2020JC016544.

Thomas, L. N., and C. M. Lee, 2005: Intensification of ocean fronts by down-front winds. *J. Phys. Oceanogr.*, **35**, 1086–1102, https://doi.org/10.1175/JPO2737.1.

—, and J. R. Taylor, 2010: Reduction of the usable wind-work on the general circulation by forced symmetric instability. *Geophys. Res. Lett.*, **37**, L18606, https://doi.org/10.1029/2010GL044680.

—, —, R. Ferrari, and T. M. Joyce, 2013: Symmetric instability in the Gulf Stream. *Deep-Sea Res. II*, **91**, 96–110, https://doi.org/10.1016/j.dsr2.2013.02.025.

Uchida, T., 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.

Vallis, G. K., 2017: Atmospheric and Oceanic Fluid Dynamics: Fundamentals and Large-Scale Circulation. 2nd ed. Cambridge University Press, 946 pp., https://doi.org/10.1017/9781107588417.