## Abstract

A large fraction of the kinetic energy in the ocean is stored in the “quasigeostrophic” eddy field. This “balanced” eddy field is expected, according to geostrophic turbulence theory, to transfer energy to larger scales. In order for the general circulation to remain approximately steady, instability mechanisms leading to loss of balance (LOB) have been hypothesized to take place so that the eddy kinetic energy (EKE) may be transferred to small scales where it can be dissipated. This study examines the kinetic energy pathways in fully resolved direct numerical simulations of flow in a flat-bottomed reentrant channel, externally forced by surface buoyancy fluxes and wind stress in a configuration that resembles the Antarctic Circumpolar Current. The flow is allowed to reach a statistical steady state at which point it exhibits both a forward and an inverse energy cascade. Flow interactions with irregular bathymetry are excluded so that bottom drag is the sole mechanism available to dissipate the upscale EKE transfer. The authors show that EKE is dissipated preferentially at small scales near the surface via frontal instabilities associated with LOB and a forward energy cascade rather than by bottom drag after an inverse energy cascade. This is true both with and without forcing by the wind. These results suggest that LOB caused by frontal instabilities near the ocean surface could provide an efficient mechanism, independent of boundary effects, by which EKE is dissipated. Ageostrophic anticyclonic instability is the dominant frontal instability mechanism in these simulations. Symmetric instability is also important in a “deep convection” region, where it can be sustained by buoyancy loss.

## 1. Introduction and motivation

The general circulation of the ocean is forced by surface fluxes of momentum, heat, and freshwater at basin scales. A large fraction of the kinetic energy *E*_{k} associated with the large-scale forcing must be dissipated at molecular scales in order for the circulation to remain approximately steady. The *E*_{k} pathways across this wide range of scales remain poorly understood (Ferrari and Wunsch 2009). Possible routes to dissipation include nonlinear internal wave interactions in the interior, instability of geostrophic motions, and direct interactions with side and bottom boundaries (Müller et al. 2005; Ferrari and Wunsch 2010). About 90% of the *E*_{k} in the oceans is stored in the geostrophic eddy field (Ferrari and Wunsch 2009) and, given this large fraction, we focus here on the “instability” pathway to dissipation.

Geostrophic eddies are formed, primarily, via baroclinic instability of large-scale ocean currents that are in approximate thermal wind balance (Gill et al. 1974). These geostrophic currents and eddies are considered “balanced” in the sense that they are in geostrophic and hydrostatic balance, and the associated velocity and buoyancy fields are entirely determined by the potential vorticity (Ford 1993). The quasi-two-dimensional geostrophic eddy field is expected, based on geostrophic turbulence theory, to transfer energy to larger scales (Salmon 1980). Kinetic energy dissipation, however, is observed to take place at small scales (Wunsch and Ferrari 2004 and references therein). Instability mechanisms leading to loss of balance (LOB) are therefore hypothesized to take place in order for some of the kinetic energy of the geostrophic eddies (EKE) to be transferred to small scales where it can be dissipated.

Recent numerical studies (Capet et al. 2008c; Molemaker et al. 2010) have suggested that loss of balance due to submesoscale instabilities can provide an efficient route to small-scale dissipation. To explore this hypothesis, we pose an idealized problem of flow in a flat-bottomed reentrant channel, externally forced by surface buoyancy fluxes and wind stress, in a configuration that resembles the Antarctic Circumpolar Channel (ACC). The flow is allowed to reach a statistical steady state at which point it exhibits both forward and inverse energy cascades. Flow interactions with irregular bathymetry are excluded so that bottom drag is the sole mechanism available to dissipate the upscale energy cascade. This simple setting allows us to compare the *E*_{k} energy pathways along two main routes: the “inverse” route—in which EKE is transferred to larger scales where it is dissipated via bottom drag—and the “forward” route—in which EKE cascades to small dissipative scales, presumably via loss of balance. We analyze these competing routes to dissipation for two different forcing scenarios and evaluate the relative efficiency of each pathway.

Capturing both the inverse and forward pathways while computing kinetic energy dissipation accurately requires a wide range of scales to be resolved and is a challenging numerical task. We are therefore forced to make compromises in the problem setup. We evaluate the validity of our numerical simulations in representing realistic oceanic flows based on five nondimensional parameters (see section 2 for details) and find them to be in quantitative agreement with the oceanic dynamical regime. In addition, we diagnose the mean buoyancy structure and the mean and eddy-driven circulations and find them to be in qualitative agreement with zonal-mean theories of the ACC (see section 3 for details). This suggests that understanding the physical mechanisms illuminated by this process study is helpful in improving our dynamical understanding of real ocean phenomena.

The paper is organized as follows: in section 2 we discuss the problem setup; in section 3 we describe the overall flow features; in section 4 we show the *E*_{k} balance; in section 5 we provide spectral analyses of velocity variance, *E*_{k} equation, and *E*_{k} fluxes; in section 6 we evaluate loss of balance in relation to EKE dissipation; and finally, in section 7 we provide a summary and discussion.

## 2. Problem setup

We consider an idealized ocean basin on an *f* plane in a rectangular channel of volume *V* and depth *H*, with zonal and meridional dimensions *L*_{x} = *L*_{y}. The vertical coordinate is −*H* ≤ *z* ≤ 0, and density is expressed as *ρ* = *ρ*_{0}(1 − *g*^{−1}*b*), where *b* is the buoyancy. The corresponding Cartesian and Boussinesq equations are

The pressure is *ρ*_{0}*p*, *f* > 0 is the Coriolis frequency, *κ* is the diffusivity, *ν* is the kinematic viscosity, *r* (s^{−1}) is the coefficient of bottom drag, *f*_{bot}(*z*) is a near-bottom localization function discussed in the appendix, and are the unit vectors in the zonal (*x*), meridional (*y*), and vertical (*z*) directions. No penetration conditions are imposed on the top, bottom, and meridional sides, where **u** = (*u*, *υ*, *w*) and is the outward normal to the surface of *V*. The flow is taken to be periodic in the zonal direction. Free-slip and no-flux scalar conditions are prescribed on the bottom and meridional sides.

Our analysis is based on two direct numerical simulations (DNS) in the configuration described above. The first simulation, hereinafter referred to as buoyancy forced (BF), is forced at the surface solely by a buoyancy flux of the form

where *B*_{max} (m^{2} s^{−3}) is the magnitude of the maximal buoyancy flux applied at the surface and *λ*_{y} = *L*_{y}/8; (2) corresponds to supplying positive buoyancy anomalies near *y* = 0 and negative buoyancy anomalies near *y* = *L*_{y}, as is shown in Fig. 1.

The second simulation, hereinafter referred to as wind and buoyancy forced (WBF), is forced by the same buoyancy flux, but in addition a surface stress

is applied at the top, where *μ* is the dynamic viscosity and *ρ*_{0}*τ*_{max} (N m^{−2}) is the magnitude of the maximal surface stress. The Ekman pumping −*f*^{−1}∂*τ*/∂*y* associated with (3) corresponds to the downwelling positively buoyant fluid near *y* = 0 and the upwelling negatively buoyant fluid near *y* = *L*_{y} (see Fig. 1). This configuration resembles the ACC in the sense that the wind stress curl acts to increase available potential energy (APE) and is therefore “eddy favoring” (Cessi 2007).

Both simulations were run for ~0.25 of a diffusive time (*H*^{2}/*κ*) by which [(8a, below)] had reached a statistical steady state (as shown in Fig. 5, described later). Time averages and all analyses presented in this paper are based on the last 0.05 diffusive time, during which all physical variables were stored every 6.5 inertial periods (*f*^{−1}) for a total of 58 snapshots. Throughout the manuscript we present probability distribution functions (PDFs) and cumulative distribution functions (CDFs) of various diagnostic parameters computed for mean and total quantities at selected regions of the domain. The PDFs and CDFs are computed using all grid points in the specified domain for all snapshots. For example, if a PDF is computed based on the entire volume *V* then the total number of samples is 65 × 1025 × 1024 × 58 ~ *O*(10^{9}). The PDFs in Figs. 2, 4, and 17 (as seen below) are normalized to lie between zero and one.

The model Eqs. (1) for both simulations are solved using the nonhydrostatic three-dimensional pseudospectral model flow_solve. See the appendix for more details on the model and problem setup.

Throughout the paper 〈( )〉 denotes a volume average, denotes a spatial average over the horizontal area *A* = *L*_{x} × *L*_{y}, denotes an average over the zonal length *L*_{x}, denotes a time average, and eddies (denoted by ′) are defined as perturbations from both time and zonal averages (denoted by ). Note that at statistical steady state .

To compare our simulations to the ocean we compute, at steady state, the resulting Rossby number Ro, the Ekman number Ek, and the Burger number Bu, defined as

where *ζ* is the vertical component of vorticity, subscripts denote derivatives, and *R*_{d} is the Rossby deformation radius. Two additional nondimensional numbers are required to completely define our problem, and we choose the externally determined parameters

where Pr is the Prandtl number and *α* is the aspect ratio.

The value for Ek represents the mean computed in the interior (0.8 < *y*/*H* < 14.2 and 0.06 < *z*/*H* < 0.94) away from the boundaries where frictional effects are expected to be important (Ek ~ 0.1 there). Of the Ek values, 90% are below 0.02. The value for Ro represents the mean computed over the entire domain; 90% of the values are below 0.05. The analyses shown in Figs. 8 and 10 below help us estimate the ratio *R*_{d}/*L*_{y} as *O*(10^{−1}) (see sections 5b and 5c for details) and the corresponding Bu as *O*(10^{−2}). Our simulations are thus well within the dynamical regime of the ocean with respect to all of the nondimensional numbers aside from the aspect ratio that is one to two orders of magnitude larger than that of typical ocean basins.

To verify that our simulations are in hydrostatic balance at basin scales as well as mesoscales, we follow Molemaker et al. (2010) who have suggested a measure to assess the degree of which the solutions are in hydrostatic balance using the parameter

where *μ* = (∂_{z}*p*)_{rms} + *b*_{rms} is added to the denominator of (7) to ensure that situations with locally weak vertical pressure gradients are not being identified as significantly unbalanced: regions where are thus considered to be hydrostatically balanced and regions where are considered highly nonhydrostatic. Figure 2 shows the PDF of computed over the entire volume. Of the values, 99% are smaller than 0.015 and 0.04 for BF and WBF, respectively, indicating that the flow is hydrostatically balanced at all scales. A similar analysis using the time- and zonal-mean fields (not shown) indicates that 99% of the mean values are below 0.005 and 0.009. These results ensure that the flow in our simulations is hydrostatically balanced at nearly all scales of motion. Figures 13e and 14e (see sections 6 and 7 for discussion) show that hydrostatic balance is partially lost at the sharp, near-surface, submesoscale fronts that develop at the edges of the mesoscale eddies.

The volume-averaged kinetic (*E*_{k}) and potential (*E*_{p}) energy equations for the model (1) take the form

where *ϵ* = *ν*〈|**∇ u**|

^{2}〉,

*ϵ*

_{d}=

*r*〈|

**u**

_{bot}|

^{2}〉, ,

**u**

_{bot}and

**u**

_{s}are the horizontal velocities (

*u*,

*υ*) evaluated at the bottom and surface, respectively, and is the reversible conversion term between kinetic and potential energies associated with the mean and the eddies, respectively. Similarly the Laplacian EKE dissipation is

*ϵ*′ =

*ν*〈|

**∇u**′|

^{2}〉, and EKE dissipation due to the bottom drag is . In the following sections, we will evaluate the relative efficiency of the “inverse” to “forward” routes by comparing the reservoirs of

*ϵ*,

*ϵ*′, and .

The spectral representation of the *E*_{k} equation is derived by taking the cospectra of the dot product of **u** and (1a) to give

where { }_{k} denotes a Fourier transform, ^{⋄} denotes a complex conjugate, and denotes the real part. The terms denote the dissipation (both Laplacian and bottom drag) and wind input, respectively.

## 3. Overall flow features

Figure 1 (left) shows the time-averaged buoyancy field for BF. The specified flux boundary condition [(2)] corresponds to buoyancy gain near *y* = 0 and buoyancy loss near *y* = *L*_{y}. Figure 3 shows the time- and zonal-mean isopycnal height (top), Eulerian-mean streamfunction (middle), and residual streamfunction *ψ*_{res} (streamfunction computed in buoyancy space and interpolated back to ** z** space). The residual streamfunction provides a useful framework for understanding buoyancy transport associated with the mean flow and the eddies (

*ψ**) (Wolfe 2014). For BF, the solutions consist of a plume (deep convection region) near

*y*/

*H*= 15 and a thermocline near

*y*/

*H*= 0. The residual streamfunction

*ψ*

_{res}is much larger than the Eulerian-mean streamfunction , indicating that buoyancy is mainly transported by the eddies (denoted by wavy arrows in Fig. 1, left). The depth of the thermocline is set by the interaction between the deep convection and baroclinic eddies (Barkan et al. 2013).

Figure 1 (right) shows the time-averaged zonal velocity field *u* for WBF. The spatial structure of the zonal jet is shown as well as the Ekman upwelling (*y* ≈ *L*_{y}) and downwelling (*y* ≈ 0) regions associated with the curl of the wind stress [(3)]. The thermocline depth in this case is determined by the balance between the wind-induced Ekman flow and the baroclinic eddies (Marshall and Radko 2003). Figure 3 shows that the thermocline, near *y*/*H* = 0, is deeper than in BF. Near *y*/*H* = 15 the interaction between the deep convection, wind-induced upwelling, and baroclinic eddies sets the isopycnal structure. The Eulerian-mean streamfunction , associated in WBF with the wind-induced Ekman flow (thermally indirect), is much stronger than the viscously driven one in BF. The magnitudes of *ψ*_{res} and are similar for WBF, but the circulation is reversed. This illustrates the eddies’ tendency to counteract the Ekman flow and produce a thermally direct flow. For nearly adiabatic flow in a channel configuration, *ψ*_{res} is expected to decrease as *κ* decreases (hence the name “residual”), a well-known result from GCM simulations of the ACC (Marshall and Speer 2012 and references therein). Note that in the thermocline region (0 ≤ *y*/*H* ≤ 6) *ψ*_{res} ≈ 0, suggesting that the flow is approximately adiabatic there. To assess how adiabatic our interior flow is for WBF, we examine the departure of the time and zonally averaged buoyancy [(1b)] from

using the parameter

where . Analogous to (7), regions where are completely adiabatic and regions where are completely diabatic. The PDF of the interior (e.g., in the domain 0.8 < *y*/*H* < 14.2 and 0.06 < *z*/*H* < 0.94) is plotted in Fig. 4. Of the values, 90% are below 0.05, indicating that diffusivity is negligible in the interior of our domain in agreement with zonally averaged theories of the ACC. Note that typical stratification values produced in the WBF simulation are lower than in the ACC, *N*/*f* ~ 5 compared with ~25 (Garabato et al. 2004). In reality, the thermal structure of the ACC depends on bathymetric detail, lateral exchanges with the basin to the north, and the shelf dynamics at the Antarctic margin. None of these processes are included in our model, while they are presumably necessary to quantitatively match observations.

## 4. Kinetic energy balance

Figure 5 shows the terms in (8a) demonstrating that, in both cases, the *E*_{k} reservoir has reached a statistical steady state. The ratio illustrates that the *E*_{k} reservoir is much larger for WBF. The ratio is 1.8 for BF and 0.27 for WBF. Bottom drag thus provides a much larger sink for *E*_{k} in WBF.

Figure 6 shows the decomposition of 〈*wb*〉, *ϵ*, and *ϵ*_{d} in (8a) into the mean and eddy components. For BF, the vertical buoyancy flux is accomplished by the eddies , and 80% of the total *E*_{k} is EKE. For WBF, . The wind generates APE by tilting the isopycnals , and the eddies release APE via baroclinic instability and restratify the fluid . Note that both of these terms are equal in magnitude to the total dissipation of EKE and to the vertical buoyancy flux by the Ekman velocity *w*_{E} = −(2*f*)^{−1 }*dτ*_{s}/*dy*. In this simulation, most of the wind work is used to form a barotropic jet in the zonal direction (Fig. 1). The energy associated with this jet (~75% of the total *E*_{k}) is dissipated by bottom drag . Only a small fraction of the wind work , sometimes referred to as “useful wind work” (Cessi et al. 2006), is expended in generating APE and is thus available to drive baroclinic eddies. Accordingly, the ratio , which is smaller than for the total fields. This result suggests a way to evaluate the relative importance of wind to surface buoyancy forcing by comparing their role in generating EKE. For BF, Fig. 6 and (8b) suggest that at steady state, , where *b*_{max} is the maximal specific buoyancy difference. We therefore propose the nondimensional number

as a measure for the relative importance of wind to buoyancy forcing in generating EKE in the ocean. For a typical ocean basin with *α* = 10^{−3}, *κ* = 10^{−5} m^{2} s^{−1}, *f* = 10^{−4} s^{−1}, and *τ* = 10^{−4} m^{2} s^{−2} (corresponding to a wind stress of 0.1 N m^{−2}) *S*_{EKE} ~ 100, illustrating that buoyancy fluxes provide only 1% of EKE compared with the useful wind work. In the WBF simulation *S*_{EKE} ~ 6, suggesting that our wind forcing is much weaker than that of the ACC (see section 7 for discussion).

The ratio for BF and 1.67 for WBF. This is the main result of this paper, demonstrating that both with and without mechanical forcing a large fraction of EKE dissipation is established independent of bottom drag. This suggests that in both cases the “forward route” to dissipation is more effective than the “inverse route” in removing eddy kinetic energy once the flow has reached a quasi-steady state. It remains to be shown that the Laplacian EKE dissipation *ϵ*′ is directly linked to a forward energy cascade and loss of balance.

## 5. Spectral analysis

### a. Horizontal velocity wavenumber spectra

Figure 7 shows the time-averaged horizontal velocity (*u* and *υ*) wavenumber spectra in the *k*_{x} (black) and *k*_{y} (gray) directions, computed both near the surface (solid lines) and at depth (dashed line). For both BF and WBF, the surface spectral shapes are nearly flat at low wavenumbers, exhibit a rolloff, and become linear over more than a decade of wavenumbers and then steepen further at higher wavenumber where horizontal isotropy is lost and the flow becomes more three-dimensional (the horizontally isotropic regions are gray shaded). The spectral shapes at depth closely resemble those at the surface, although the spectral level is reduced by about a decade over the entire wavenumber range (in WBF the amount of energy at the lowest meridional wavenumber is independent of depth, a signature of the barotropic zonal jet).

For BF, the spectral slope after the rolloff is close to *k*^{−3}, whereas in WBF it is closer to *k*^{−2}. A *k*^{−3} spectral slope is typically associated with an enstrophy inertial subrange in quasigeostrophic turbulence theory (Salmon 1980). A *k*^{−2} spectral slope has previously been observed in submesoscale-resolving simulations (Capet et al. 2008a, their Fig. 6) and was interpreted as an indication of a significant amount of kinetic energy in the submesoscale range. This suggests that the energy reservoir in the wavenumber band associated with submesoscale features is larger for WBF.

### b. Kinetic energy spectral balance

Figure 8 shows the individual terms in the spectral energy balance (9). The cospectra are vertically integrated, time averaged, and are computed as a function of the horizontal wavenumber . Solid lines denote total fields and dashed lines denote eddy fields.

For BF, there is very little difference between total and eddy fields. Positive (blue line) illustrates the release of APE by baroclinic instability. The peak around *k*_{h}*L*_{y} ~ 5 corresponds to the most unstable baroclinic mode *k*_{bc} (vertical dashed lines on the right panels) and can be used to approximate the effective Rossby radius of deformation wavenumber *k*_{d}. Note that, in effect, *k*_{d} is probably a few times larger than *k*_{bc} [for the classical Eady (1949) baroclinic instability analysis *k*_{d} ≈ 3.9*k*_{bc}]. The green line is the momentum advection term that, by definition, is identically zero when integrated over the entire wavenumber range. It thus only describes a redistribution of energy across wavenumbers. Note that the transition between negative and positive energy flux divergence occurs around *k*_{bc}. The Laplacian EKE dissipation term (black line) and the EKE dissipation due to bottom drag (red line) extract roughly the same amount of energy for *k*_{h}*L*_{y} < 7. However, for higher wavenumbers Laplacian EKE dissipation exceeds bottom drag EKE dissipation.

In contrast, for WBF there is a clear difference between total and eddy fields at low wavenumbers (*k*_{h}*L*_{y} < 2.5), where most of the energy is in the mean flow. Wind work (magenta line) supplies energy at large scales, primarily to the barotropic zonal jet (Fig. 1) whose energy is dissipated by bottom drag (red line). Negative values of (blue line) indicate the generation of APE by the wind, corresponding to negative values in Fig. 6. For *k*_{h}*L*_{y} > 2.5, the energy is mainly in the eddy field and the picture resembles qualitatively that of BF. Quantitatively, the reservoir of *E*_{k} is larger (note the different *y* axis range between top and bottom panels), and *k*_{bc} is slightly larger as well. Assuming *k*_{d} ≈ 2*k*_{bc}, which is less than the 3.9 value in the classical Eady (1949) problem, results in *k*_{d}*L*_{y} ≈ 10. This means that our domain contains roughly 10 Rossby radii of deformation, leading to our estimate of the Bu in (4) as being *O*(10^{−2}). Bottom drag EKE dissipation exceeds Laplacian dissipation in the range *k*_{h}*L*_{y} < 10. However, as in BF, for higher wavenumbers the EKE dissipation by bottom drag is smaller. Note that the observed elevated Laplacian dissipation at small horizontal wavenumbers (large horizontal scales) implies that anisotropic motions with small vertical scales are being dissipated.

### c. Kinetic energy spectral fluxes: Forward and inverse cascades

Forward and inverse energy transfers can be expressed explicitly in terms of the *E*_{k} spectral flux:

As before, the cospectra are depth integrated, time averaged, and are computed as a function of the horizontal wavenumber *k*_{h}. Figures 9 and 10 show Π(*k*_{h}) for BF and WBF respectively (top-right panels). By definition, . Positive/negative slopes (∂Π/∂*k*_{h}) indicate wavenumber bands where the flux is divergent/convergent. As a result, regions where Π(*k*_{h}) > 0 denote regions of forward energy transfer, and regions where Π(*k*_{h}) < 0 denote regions of inverse energy transfer. Wavenumber bands denoted by *a*, *b*, and *c* thus correspond to scales of motion that exhibit forward, inverse, and forward energy transfers, respectively. Also displayed are the depth-averaged enstrophy *ζ*^{2} (top-left panels) and bandpass-filtered enstrophies corresponding to wavenumber bands *a*, *b*, and *c* (the applied filters excluded narrow wavenumber bands around the zero crossings).

For BF, band *a* is associated with large-scale motions in the deep convection region. Band *b* corresponds to “mesoscale” geostrophic eddies that, in agreement with geostrophic turbulence theory, are seen to transfer energy to larger scales. Note that EKE dissipation due to both bottom drag and Laplacian dissipation have maxima (*k*_{h}*L*_{y} = 4.2) near the transition wavenumber (*k*_{h}*L*_{y} ~ 3.5) between regions *a* and *b* (Fig. 8, top). Our interpretation is that the maximum in Laplacian EKE dissipation is an energy sink associated with the forward transfer in the deep convection region, and the maximum in bottom drag EKE dissipation corresponds to the arrest of the inverse energy cascade associated with the geostrophic eddies. Band *c* corresponds to the sharp fronts that develop at the edges of the eddies and deep convection filaments, and we argue in section 6 that the flow features in this wavenumber range are associated with the loss of balance.

For WBF (Fig. 10), the dynamics associated with wavenumber bands *b* and *c* and the corresponding bandpass-filtered enstrophies are qualitatively similar to BF. Quantitatively, the ratio of forward to inverse energy fluxes is larger than in BF, suggesting that loss of balance at these scales is more vigorous in WBF. This is in agreement with the shallower spectral slope of WBF compared with that of BF (Fig. 7). The dynamics in band *a* are now dominated by the large-scale wind-driven jet. As indicated by Fig. 8 (bottom left), the total kinetic energy at these scales (associated with the zonal jet) is removed by bottom drag (solid red line). As in BF, we interpret the maximum bottom drag EKE dissipation (dashed red line in bottom-right panel) observed in the transition wavenumber between bands *a* and *b* (*k*_{h}*L*_{y} = 4.2) to be associated with the halting of the inverse energy cascade. Laplacian dissipation (dashed black line) is more evenly distributed across the wavenumber with a much broader maximum than in BF. We note, however, that inferring dynamics solely from spectral fluxes can be misleading, particularly in situations where the energy transfers are nonlocal in wavenumber space.

## 6. Loss of balance

Figures 11a,b and 12a,b show snapshots of the Richardson number and the Rossby number Ro = *ζ*/*f*. Positive *O*(1) values of Ro and Ri are observed in the fronts that develop at the edges of the larger-scale eddies and in the filaments within the deep convection zone. These positive *O*(1) values indicate regions of ageostrophic motion (Thomas et al. 2008) and correspond spatially to elevated levels of EKE and Laplacian *E*_{k} dissipation (Figs. 11c,d and 12c,d). The correlation between small-scale dissipation and ageostrophic motions suggests that LOB is related to a forward energy cascade at frontal scales corresponding to the wavenumber band *c* in Figs. 9 and 10.

McWilliams and Yavneh (1998), McWilliams (2003), and Molemaker et al. (2005) suggested criteria for LOB based on the change of type of the balanced equations, from elliptic to hyperbolic, which occurs under any of the following conditions:

*N*^{2}changes sign (becomes negative);the absolute vorticity in isopycnal coordinates

*A*=*f*+*ζ*changes sign (becomes negative for*f*> 0);change of sign of

*A*− |*S*| (becomes negative for*f*> 0), where*S*^{2}= (*u*_{X}−*υ*_{Y})^{2}+ (*υ*_{X}+*u*_{Y})^{2}is the variance of the horizontal strain rate in isopycnal coordinates; orthe Richardson number Ri < ¼.

The physical mechanisms suggested to be associated with these LOB conditions are convective instability (CI) for condition i, inertial instability (INI) or symmetric instability (SI) for condition ii, anticyclonic–ageostrophic instability (AAI) for condition iii (Müller et al. 2005), and classical Kelvin–Helmholtz instability (KHI) for condition iv (Miles 1961; Howard 1961). Note that for flow in hydrostatic balance

where PV is the potential vorticity and the subscript *H* denotes horizontal components. Condition ii can thus be expressed as the change in sign of PV (becomes negative for *f* > 0). As discussed in Hoskins (1974) and Bennetts and Hoskins (1979), if the PV is negative because of the horizontal components PV_{H} [first term on the right-hand side of (14)], the instability is of type SI, and if the PV is negative because of the vertical components PV_{V} [second term on the right-hand side of (14)], the instability is of type INI. We will thus distinguish between INI and SI when examining condition ii. Note that according to ii and iii above if the condition for SI or INI is satisfied, the condition for AAI must also be satisfied because *A* − |*S*| ≤ *A*.

Other forms of submesoscale instabilities are clearly present in our numerical simulations. Mixed layer instability (MLI) (Fox-Kemper et al. 2008) is an obvious one, as indicated by the *O*(1) Richardson number values in Figs. 11a and 12a. MLI is important for the restratification of the mixed layer (Boccaletti et al. 2007) and may develop into finite-amplitude mixed layer eddies that can lead to intense frontogenesis, a precondition for LOB. Nevertheless, it is a balanced instability that does not lead directly to LOB (Molemaker et al. 2005; Thomas et al. 2008); that is, additional instability mechanisms, such as the ones described above, must take place for LOB to occur. The indirect pathways to LOB due to MLI and other balanced phenomena are therefore assumed to be accounted for by the instability mechanisms i–iv above.

Capet et al. (2008b) have suggested an assessment for the degree of balance based on the departure of the divergence of the horizontal momentum equation (McWilliams 1985) from

using the parameter

where the *H* subscripts denote horizontal velocities and gradients. The term is added to the denominator of (16) to exclude situations with locally weak force divergences from being identified as significantly unbalanced. Regions where are thus considered highly balanced and regions where are considered highly unbalanced.

Figures 13 and 14 show surface slices of the conditions for LOB (i–iii) in BF and WBF, respectively. The contour plot in Figs. 11a and 12a shows condition iv for LOB. Regions where *N*^{2} < 0 are not plotted in conditions ii–iv. For SI, only the regions where PV < 0 and PV_{V} > 0 are plotted. Similarly, for INI, only regions where the PV < 0 and PV_{H} > 0 are plotted.

Qualitatively, in BF, on large scales, LOB occurs because of CI in the deep convection region near *y*/*H* > 12 (Fig. 13a). At smaller scales, within the deep convection region, LOB occurs in some locations due to SI (Fig. 13b) and more extensively due to AAI, as indicated by the fine, blue filaments in Fig. 13c. The condition for INI is never satisfied (Fig. 13d). The condition KHI is rarely satisfied, as indicated by the hardly visible black contours in Fig. 11a.

In WBF, LOB due to SI and more so due to AAI occurs over a much greater portion of the domain, extending throughout the thermocline region (*y*/*H* < 5). The condition for INI is, again, never satisfied. LOB due to CI is still observed in the deep convection region (now *y*/*H* > 10) and, in addition, across the jet (6 < *y*/*H* < 9) at much smaller scales. In both flows, there are regions, especially near fronts, where PV ~ *A* > 0 but *A* − |*S*| < 0, emphasizing the importance of elevated strain variance in triggering LOB in these locations. The condition for KHI is satisfied more extensively as well, primarily in the jet region (contour lines in Fig. 12a). Note, however, that the level of EKE in the jet region is low compared with levels at the thermocline and deep convection regions (Fig. 12c). Visually, regions of unbalanced, nonhydrostatic flow (Figs. 14e and 14f) correspond spatially to the regions of elevated Laplacian dissipation shown in Fig. 12d. This is consistent with the notion that the forward energy cascade at frontal scales (wavenumber band *c* in Figs. 9 and 10) and the corresponding Laplacian dissipation are dynamically linked with LOB.

Figure 15 shows the CDF of (16), denoted by , computed based on the full zonal and meridional extents in the range 0.92 ≤ *z*/*H* ≤ 1 (black) and 0.3 ≤ *z*/*H* ≤ 0.6 (gray). For both BF and WBF, higher values of are more probable near the surface, indicating that LOB is more likely to occur there than in the interior.

Figure 16 shows the time- and zonal-mean Laplacian EKE dissipation . For BF, occurs primarily near the surface and is concentrated in the deep convection region. This agrees with Fig. 15 and the fact that LOB primarily occurs within the deep convection region (Fig. 13).

For WBF, is most intense near the surface (dashed line) in the thermocline region and to a lesser extent in the deep convection region with a weak signature at depth. Interestingly, little is observed in the jet region in agreement with Fig. 12c. This observation, together with the spatial structure of Laplacian dissipation and LOB instability mechanisms (Figs. 12d, 14), suggests that LOB leads to the observed forward energy cascade at frontal scales and the consequent small-scale EKE dissipation.

To quantify which of the LOB instabilities is most probable we have computed the CDFs of each of the LOB conditions (i–iv). In both simulations, the CDFs were computed near the top (0.92 ≤ *z*/*H* ≤ 1), extending over the full zonal and meridional extents. For WBF, we have also computed the CDFs in the same vertical and zonal ranges, however, limiting the meridional extent to the thermocline region. Table 1 summarizes the results of the CDF analysis. Figure 17 shows the CDFs for the WBF simulations.

For BF, in agreement with Fig. 13, CI is the most probable instability condition followed by AAI and SI. The conditions for INI and KHI are essentially never satisfied (for INI the probability is, in fact, identically zero).

For WBF, the probability for all conditions increases by an order of magnitude, however, the condition for INI is rarely satisfied. CI is the most probable instability condition followed by AAI, SI, and KHI. In the thermocline region, however, the probability for AAI exceeds all other conditions by an order of magnitude. The reduction in CI is clear from Fig. 14a that shows that *N*^{2} > 0 in the thermocline region. CI is thus most probable in the jet and convection regions. SI is a “forced” instability and can only be sustained if it is maintained by winds and/or buoyancy fluxes (Thomas and Lee 2005). Otherwise, SI will tend to restratify the surface layer and adjust the background flow to a state of marginal stability with PV ≈ 0, thus damping itself out (Taylor and Ferrari 2009). This is clearly seen in the PDFs of SI (inset in Fig. 17), where the most probable PV values are positive and very close to zero. The reduction in the probability of SI in the thermocline region compared with the entire domain suggests that it is sustained, primarily, by buoyancy loss in the deep convection region. Finally, as depicted in Fig. 12a (contours), the Richardson number is less than a quarter, primarily in the jet region. This explains the reduction in the probability for KHI within the thermocline region.

## 7. Summary and discussion

We have investigated two idealized DNS of channel flow with different external forcing: one forced solely by surface buoyancy fluxes (BF) and the other by both surface buoyancy fluxes and wind stress (WBF). In both experiments there are two possible routes to *E*_{k} dissipation: upscale energy cascade leading to dissipation due to bottom drag (the inverse route) and downscale energy cascade, catalyzed by loss of balance (LOB), leading to Laplacian dissipation at small scales (the forward route).

In WBF, a strong barotropic zonal jet is directly forced by the wind stress. As a result, the total kinetic energy reservoir is 7 times larger than in BF. Only a fraction of that wind work, however, that is, useful wind work (Cessi et al. 2006), is converted to available potential energy and then drained to excite geostrophic eddies. Consequently, the eddy kinetic energy (EKE) reservoir in WBF is only twice that of BF. Furthermore, in both simulations, the amount of EKE dissipation due to Laplacian dissipation at small scales (*ϵ*′) exceeds that due to bottom drag , suggesting that the forward route dominates. Interestingly, the ratio between *ϵ*′ and is about 1.6 in both cases.

Analyses of the Rossby number, Richardson number, violation of the different balance constraints, and of *E*_{k} spectral fluxes indicate that *ϵ*′ is directly linked to ageostrophic motions (Ro ≥ 1, Ri ≤ 1) that undergo LOB instabilities and exhibit a forward energy cascade. Spatial maps of time- and zonal-mean EKE dissipation are correlated with regions of LOB at frontal (submeso) scales. As a result, for BF, EKE dissipation is strongest in the upper domain within the deep convection region. For WBF, where the forward energy cascade at submesoscales and the degree of LOB are much larger, EKE dissipation is still strongest near the surface (e.g., in the “surface layer”) and particularly within the thermocline region. In both cases a CDF of the parameter quantifying the degree of LOB [(16)] indicates that LOB occurs much more frequently in the surface layer than in the interior.

Analyses of the CDFs and spatial maps of the different LOB conditions (section 6) that lead to the forward energy transfers indicate that, in BF, convective instability (CI) is most probable on the large scales followed by ageostrophic anticyclonic instability (AAI) and to a lesser extent symmetric instability (SI), both of which act at submesoscales. Inertial instability (INI) and Kelvin Helmholts instability (KHI) are essentially never found. All instability mechanisms in this case take place within the deep convection region.

For WBF, the probability for all conditions is increased by an order of magnitude compared with BF. In addition, different instabilities are more/less probable at different meridional regions. In the deep convection region (10 < *y*/*H*), CI is most probable at large scales, followed by AAI and SI that act at submesoscales. In the jet region (6 < *y*/*H* < 9), CI, AAI, and KHI are all diagnosable, although the amount of EKE there is very low. In the thermocline region (*y*/*H* < 5), where most of the EKE is found and where most of the EKE dissipation takes place, AAI is the most probable instability, exceeding all others by an order of magnitude. This suggests that AAI is a more prevalent mechanism for LOB at these submesoscales. It is important to note that because of our model output frequency (≈6.5 × *f*^{−1}), we may underestimate the importance of SI that equilibrates on much faster time scales than AAI, thereby wiping itself out as it resets PV to zero. Nevertheless, we argue that it is the lack of external forcing needed to sustain symmetric instability that causes the probability reduction in the thermocline region. The fact that in the deep convection region where buoyancy loss can sustain SI the probabilities are high supports this argument.

It is noteworthy that a variety of balanced processes, such as mixed layer instability, are important to catalyze frontogenesis and thus precondition the flow to be more susceptible to LOB. Nevertheless, additional instability mechanisms, such as the ones described above, must take place for LOB to occur. Therefore, we have not explicitly analyzed all possible preconditioning processes and assumed LOB to be accounted for by the specific instability mechanisms listed in section 6.

### Oceanographic implications

The numerical challenge of resolving the wide range of scales necessary to capture both forward and inverse energy cascades while computing dissipation accurately is considerable, particularly for an externally forced problem that requires a long time to equilibrate. We therefore had to make compromises in the problem setup. The results of such an exercise are intended to provide insight into the physics of nearly balanced flows that are dynamically similar to oceanic flows and that drive both up- and down-scale energy cascades. The similarity to oceanic flows has been demonstrated by diagnosing the quasi-steady simulations to be nearly hydrostatically balanced and with small characteristic values of the Rossby, Burger, and Ekman numbers. Though the nearly adiabatic nature of the WBF simulation and the overall thermal structure in the thermocline region are in reasonable qualitative agreement with flow in the ACC (Marshall and Speer 2012, and references therein), the strength of the stratification produced in the simulation is weaker than that observed. This is partially because of processes associated with bathymetric detail and lateral exchange with other basins that are not included in our model. Whether or not higher stratification would alter the relative transfer rates of the forward and inverse energy cascades by suppressing certain frontal instability mechanisms is not clear. In the deep convection zone the dynamics are perhaps more relevant to weakly stratified regions subject to deep convection such as the Labrador Sea (Jones and Marshall 1997).

The useful wind work concept lead us to propose the nondimensional parameter *S*_{EKE} [(12)] to account for the relative contribution of wind to buoyancy forcing to EKE generation. The parameter *S*_{EKE} ~ 100 for representative ACC parameters, suggesting that wind forcing is responsible for most of the EKE generation, with buoyancy forcing providing only 1%. In our WBF simulation, *S*_{EKE} is much smaller (~6). Figures 12a and 12b show that wind forcing greatly increases/decreases the Ro/Ri values at the fronts. This suggests that higher *S*_{EKE} values may, in fact, lead to higher/lower values of Ro/Ri, that is, stronger ageostrophic motions, which will have even greater susceptibility to LOB and forward energy cascade.

One can better relate our simulations to ocean scales based on our approximation *k*_{d}*L*_{y} ≈ 10 (see section 5b for detail), meaning that the domain in our simulations contains roughly 10 Rossby radii of deformation. This estimate is in agreement with the enstrophy structures associated with bandpass filter (*b*) in Figs. 9 and 10. Our horizontal grid spacing is thus ≈*R*_{d}/100 and that for high latitudes (*R*_{d} = 15–25 km) would be about 150–250 m. The vertical resolution needs to be sufficient to resolve the elevated EKE dissipation in the surface “mixed” layer (dashed line in Fig. 16). In our case, seven grid points were required to properly resolve the surface layer. Note that a low-order finite difference model would require many more grid points than a spectral model to accurately resolve the same flow features. Finally, Fig. 14e shows that hydrostatic balance is partially lost in the sharp frontal features that develop at the edge of the baroclinic eddies. This suggests that nonhydrostatic dynamics may be important at these small scales and can affect the energy pathways and LOB instabilities, as has previously been noted by Mahadevan (2006). Magaldi and Haine (2015), for example, have found the amount of strain variance to be much larger in nonhydrostatic versus hydrostatic simulations, a feature that would undoubtedly increase the likelihood for AAI.

Figure 18 shows horizontal and vertical slices of a local region around a representative front in WBF (rectangular box in Fig. 12). Figure 18a shows the buoyancy field at the front with the predominantly geostrophic flow (solid arrows) in the vicinity of the high pressure, positively buoyant eddy. Figures 18c and 18e show the vertical velocity and EKE dissipation (log scale) in the same region. EKE dissipation is strongest at the front and correspond well to downwelling regions. Vertical slices of the same fields at *x*/*H* = 5.2 are shown in Figs. 18b, 18d, and 18f. As expected from classical theory (Hoskins and Bretherton 1972), ageostrophic circulation with downwelling on the cold side of the front and upwelling on the warm side is observed. The correlation between elevated dissipation and downwelling is found all through the thermocline region. The largest values of dissipation, however, are confined to the surface layer (dashed line), as expected from Fig. 16.

Although we do not attempt to perform a thorough analysis of all possible mechanisms that may lead to downwelling and small-scale dissipation, our results suggest that AAI may play an important role in this process. Understanding the dynamics associated with AAI is an active area of research (Ménesguen et al. 2012). It has been previously identified in high-resolution numerical studies of variable configurations and complexities (Mahadevan and Tandon 2006; Capet et al. 2008b; Molemaker et al. 2010), although its significance as a mechanism leading to LOB was not quantified.

The fact that EKE is dissipated preferentially at small scales near the surface, rather than at large scales via bottom drag suggests that the forward route associated with LOB may provide an important route to dissipation for the kinetic energy stored in the geostrophic eddy field. Although we do not incorporate LOB mechanisms associated with topographic interactions, except as a parameterized bottom drag, the spatial structure of small-scale EKE dissipation suggests that LOB caused by frontal instabilities near the ocean surface can be very efficient, independent of boundary effects.

It has been shown in models (Chen and Kamenkovich 2013) and in observations (Zajaczkovski and Gille 2014, manuscript submitted to *J. Geophys. Res.*) that baroclinic eddies in the ACC are typically generated in localized regions downstream of topography. Whether this localization would affect the degree of LOB or the instability mechanisms leading to it is unclear. Furthermore, Nikurashin et al. (2013) have demonstrated in numerical simulations that include realistic Southern Ocean topography that topographically generated internal waves play an important role in the total *E*_{k} dissipation. However, because they did not distinguish between mean and eddy fields it is hard to evaluate what fraction of the internal waves was excited by the mean flow and what fraction was excited by the eddies. In the ocean, LOB due to both frontal instability near the surface and internal wave generation near topography may provide a route leading to EKE dissipation. To our knowledge, no study has yet quantified which of the routes dominates.

Other phenomena such as downfront winds (Thomas and Lee 2005) and convection due to the diurnal cycle (Taylor and Ferrari 2010) may trigger or suppress certain frontal instability mechanisms. In addition, inertial gravity waves (Thomas 2012) and Stokes drift associated with the surface gravity wave field (McWilliams and Fox-Kemper 2013) may affect the frontal processes and complicate the picture further. These complexities, however, only emphasize the importance of understanding the physics of fronts, which lie at the cusp between balanced and unbalanced motions near the ocean surface.

## Acknowledgments

RB thanks George Carnevale for useful discussion on spectral fluxes. Comments from two anonymous referees greatly impoved this manuscript. This work was supported by the National Science Foundation, Grant Numbers OCE-0926481 and OCE-1259580. XSEDE computing resources at the National Institute of Computational Sciences, University of Tennessee were made available under Grant TG-OCE120004.

### APPENDIX

#### Details of the Numerical Simulations

The model Eqs. (1) are solved using the nonhydrostatic, three-dimensional, pseudospectral model flow_solve (Winters and de la Fuente 2012) on an equally spaced grid with (*nz* × *ny* × *nx*) = (65 × 1025 × 1024) grid points. The grid spacing is *dx* = *dy* = *dz*, ensuring that nonhydrostatic motions are completely resolved.

where *σ*_{b} = *dz*, and is defined in (2). In the limit of infinite resolution (*dz* → 0), the inhomogeneity in the boundary condition is exactly exchanged for inhomogeneity in the governing equation (Winters and de la Fuente 2012) so that (A1) and (2) are identical. The corresponding maximal surface buoyancy flux . Note that to ensure mass conservation.

Similarly, the surface stress [(3)] is applied as a body force in the zonal momentum equation

where is defined in (3). The width *σ*_{τ} = 6*dz* is chosen to distribute the stress-induced momentum over a resolvable thin layer near the top of thickness *O*(6*dz*), effectively specifying the thickness of the surface boundary layer. The corresponding maximal wind stress .

The form of the bottom drag *f*_{b}(*z*) in (1a) is

where *σ*_{d} = 6 *dz*, which smoothly confines the action of the drag term to a thin but well-resolved near-bottom layer, thus imposing the thickness of the bottom boundary layer. The magnitude of the bottom drag was set to be , where the constant was chosen such that the associated wavenumber , with as a representative kinetic energy flux in the simulations, was smaller than the deformation wavenumber *k*_{d} and larger than the domain wavenumber . This ensures that the simulations reach a steady state (Fig. 5) while allowing for a reasonable wavenumber band over which an inverse energy cascade takes place (see Figs. 9, 10). In oceanic simulations, typical linear drag values are *O*(10^{−7} − 10^{−6}) (Cessi et al. 2006) and the ratio *r*/*f* is *O*(10^{−3} − 10^{−2}). In our simulations, *r*/*f* ~ *O*(10^{−1}). This is because our domain is relatively small, which means we must have a larger magnitude of bottom drag in order to halt the inverse cascade before the eddies reach the domain size. Our results should therefore not be sensitive to the choice of bottom drag magnitude as long as the associated wavenumber satisfies *k*_{domain} < *k*_{r} < *k*_{d}.

To verify that all simulated fields were resolved down to the buoyancy variance dissipation scale, the spectra of the second derivative (highest derivative in the equations of motion) of buoyancy in each direction were analyzed. The value of Pr was 7 for both simulations, so the buoyancy field exhibits the smallest scales in these simulations and the second derivative of buoyancy is the hardest computed field to resolve. The decay in the spectra of the second derivative down to the highest wavenumber was the criterion used to ensure sufficient resolution. Table 2 provides the dimensional scales of the variables used to set up the simulations. Finally, the Kolmogorov microscale *l*_{ko} = (*ν*^{3}/EKE_{flux})^{¼} based on the EKE flux measured in our WBF simulation (Fig. 6) is *l*_{ko} ≈ 3*dy*. This is an indication that dissipative scales are adequately resolved. Furthermore, it illustrates that the forward cascade is reasonably inertial with peak values being ≈13*l*_{ko} (Fig. 10).

## REFERENCES

*Nonlinear Processes in Geophysical Fluid Dynamics,*O. V. Fuentes, J. Sheinbaum, and J. Ochoa, Eds., Kluwer, 287–304.

*Marine Turbulence,*H. B. J. Simpson and J. Sündermann, Eds., Cambridge University Press, 397–405.

*Ocean Modeling in and Eddying Regime, Geophys. Monogr.,*Vol. 177, Amer. Geophys. Union, 17–38.

## Footnotes

This article is included in the Ocean Turbulence Special Collection.