# Abyssal Slope Currents

Esther Capó, James Mcwilliams, Jonathan Gula, M. Jeroen Molemaker, Pierre Damien, René Schubert

► To cite this version:

Esther Capó, James Mcwilliams, Jonathan Gula, M. Jeroen Molemaker, Pierre Damien, et al. . Abyssal Slope Currents. Journal of Physical Oceanography, 2024, 54 (11), pp.2373-2392. 10.1175/JPO-D-24-0028.1 . hal-04841641

HAL Id: hal-04841641

https://inria.hal.science/hal-04841641v1

Submitted on 17 Dec 2024

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

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

<figure><img src="image_1.png" /></figure>
Abyssal Slope Currents

Esther Capó<sup>a</sup>, James C. McWilliams<sup>a</sup>, Jonathan Gula<sup>a,b</sup>, M. Jeroen Molemaker<sup>a</sup>, Pierre

Damien<sup>a</sup>, René Schubert<sup>b</sup>

<sup>a</sup> Department of Atmospheric and Oceanic Sciences, University of California, Los Angeles, CA, USA

<sup>b</sup> Univ Brest, CNRS, IRD, Ifremer, Laboratoire d'Océanographie Physique et Spatiale (LOPS), IUEM, Brest, France

<footnote>Corresponding author: Esther Capó, estherct@ucla.edu</footnote>
ABSTRACT: Realistic computational simulations in different oceanic basins reveal prevalent prograde mean flows (in the direction of topographic Rossby wave propagation along isobaths; a.k.a. topostrophy) on topographic slopes in the deep ocean, consistent with the barotropic theory of eddy-driven mean flows. Attention is focused on the Western Mediterranean Sea with strong currents and steep topography. These prograde mean currents induce an opposing bottom drag stress and thus a turbulent boundary-layer mean flow in the downhill direction, evidenced by a near-bottom negative mean vertical velocity. The slope-normal profile of diapycnal buoyancy mixing results in down-slope mean advection near the bottom (a tendency to locally increase the mean buoyancy) and up-slope buoyancy mixing (a tendency to decrease buoyancy) with associated buoyancy fluxes across the mean isopycnal surfaces (diapycnal downwelling). In the upper part of the boundary layer and nearby interior, the diapycnal turbulent buoyancy flux divergence reverses sign (diapycnal upwelling), with upward Eulerian mean buoyancy advection across isopycnal surfaces. These near-slope tendencies abate with further distance from the boundary. An along-isobath mean momentum balance shows an advective acceleration and a bottom-drag retardation of the prograde flow. The eddy buoyancy advection is significant near the slope, and the associated eddy potential energy conversion is negative, consistent with mean vertical shear flow generation for the eddies. This cross-isobath flow structure differs from previous proposals, and a new one-dimensional model is constructed for a topostrophic, stratified, slope bottom boundary layer. The broader issue of the return pathways of the global thermohaline circulation remains open, but the abyssal slope region is likely to play a dominant role.
# 1. Introduction

The abyssal ocean below the pycnocline is little explored relative to the upper ocean. It is understood to be generally full of mesoscale eddy currents that reach deeply down from the pycnocline, inertia-gravity waves generated at both the bottom and surface, stratified shear turbulence, a turbulent bottom boundary layer with topographic vortical wakes extending into the interior, submesoscale coherent vortices, and narrow deep western boundary currents that are part of the global thermohaline circulation. This paper is about additional phenomena that as yet are even less well understood.

The first of these phenomena is topostrophy (Holloway 2008), defined as a prograde time-mean along-isobath current over a sloping bottom (i.e., in the same direction as the phase propagation of topographic Rossby waves, with shallower water to the right of the current in the Northern Hemisphere where the Coriolis frequency *f* is positive, and vice versa in the Southern Hemisphere). Its measure is defined by

$$ \mathcal{T} = \operatorname{sgn}[f] (\hat{z} \cdot \mathbf{u} \times \nabla h), \qquad (1) $$

and its time average near the bottom, denoted by $\overline{\mathcal{T}}_b$, is usually positive. The subscript *b* denotes a near-bottom value (i.e., within the Monin-Obukhov turbulent layer, outside of any viscous sublayer that would support a no-slip boundary condition). The overbar is a time average; bold-face symbols are horizontal vectors; overlying arrows indicate a 3D vector; $\hat{z}$ is a unit vertical vector in the upward direction against gravity; and $h(\mathbf{x})$ is the resting depth of the ocean. A unit horizontal vector in the topostrophic direction is

$$ \hat{\mathbf{s}} = \operatorname{sgn}[f] \frac{\nabla h}{|\nabla h|} \times \hat{z}, \qquad (2) $$

whence $\mathcal{T} = \hat{\mathbf{s}} \cdot \mathbf{u} |\nabla h|$.

The theory for this phenomenon has its origins in topographic turbulence theory for two-dimensional flow over variable bathymetry (Bretherton and Haidvogel 1976; Salmon et al. 1976; Herring 1977). It shows that random eddying flows evolve to develop a negative correlation between the vertical vorticity, $\zeta = \partial_x v - \partial_y u$ where ($u, v$) are (eastward,northward) horizontal velocities, and resting-depth variations, $h-H$, where $H$ is a mean depth. This implies the development of a persistent topostrophic flow that can be described as an eddy-driven mean flow. The rationales for this outcome have variously been given as viscous enstrophy dissipation due to the forward cascade
of $\overline{\zeta'^2}$ in two-dimensional turbulence, entropy maximization in the sense of statistical mechanical equipartition of evolutionary end-states, and potential vorticity homogenization through horizontal eddy stirring. Topostrophy is also called the Neptune effect, and Greg Holloway has been a strong advocate for its importance (e.g., Holloway 1987, 1996). For oceanic basins, the most energetic eddies are mesoscale, and the largest $|\nabla h|$ values occur over continental slopes, seamounts, mid-ocean ridges, and islands, and this is where topostrophy might be most likely to occur; such flows would be in the cyclonic direction on basin edges and anticyclonic around interior bumps. Some striking examples have been reported (Alvarez et al. 1994; Saunders and King 1995; Weijer et al. 2020; Xie et al. 2022), and a climatological data set from bottom current meters supports the frequent occurrence of topostrophy on slopes (Holloway 2008). In a companion paper (Schubert et al. 2024) and in this one, high-resolution, realistic simulations are used to demonstrate the extreme prevalence of topostrophic flows on abyssal slopes.

The second phenomenon focused on in this paper is a prevalent downhill mean flow. If the along-isobath flow is topostrophic, then there is a horizontal drag stress on the bottom, $\tau = \rho C_d |\overline{\mathbf{u}_b}| \overline{\mathbf{u}_b}$, in the same direction, acting to retard the near-bottom $\overline{\mathbf{u}_b}$. This stress could cause a turbulent bottom (Ekman) boundary layer, whose mean horizontal transport velocity — ignoring the influence of stable density stratification and buoyancy mixing (but see Sec. 7) — is perpendicular to the bottom stress in the downhill direction, i.e., $\overline{\mathbf{u}_b} \cdot \nabla h > 0$. Alternatively expressed, the unstratified Ekman layer transport is to the left of the adjacent interior flow in the Northern Hemisphere (and vice versa in the South). Furthermore, the mean kinematic bottom boundary condition of no normal flow at a solid boundary is

$$ \overline{w}_b = - \overline{\mathbf{u}_b} \cdot \nabla h, \qquad (3) $$

and this implies that the vertical velocity $\overline{w}_b < 0$ at the bottom.

Care must be taken about details of the vertical profile of the current structure within the boundary layer. In this regard, the third primary phenomenon that we focus on is the mean diapycnal buoyancy mixing associated with these slope flows. An important idealized, one-dimensional (1D) bottom boundary layer (BBL) model (Wunsch 1970; Phillips 1970; Garrett 1990) for the effect of mixing in a stratified fluid with flat buoyancy surfaces in the interior next to a topographic slope has steady solutions with an upslope mean current at the bottom (hence $\overline{w}_b > 0$), a reverse downslope current above, and a retrograde along-isobath mean current in the adjacent interior — all as a consequence
<figure><img src="image_2.png" /><figcaption>FIG. 1. Sketch of the currents and buoyancy surfaces near a sloping boundary for a topostrophic interior flow when $f > 0$. Y and Z are the coordinates parallel and perpendicular to the slope, respectively. $h_{BBL}$ is the BBL thickness, and $V(Z)$ is the cross-slope velocity profile within the BBL; see (15).</figcaption></figure>

of no buoyancy flux through the boundary (ignoring the small geothermal heat flux (Thompson and Johnson 1996)); see Fig. 1. Stable stratification is expressed as $\partial_z b > 0$, where $b = -g\rho/\rho_0$ is the buoyancy, $\rho$ is the local potential density, and $z$ is the upward coordinate aligned with gravity. The direction of these currents is in contradiction to the topostrophic and Ekman-layer expectations in the preceding two paragraphs. This contradiction is resolved in Sec. 7.

The small-scale turbulent buoyancy flux,

$$ \vec{\mathcal{B}} = \vec{u}'b', \qquad (4) $$

(where primes denote deviations from mean values that must be defined in each context; in this case, deviations include sub-grid scale fluctuations due to parameterizations and model discretization errors) is “downward” in stable stratification, i.e., toward less buoyant water (e.g., by mixing warm surface water and abyssal cold water). Extensive micro-structure measurements (St. Laurent et al. 2001; Garabato et al. 2004; Waterhouse et al. 2014; Garabato et al. 2019) show that in the abyss above “rough” topography this mixing flux is often largest at a depth somewhat above the sloping bottom due to both boundary-layer turbulent mixing and nearby breaking of bottom-generated inertia-gravity waves or topographic wakes. The diapycnal mixing velocity across a buoyancy
102 iso-surface is defined by

$$ \tilde{\omega} = - \frac{\vec{\nabla} \cdot \vec{\mathcal{B}}}{|\vec{\nabla} \vec{b}|} \qquad (5) $$

103 Integrated over the horizontal area of a buoyancy surface in the abyss, $\tilde{\omega}$ must be positive both to serve as the upward return flow compensating polar sinking in the global thermohaline circulation, 104 and to account for the water mass transformation of less buoyant into more buoyant water. Invoking 105 the slope buoyancy-mixing BBL model cited above and the observed vertical structure of diapycnal 106 mixing near topography (e.g., Ferrari et al. 2016; McDougall and Ferrari 2017; Callies 2018; 107 Holmes and McDougall 2020; Drake et al. 2022; Baker et al. 2023), it is hypothesized that most 108 of the dominant positive $\tilde{\omega}$ occurs in the BBL next to the slopes and weaker negative values occur 109 in the interior, with a spatially inhomogeneous pattern concentrated near topographic slopes. This 110 is a modern revision of the older conceptions of positive $w$ and $\tilde{\omega}$ more homogeneously spread 111 throughout the abyssal interior (Stommel 1958; Munk 1966). 112

113 The intersection of these processes occurs principally in the mean buoyancy balance,

$$ \partial_t \bar{b} = - \overline{\vec{u} \cdot \vec{\nabla} b} - \vec{\nabla} \cdot \vec{\mathcal{B}} \qquad (6) $$

$$ = - \overline{\vec{u} \cdot \vec{\nabla} b} + \tilde{\omega} |\vec{\nabla} \vec{b}| \approx 0, \qquad (7) $$

114 especially near bottom slopes. The mean and eddy buoyancy advection approximately balances 115 the small-scale mixing, i.e., their tendencies must have opposite signs in equilibrium.

116 In this paper we partially synthesize these processes guided by realistic regional simulations, but 117 leaving the full characterization of the upward branches of the thermohaline circulation to future 118 global simulations. Our findings are that topostrophic flows are very common and that $\overline{w}$ and 119 $\tilde{\omega}$ are mostly negative along the abyssal slopes, but with sign reversals toward the interior. Our 120 focus is on the abyss because there can be many competing influences on the mean currents in 121 the upper ocean (e.g., the wind-driven western boundary current in a subtropical gyre); whereas, 122 mesoscale eddies are ubiquitous and are often the dominant influence on the mean currents in the 123 abyss (e.g., Rhines and Holland 1979). We use the Western Mediterranean Sea to illustrate this 124 behavior. Schubert et al. (2024) is a companion paper showing that this is common throughout 125 the Atlantic and Pacific Oceans, and Gula et al. (2024) explains the method of diagnosing $\vec{\mathcal{B}}$ in
a terrain-following coordinate model (as used here), inclusive of both explicitly parameterized mixing and the contribution from model discretization errors.

The organization of the paper is as follows. The general circulation of the Western Mediterranean Sea is summarized in Sec. 2; the simulation specifications are in Sec. 3; the near-slope circulation results are in Sec. 4; the mean along-isobath momentum $\hat{s} \cdot \bar{u}$ and buoyancy $\bar{b}$ balances are diagnosed in Secs. 5-6; a simple 1D BBL model combining diapycnal mixing and a topostrophic interior flow is in Sec. 7; and the results are further discussed and summarized in Secs. 8-9.

## 2. Western Mediterranean circulation

Topostrophy is a widespread phenomenon (Holloway 2008; Schubert et al. 2024). The analysis here is conducted in the Western Mediterranean basin, with a particular focus in the Alboran Sea, partly taking advantage of existing simulations and previous analyses (Capó et al. 2021; Capó and McWilliams 2022; Capó et al. 2023).

The Western Mediterranean Sea is a semi-enclosed basin located between the Strait of Gibraltar, 36°N–5.5°W, and the Strait of Sicily, 37°N–12°E (Fig.2a). The region is characterized by a complex bathymetry shaped by numerous capes, seamounts, ridges and islands of different shapes and sizes. While the interior of the basin is relatively flat, with an averaged depth not exceeding 3000 m, the continental and island/ridge margins are characterized by steep slopes, often shaped with narrow, deep canyons. The Alboran Sea is the transitional area between 5.5°W–1°W that connects the Western Mediterranean Sea with the Atlantic Ocean through the narrow Strait of Gibraltar (Fig.2a,b).

In the upper 150-200 m the large scale circulation in the Alboran Sea is mainly driven by the incoming Atlantic Jet, a strong fresh, surface current with velocities of $O(1)$ m s$^{-1}$ forced by a strong surface pressure gradient between both sides of the Strait and modulated by the tides affecting the Atlantic side. This jet flows along the basin shaping the quasi-permanent Western Anticyclonic Gyre, and an intermittent, less intense Eastern Anticyclonic Gyre (Viúdez et al. 1996; Vargas-Yáñez et al. 2002). To the east, the meeting of the Atlantic fresh current with the resident salty Mediterranean water generates a strong surface density (mainly salinity) contrast, commonly named the Almeria-Oran Front (Tintore et al. 1988; Allen et al. 2001). The Atlantic Jet exits the Alboran Sea and flows eastwards along the African coast towards the Strait of Sicily as the
Algerian Current. The Mediterranean thermohaline circulation transforms the incoming fresh Atlantic water into saltier water, which leaves the basin towards the Atlantic Ocean through the Strait of Gibraltar in the form of Levantine Intermediate Waters, within a 200-600 m depth range on average, and Deep Mediterranean Water below. These currents follow a complicated path that is highly constrained by the shape of topography: a major branch of the Levantine undercurrent flows along the Spanish margin towards the Strait, while the denser Deep Mediterranean Water occupies the central (and deepest) part of the basin and finally flows towards the Strait along the African slope (Sánchez-Garrido and Nadal 2022). The influence of sub-basin scale interactions between oceanic eddies and bottom topography in the western Mediterranean Sea was first investigated in Alvarez et al. (1994). By testing a parameterization of the eddy-topography interaction through various numerical experiments, their study suggested that topographic stress could play a significant role in driving circulation at both basin and sub-basin scales.

The Western Mediterranean Sea exhibits a high degree of eddy variability, manifested in its boundary currents, surface density fronts, pycnocline depth variations, large and rapid vertical parcel displacements, topographic vorticity generation, and abyssal currents, e.g., partly as shown in the recent CALYPSO experiments (Mahadevan et al. 2020).

### 3. Mediterranean simulations

We use two realistic simulations with the Regional Oceanic Modeling System (ROMS) (Shchepetkin and McWilliams 2005) and different spatial grid resolutions: (1) the Western Mediterranean basin (WMED) with horizontal $dx = 1500$ m and 60 $\sigma$-levels in the vertical (with stretching parameters $\theta_s = 7$, $\theta_b = 5$ near the surface and the bottom respectively, following Lemarié et al. (2012)); and (2) the nested subdomain for the Alboran Sea (ALB) with $dx = 500$ m, 140 $\sigma$-levels, and $\theta_s = \theta_b = 7$. WMED was calculated using CROCO, a modeling system built upon ROMS with a numerical kernel close to the UCLA version of the ROMS model (Shchepetkin and McWilliams 2009, http://www.croco-ocean.org). ALB was calculated using the UCLA-ROMS version of the model (https://github.org/CESR-lab/ucla-roms.git). The two codes have enough of a common heritage that their solution behaviors are similar.

ROMS is a split-explicit, 3D free-surface oceanic model that solves the hydrostatic Primitive Equations for momentum, temperature, and salinity using the Boussinesq approximation and
hydrostatic vertical momentum balance (Shchepetkin and McWilliams 2005). Horizontal advection is based on a modified version of the upstream 3rd-order scheme, which is split into a 4th-order centered advection and a rotated biharmonic diffusion with grid-dependent diffusivity satisfying the Peclet constraint (Marchesiello et al. 2009; Lemarié et al. 2012). To avoid numerical instabilities, the algorithm for vertical advection of momentum combines a high-order explicit scheme with an implicit one, which mitigates the vertical Courant-Friedrichs-Lewy restriction (Shchepetkin 2015). Vertical advection of tracers is performed using a 4th-order Akima scheme.

In both configurations the top and bottom boundary layers are well resolved using a K-Profile Parameterization (KPP) for the turbulent mixing (Large et al. 1994). With this scheme the mixing coefficients reach maximum values in the middle of the layers, decreasing to minimal background values within the interior, except where the Richardson number, $Ri$, becomes small. Bottom drag is parameterized using the quadratic law, $\tau = \rho_0 C_D ||u_b|| u_b$, where $\rho_0$ is the reference density, $u_b$ is the horizontal velocity in the bottom layer, and $C_D$ is a drag coefficient. We use the logarithmic formulation for $C_D = [k/\log(\Delta z_b / z_r)]^2$, where $k = 0.41$ is the von Karman constant, $\Delta z_b$ is the thickness of the bottom layer, and $z_r$ is the roughness scale, set to 1 cm.

Bathymetry for both domains is obtained from the SRTM30_PLUS_V11 global dataset (Becker et al. 2009). A steepness factor $r < 0.2$ is applied to realistic topography to limit pressure gradient errors. An additional smoothing with a radius of 4 km half-width is also applied to the grid in the ALB solution.

Realistic atmospheric forcing for the WMED solution is derived from the NCEP-CFSR data set at 0.2° resolution. Wind velocity, temperature and humidity fluxes, and short and long wave radiation fluxes feed the model at hourly time steps, and river inflow is derived from a climatological dataset. The initial and boundary conditions for WMED are interpolated from the Copernicus Marine Service global ocean forecasting system (Lellouche et al. 2018), providing a realistic representation of the Western Mediterranean Sea in the mesoscale range (Mason et al. 2019). This simulation was run for four years, 2010-2011, with a spin-up period of approximately three months, although the first year is excluded from the analysis. Capó (2020) provides a detailed description of this solution, including validation and data analysis.

The corresponding atmospheric forcing for ALB is obtained from the ECMWF atmospheric reanalysis of the global climate product (ERA5) at hourly time steps (Hersbach et al. 2020). The
initial state and boundary conditions for ALB are obtained from WMED, and the simulation is
initialized in December 2010, eleven months after WMED initialization to ensure sufficient spin-up
of the latter.

Boundary conditions consist of a Flather-type scheme for the barotropic mode, and Orlansky-type
for baroclinic velocity and for tracers, following the nesting procedure described in Mason et al.
(2010).

Neither of the simulations analized here is forced with tides. Although internal tides play a crucial
role in enhancing mixing over abyssal slopes in the oceans, these are not the sole contributors in our
study region and are not so essential for our analysis. Lee waves, topographic wakes, submesoscale
instabilities and interactions with mesoscale eddies also play important roles in enhancing mixing
along abyssal slopes (Legg 2021; Mashayek et al. 2024), and these processes occur in our solution
at submesoscale resolution despite the absence of tidal forcing. While the tidal effect is important
as part of the Atlantic inflow through the Strait of Gibraltar, tidal amplitudes decay towards the
interior of the basin, and the effect of topographic wakes, lee waves and other instability sources
can be as relevant as internal tides in enhancing mixing near abyssal slopes in the Mediterranean
Sea. Comparisons with tidal solutions show modest differences in the topostrophy and downhill
mean flows that are the primary focus in this paper (Sec. 4).

The ALB simulation is based on a previous solution with identical configuration that included
tidal forcing, ALB_T, and is extensively validated and analyzed in Capó (2020). The new ALB
solution presented here shows consistency with the parent WMED and with ALB_T. ALB validation
entailed comparison of the 3D flow structures, temperature and salinity fields, and the analysis of
time series of area averaged kinetic energy and local sea surface height throughout the domain.
The latter are consistent with those obtained from the parent WMED and serve to determine the
spin up period required until the model stabilizes. The absence of tidal forcing facilitates rapid
stabilization, which is achieved after two weeks. For greater reliability, the first two months of the
ALB simulation are excluded from the analyzed period.

The results presented in this paper, if not otherwise indicated, are obtained from daily-averaged
fields of the last 3 years of the WMED solution (2011–2013) and from the last year of ALB (Feb
2011-Jan 2012).
<figure><img src="image_3.png" /><figcaption>FIG. 2. Time-averaged $\overline{T}_b$ (a,b) and $\overline{w}_b$ (c,d) in the Western Mediterranean (a,c) and Alboran (b,d) simulations. Bathymetry is contoured in black every 500 m. The black box in (b) frames the Cartagena subdomain analyzed in Fig. 4 *et seq.*, and the dashed line indicates the section shown in Fig. 10. Locations of interest are labeled in (a) and (b).</figcaption></figure>

4. Abyssal topostrophy and downhill flow

In this section we examine the first two phenomena exposed in Sec. 1, topostrophy, $\overline{T}$, and downhill bottom flow, $\overline{w}_b$. Using our two simulations in the Western Mediterranean, WMED and ALB, and the definitions in (1)-(3). Figure 2 shows these fields over the whole domains of WMED and ALB. As anticipated, $\overline{T}_b$ is primarily positive, and $\overline{w}_b$ is primarily negative. Recall that both quantities are weighted by $\nabla h$ and thus tend to be largest over slopes. Not all slope regions have these signatures, but there is very little evidence of opposite signs with significant magnitude.

A more complete view is provided by the Probability Distribution Functions (PDFs) for these fields at locations with resting depths $\ge$ 500 m (Fig. 3). There are many small values due to locations with small $|\nabla h|$, but the PDFs are highly skewed in both their bulk and tails.
<figure><img src="image_4.png" /><figcaption>FIG. 3. Probability Density Functions for the bottom fields, $\bar{T}_b$ [m/s] (a,c) and $\bar{w}_b$ [m/s] (b,d), for time-averaged values at each location in the Western Mediterranean (a,b) and Alboran (c,d) simulations (in m s⁻¹). Only locations with resting depths $h \ge 500$ m are included. The predominance of negative $\bar{w}_b$ and positive $\bar{T}_b$ is evident, e.g., in the indicated values of skewness.</figcaption></figure>

The ALB simulation has a significantly increased spatial bathymetry and grid resolution compared to WMED (Sec. 3). It therefore exhibits greater variability in these measures, especially in $\bar{w}_b$. Partly, this is expressed in dipole patterns over bumps and ridges and along the sides of submarine canyons that interrupt the primary continental slopes. In particular, abundant submarine canyons interrupt the Spanish slope between 1.2°W - 0.2°E (the Cartagena subregion delimited in black in Fig. 2) and the African margin between 2°W - 1.5°E (the path of the eastward Alboran Current). It has been widely shown from both numerical simulations and observations that currents flowing across canyons induce opposite-signed $\bar{w}_b$ values on the upstream (negative $w$) and downstream (positive $w$) sides (Klinck 1996; Flexas et al. 2006; Spurgin and Allen 2014) with reversed $w$ signs for flow over bumps. In combination with prograde slope flow, this local canyon effect results in a net downhill transport. These features are shown in detail in Fig. 4 for
the Cartagena subregion. The prevailing topostrophic, downhill flow in this regions is part of the Levantine and Deep Mediterranean undercurrent that crosses the Alboran Sea westward towards the Strait of Gibraltar. We have chosen to focus on the Cartagena region, both here and below, as

<figure><img src="image_5.png" /><figcaption>FIG. 4. A zoom into the Cartagena region depicted within the black box in Fig. 2b for time-averaged $\bar{T}_b$ (a), $\bar{w}_b$ (b; colors), and $\bar{u}_b$ (b; vectors) for resting depths ≥ 500 m in the ALB simulation. Bathymetry is contoured in black every 100 m. The thick black line defines the sub-region for the along-isobath averages in Fig. 5 *et seq.*. The points C1 and C2 in (b) mark two representative locations on opposite sides along a ridge, used in Fig. 7 to illustrate the local flow vertical profiles.</figcaption></figure>

one that seems at least not atypical from the perspective of Fig. 2. No doubt there are meaningful
regional differences, but based on our present level of exploration, the behavior off Cartagena is representative of the structure and dynamics of abyssal slope flows.

To better discern the behavior associated with the primary abyssal slope, we average over the local flow structures evident in Fig. 4 by averaging in the $\hat{s}$ direction (i.e., along isobaths) within the thick black sub-region there. This is done by binning all grid points along the slope with a discretization interval of 40 m in the local depth.

The result is denoted by angle brackets, $\langle \cdot \rangle$, and we show cross-sections in vertical height $z$ and horizontal distance across the mean slope $\eta$ (increasing toward deeper water), e.g., $\langle \overline{T} \rangle(\eta, z)$ and $\langle \overline{w} \rangle(\eta, z)$ in Fig. 5. Topostrophy is positive for hundreds of meters above the slope, with maximum values at ~1200 m depth, and the averaged $\langle \overline{w} \rangle$ profile shows predominance of strong downwelling within the lowest ~80 m with weaker upwelling above and even weaker downwelling further into the interior.

Topostrophy and vertical velocities near the bottom have been also examined in a solution similar to ALB but forced with tides (Capó et al. 2021, 2023), with very similar results in the distribution of positive $\overline{T}$ and predominantly negative $\overline{w}_b$, with slightly increased values due to tidal effects on bottom velocities. Further ROMS simulations in the Atlantic and Pacific Oceans examined in Schubert et al. (2024) — with different topographic smoothing, horizontal resolution, and time averaging over several time scales — show the same overall predominance of prograde flows, with positive topostrophy and negative vertical velocity near the abyssal slopes. Their work includes diagnostics from a multitude of current-meter measurements that support the hypothesis of prevalent, widespread prograde currents and downhill flows.

## 5. Mean momentum

To investigate the various forces sustaining the topostrophic flow, in this section we analyze the momentum balance in our solution, with a focus on the Cartagena slope. Then we explore the eddy-driven nature of the widespread mean prograde flow by decomposing the advective forcing term into contributions from the mean flow, transient eddies, and standing eddies. Finally, an analogous decomposition of the time-averaged kinetic energy (KE) highlights the prominence of eddies as a principal characteristic.
<figure><img src="image_6.png" /><figcaption>FIG. 5. Time- and along-isobath-averaged cross-sections on the Cartagena slope (i.e., within the thick black line in Fig. 4, between 0.25-1.1°W) in the ALB simulation: (a) topostrophy $\langle \overline{\tau} \rangle (\eta, z) (\times 10^{-2})$ and (b) vertical velocity $\langle \overline{w} \rangle (\eta, z) (\times 10^{-3})$, both referenced to the color bar values on the right side. $\eta$ is a horizontal coordinate in the cross-isobath direction, increasing toward deeper water. Isopycnals are contoured in black in (a), and the inset box shows in detail the tilting of isopycnals perpendicular to the sloping floor, consistent with zero buoyancy flux at the boundary. Also in the inset, the dashed black line represents the mean bottom boundary-layer depth (as diagnosed in the model from KPP), and the local configuration of the sigma levels is contoured in gray, showing that the vertical resolution of the model increases significantly towards the bottom.</figcaption></figure>

a. Along-isobath momentum balance and bottom drag

The time-averaged momentum balance in the along-isobath direction is

$$ \hat{\mathbf{s}} \cdot \partial_t \bar{\mathbf{u}} \approx 0 = + \underbrace{\hat{\mathbf{s}} \cdot (\vec{u} \cdot \vec{\nabla}) \bar{\mathbf{u}}}_{\text{ADV}} - \underbrace{\hat{\mathbf{s}} \cdot \nabla \bar{\phi}}_{\text{PGF}} - \underbrace{f \hat{\mathbf{s}} \cdot \hat{\vec{z}} \times \bar{\mathbf{u}}}_{\text{COR}} + \underbrace{\hat{\mathbf{s}} \cdot \bar{\mathbf{F}}}_{\text{VMIX}}. \quad (8) $$
VMIX is the parameterized vertical mixing in ROMS described in Sec. 3. The other right-side tendency terms are advection of along-isobath momentum (ADV) and along-isobath pressure-gradient (PGF) and Coriolis (COR) forces. This balance is diagnosed on-line and then averaged in time and along-isobath. Additional terms in the discretized ROMS momentum equations are negligible in our time averaged momentum balance (8); these include numerical dissipation by the horizontal advection scheme and a rectification that results from the coupling of the 2D-3D barotropic-baroclinic, split-explicit time-stepping scheme (Shchepetkin and McWilliams 2005).

Figure 6 shows the momentum balance terms in various combinations. On the left panels, geostrophic balance is hidden in COR + PGF, leaving behind a residual that is mostly balanced by ADV, especially in the interior away from the slope. VMIX is small except within the BBL, where it reflects the retardation of the topostrophic current by bottom drag stress and its vertical redistribution by momentum mixing. Because the force balances exhibit continuity along the slope, the slope-averaged profiles shown in the right panels as a function of height above the bottom are representative of the balances in the BBL and adjacent interior. These balances broadly support the ideas of conservative dynamics in the interior and an Ekman-like balance in the BBL (but see Sec. 7).

b. Eddy-driven topostrophic flow

To explore the influence of the eddies in driving the topostrophic flow, we now decompose the advective term ADV into the contributions of the mean flow and eddy components. Given that the flow structure near the slope in this region is highly modulated by topography (Fig. 4), we chose to make a decomposition of the flow as follows: a flow component $u_i(\vec{x}, t)$ at a given location $\vec{x}$ and time $t$ can be written as

$$u_i(\vec{x}, t) = \langle \overline{u_i} \rangle + \tilde{u}_i + u'_i, \qquad (9)$$

where $\langle \overline{u_i} \rangle$ is the mean flow, defined here as the time- and along-isobath average; $\tilde{u}_i = \overline{u_i} - \langle \overline{u_i} \rangle$ is the standing (or stationary) eddies, and $u'_i = u_i - \overline{u_i}$ is the transient eddies. The standing eddy component accounts for the persistent circulation patterns that, in this particular case, are mainly caused by the interaction of the flow with topographic irregularities, such as canyons and bumps; the transient eddy component includes temporary eddies that are typically generated by instabilities
<figure><img src="image_7.png" /><figcaption>FIG. 6. Left panels show time- and along-isobath-averaged cross-sections on the Cartagena slope for the along-isobath momentum balance (8) in the ALB simulation: (a) $\langle \text{ADV} \rangle$, (b) $\langle \text{PGF} + \text{COR} \rangle$, (c) $\langle \text{ADV} + \text{PGF} + \text{COR} \rangle$, and (d) $\langle \text{VMIX} \rangle$. The 1D profiles on the right panels show the corresponding area-averaged terms for the lower 200 m, as a function of height above the bottom ($Hab = h+z$). While in the interior (100 m and above) the balance is between advection, Coriolis, and pressure-gradient, in the lower tens of meters (including the averaged BBL delimited by the green line), the accelerating combination of Coriolis, advection, and pressure-gradient is compensated by the retarding effect of bottom drag expressed in the vertical mixing term.</figcaption></figure>
in the flow field, such as ageostrophic shear instabilities and geostrophic barotropic or baroclinic instabilities.

Using (9), the time-averaged advective term in the momentum equation can be written as

$$ \overline{(\vec{u} \cdot \vec{\nabla})\mathbf{u}} = \underbrace{(\langle\langle \vec{u} \rangle\rangle \cdot \vec{\nabla})\langle\langle \mathbf{u} \rangle\rangle}_{\text{Mean}} + \underbrace{(\langle\langle \vec{u} \rangle\rangle \cdot \vec{\nabla})\tilde{\mathbf{u}} + (\tilde{\vec{u}} \cdot \vec{\nabla})\langle\langle \mathbf{u} \rangle\rangle + (\tilde{\vec{u}} \cdot \vec{\nabla})\tilde{\mathbf{u}}}_{\text{Standing eddies}} + \underbrace{(\tilde{\vec{u}}' \cdot \vec{\nabla})\mathbf{u}'}_{\text{Transient eddies}}. \quad (10) $$

The first term on the right side is the contribution of the mean flow, the next three terms represent the contribution of the standing eddies, and the last term represents the advection by transient eddies. By projecting (10) onto the along-isobath direction $\hat{s}$ we obtain the corresponding mean, standing eddy, and transient eddy contributions to the ADV term in (8).

Using an analogous decomposition, the total time-averaged KE is

$$ \frac{1}{2} \overline{\vec{u}^2} = \underbrace{\frac{1}{2} \langle \langle \vec{u} \rangle \rangle^2}_{\text{Mean}} + \underbrace{\tilde{\vec{u}} \cdot \langle \langle \vec{u} \rangle \rangle}_{\text{Standing eddies}} + \frac{1}{2} \tilde{\vec{u}}^2 + \underbrace{\frac{1}{2} \overline{\tilde{\vec{u}}'^2}}_{\text{Transient eddies}}, \quad (11) $$

where the first term on the right side is the KE of the mean flow; the next two terms constitute the eddy kinetic energy (EKE) stored in standing eddies; and the last term accounts for the EKE of the transient eddy flow.

The advection and KE results of this flow decomposition are summarized in Fig. 7. To highlight the relevant influence of topographic interaction in driving the abyssal sloping currents and to provide an idea of its vertical extent, in the first three panels we compare the time-averaged standing eddy component of the along-isobath $\tilde{U}_s$, the cross-isobath $\tilde{U}_n$, and the vertical $\tilde{w}$ components of the flow for three cases: the area-averaged flow (dark gray line), and at two representative locations at both sides of a bump, C1 and C2 (pink and blue dashed lines, respectively; see Fig. 4). In all situations, the along-isobath component is topostrophic, while the cross-isobath and vertical components are in opposite directions at both sides of the bump, and locally much more intense than the widespread downhill flow.

Figure 7 also shows that standing eddies dominate the mean advection of momentum in the lower $\sim 80$ m, while transient eddies are the dominant contribution above. Transient eddies dominate the mean KE at all depths near the slope. These results provide evidence of the eddy-driven nature
<figure><img src="image_8.png" /><figcaption>FIG. 7. The left three panels show time-averaged profiles of the along-slope ($\tilde{U}_s$), cross-slope ($\tilde{U}_n$), and vertical ($\tilde{w}$) velocities, both for the area-averaged flow in the Cartagena region (dark gray) and for two specific canyon locations shown in Fig. 4, C1 and C2 (pink and blue shadded lines, respectively). While the flow is generally topostrophic ($\tilde{U}_s > 0$, it is oriented uphill in the upstream side of a bump (C2, with $\tilde{U}_n < 0$ and $\tilde{w} > 0$) and downhill in the downstream side (C1, with $\tilde{U}_n > 0$ and $\tilde{w} < 0$). The right two panels show the time- and area-averaged profiles of the decomposed advective term in the momentum balance (ADV) and Kinetic energy (KE) separated into the total (black), mean (gray), standing eddy (dashed gray), and transient eddy (dotted gray) components. In the lower 100 m ADV is dominated by standing eddies, and by the transient eddy component above, while transient eddies are also the principal source of kinetic energy at all levels.</figcaption></figure>

of the topostrophic flow in our abyssal sloping currents. However, a more complete analysis of abyssal eddy dynamics, perhaps as a function of spatial scale, is needed to explain the widespread topostrophic outcome. Because the statistical stationarity of the solution is imperative in such analysis, we choose to defer this until we can achieve a higher-resolution simulation with a longer time-averaging period. In lieu thereof, the along-isobath averaged momentum balance that we present here is broadly in line with previous studies that show the turbulent character of the forces driving the mean along-isobath flow within the BBL (Molemaker et al. 2015; Chen et al. 2021; Stewart et al. 2024). Garabato et al. (2019) and Spingys et al. (2021) present observational evidence for the manner in which topographically localized turbulence is generated through submesoscale overturning instabilities in the near-boundary flow (e.g., symmetric and centrifugal instabilities;
Wenegrat and Thomas (2020)) that arise in the measured topostrophic and downhill mean slope flow. In Sec. 1 some further comments are made about present theoretical ideas about the cause of topostrophy. Such arguments are clearest in barotropic or isopycnal-layer models, whereas the ocean (and ROMS) is continuously stratified. Some attempts have been made to extend the theory to three dimensions (Merryfield and Holloway 1999; Venaile 2021), but not yet with much generality.

From a horizontal momentum-balance perspective, an important distinction is what averaging framework is chosen to encompass a region of complex topography. Beyond time averaging, one can choose between a mean horizontal direction (e.g., an average $\langle \hat{s} \rangle$ over the broad slope area) or a local-along slope direction ($\hat{s}$, as here). Stewart et al. (2024) analyzes idealized simulations of corrugated slope flows and finds that the near-bottom topostrophy-accelerating force is *PGF* with a mean along-shore average, whereas the same simulations show *ADV* acceleration using an along-isobath average, consistent with our finding here. In more realistic simulations of the California Undercurrent, another topostrophic slope flow, averaging more broadly in the along-isobath direction implicates a *PGF* acceleration cause (Molemaker et al. 2015; Chen et al. 2021). Because *ADV* and *PGF* are closely related quantities in incompressible fluid dynamics that manifest turbulent behaviors, there is no contradiction in these different diagnostic perspectives.

The choice of an averaged direction for the momentum balance might be simple for an unidirectional flow past an isolated bathymetric feature, but it can become ambiguous for complex currents and multi-scale topography as in these simulations. This led us to choose the local $\hat{s}$. Another important issue is whether the direction should be spatially smoothed; i.e., on what spatial scales is topostrophy more likely to occur? We have not yet investigated this, but notice that the general pattern of $\overline{T}$ is quite similar in the WMED and ALB simulations with different grid resolutions by a factor of three (Fig. 2).

# 6. Near-slope buoyancy balance and diapycnal mixing

The widespread negative $\overline{w}_b$ fields observed in the Western Mediterranean simulations, as well as in the Atlantic and Pacific solutions presented in Schubert et al. (2024), do not match the behavior expected by the 1D BBL model by Garrett (1990), which predicts up-slope mean current at the bottom ($\overline{w}_b > 0$) and negative topostrophy. In this section we analyze the mean BBL buoyancy
balance (7) for the ALB solution to show the vertical structure near the slope. The contradiction between our results and the predictions in Garrett (1990) is further examined in an alternative 1D BBL model in Sec. 7.

The buoyancy balance terms in (7) are obtained in ROMS from the online-diagnosed advection and mixing terms in the balances for potential temperature $T$ and salinity $S$. They are combined into buoyancy fluxes with a local linearization of the equation of state (EOS),

$$ \delta b = - \frac{g}{\rho_0} (-\alpha \delta T + \beta \delta S), \qquad (12) $$

where $\delta$ denotes the local differential change in the balance terms and $\alpha = \frac{\partial \rho}{\partial T}$, and $\beta = \frac{\partial \rho}{\partial S}$ are the local density expansion coefficients. For example, for advection,

$$ \overline{\vec{u} \cdot \vec{\nabla} b} = - \frac{g}{\rho_0} (-\alpha \overline{\vec{u} \cdot \vec{\nabla} T} + \beta \overline{\vec{u} \cdot \vec{\nabla} S}), \qquad (13) $$

with an analogous combination in $\overline{\vec{\nabla} \cdot \vec{\mathcal{B}}}$ for the mean $T$ and $S$ mixing terms. ROMS uses a nonlinear EOS that corresponds to the UNESCO formulation as derived by Jackett and McDougall (1995) that computes *in situ* density as a function of potential temperature, salinity and pressure, the latter approximated as a function of depth. This EOS is then “stiffened” by a split into an adiabatic and a compressible part to reduce pressure-gradient errors associated with nonlinear compressibility effect (Shchepetkin and McWilliams 2011).

The method of analyzing in ROMS the total buoyancy mixing term, $-\overline{\vec{\nabla} \cdot \vec{\mathcal{B}}}$, due both to the explicit subgrid-scale parameterization and to model discretization errors, is presented in Gula et al. (2024).

The occurrence of a mean downhill slope flow gives an expectation that the advective tendency in (7) will be positive, i.e., acting to increase the buoyancy. Of course, advection is three dimensional and time-varying, including by the mean current and by the standing and transient eddies. Assuming a nearly stationary state with negligible net buoyancy tendency, if the advective tendency is positive, then the necessary balancing mixing tendency will be negative, acting to decrease the buoyancy though a divergent $\vec{\mathcal{B}}$ flux. Because $\overline{\omega}$ in (5) has the same sign as the mixing tendency term, this is a wide layer with diapycnal downwelling along the slope, contrary to the arguments in Ferrari
et al. (2016) etc. This difference is explained through the use of an idealized boundary layer model
in Sec. 7.

Figure 8 shows the time- and along-isobath-averaged buoyancy balance for the Cartagena region, and a detail of the lowest 100 m over the slope is shown in Figure 9. As anticipated from $\overline{T}_b > 0$ and $\overline{w}_b < 0$ in Sec. 1, the buoyancy advection tendency is positive and the mixing tendency is negative in the middle and upper part of the boundary layer, with sign reversals in the adjacent interior.

<figure><img src="image_9.png" /><figcaption>FIG. 8. Time- and along-isobath-averaged cross-sections on the Cartagena slope for the buoyancy balance (7) in the ALB simulation: (a) advection and (b) mixing.</figcaption></figure>

However, in a thin layer several meters thick right at the slope, there is a negative advective
tendency, contrary to the $\overline{w}_b < 0$ expectation. This is a subtle effect of eddy buoyancy advection
as demonstrated below in Fig. 12 and further explained in Sec. 7. The expanded scale in the
BBL close to the slope in Fig. 9 shows that the vertical sign reversals occurs well inside the BBL,
accompanied by a further sign change higher in the interior. Here we are using a transformed
vertical coordinate, the height above the bottom, $Hab = z+h$.
<figure><img src="image_10.png" /></figure>

FIG. 9. A zoom of the fields in Fig. 8 in coordinates of height above the bottom $Hab = z+h$ and cross-slope horizontal distance $\eta$. The mean boundary-layer height (as diagnosed in the model from KPP) is the black thick line.

The positive sign of $\tilde{\omega}$ (5) in the outer part of the BBL and above means diapycnal upwelling. Its sign reversal is consistent with the sign reversal of the cross-isobath velocity profile in the 1D topostrophic boundary layer model (i.e., $V(Z)$ in Fig. 1 and $V'(Z)$ in Fig. 13 below). Thus, our results are consistent with the idea that much, maybe most, of the mixing-induced abyssal upelling in the global thermohaline circulation may be happening near the major topographic slopes, even if not is confined to a thin bottom boundary layer immediately above the slope.

While the tracer mixing in our simulations is prescribed by the turbulence parameterization (with the vertical eddy diffusivity $\kappa_v$ derived from the KPP scheme; Sec. 3), in practice mixing of temperature and salinity in the model can additionally arise from implicit contributions from the horizontal and vertical advection schemes as well as other discretization errors, in particular where currents and $T$, $S$ gradients are large. The resulting non-advection $T$, $S$ fluxes will result in an “effective” diffusivity, $\kappa_{eff}$, that can differ from the prescribed vertical diffusivity $\kappa_v$ (Gula et al. 2024) and that is not confined only to a thin bottom boundary layer immediately above the slope. In the context of abyssal slope currents, bottom-intensified mixing due to breaking inertia-gravity
waves and topographic wakes is an important phenomenon (Sec. 1) even beyond the turbulent mixing in the BBL. With the horizontal grid resolution of $dx = 500$ m in the ALB simulation (Sec. 3), we can expect an incomplete representation of this adjacent interior mixing apart from what occurs in $\kappa_v$ due to small Richardson number in the model-resolved shear and stratification (Large et al. 1994). Nevertheless, we can assess both it and the effective diapycnal mixing in our solution by evaluating the prescribed $\kappa_v$ with the diagnosed $\kappa_{eff}$. The diagnostic method for the latter is presented in Gula et al. (2024) and illustrated there with ROMS simulations for the Atlantic Ocean.

<figure><img src="image_11.png" /><figcaption>FIG. 10. The time-averaged eddy diffusivity $\kappa(y, z)$ [$m^2 s^{-1}$] in the ALB simulation along a section at 0.4°W, approximately between Oran and the coast near Cartagena (see Fig. 2b): (a) vertical $\kappa_v$ as diagnosed with the KPP parameterization and (b) diapycnal $\kappa_{eff}$ as diagnosed from the total diapycnal mixing (Gula et al. 2024). The black lines are adiabatic density surfaces.</figcaption></figure>
By definition, $\kappa_{eff}$ is obtained by projecting $\vec{\mathcal{B}}$ onto the direction perpendicular to neutral density surfaces, $\vec{\nabla}b/|\vec{\nabla}b|$, where $\vec{\nabla}b$ is the adiabatic buoyancy gradient:

$$ \kappa_{eff} = - \frac{\vec{\mathcal{B}} \cdot \vec{\nabla} b}{|\vec{\nabla} b|^2}. \qquad (14) $$

In the absence of an available online $\kappa_{eff}$ code in our ALB solution, we make an offline estimation by computing the net $\vec{\mathcal{B}}$ from the diagnosed T and S mixing fluxes using the locally linearized EOS as in (12) and the adiabatic buoyancy gradient via the stiffened nonlinear EOS introduced above. The profile in Fig. 10 across the Mediterranean Sea shows bottom-intensified mixing near the bottom reaching above the BBL height by hundreds of meters, especially over the slopes, with magnitudes much larger than in the interior, as well as elevated values in the surface boundary layer. We plot $|\kappa_{eff}|$ because there are occasional small negative values of $\kappa_{eff}$ due to the uncontrolled sign of model discretization errors. These distributions are at least qualitatively consistent with the abyssal mixing measurements (Waterhouse et al. 2014), from which interior values of $O(10^{-4} m^2 s^{-1})$ are estimated in the deep ocean, but increase to $O(10^{-2} m^2 s^{-1})$ near topography. Even if our model's mixing might be quantitatively inaccurate, we do expect that the simulated structure of the slope currents, buoyancy balance, and momentum balance (Sec. 5) to be robust, partly based on extensive studies of mixing sensitivity in 1D models (Sec. 7-8).

To show the typicality of the topostrophic regime, area-averaged profiles made over steep slopes in the whole of the abyssal Alboran Sea in the ALB simulation are in Fig. 11. These averages are conditioned on a minimum slope magnitude, $|\nabla h|$; as this minimum is decreased (from $log_{10}|\nabla h| > -0.8$ to $log_{10}|\nabla h| > -1.3$), the amplitudes of all the plotted quantities decrease, but the shapes remain nearly the same. (Because these are averages over many locations of quantities that are nonlinearly related to each other, these averaged profiles should not be viewed as relevant to any particular place.) Topostrophic currents are the norm, reaching well above the BBL; downhill flow occurs adjacent to the slope; and diapycnal downwelling ($\tilde{\omega} > 0$ in (5)) occurs in the bulk of the boundary layer. The cross-isobath mean flow turns uphill above the boundary layer along with diapycnal upwelling. As mentioned above, the mean buoyancy balance again reverses sign in a thin layer very close to the slope. The prevalence of net diapycnal downwelling within the BBL is reflected by the cumulative integral of the mixing term shown in (h), which is negative well into
the nearby interior before turning positive above. Finally, the effective diffusivity profiles reflect bottom-intensified mixing reaching far above the BBL.

<figure><img src="image_12.png" /><figcaption>FIG. 11. Area- and time-averaged, near-bottom vertical profiles for (left to right and top to bottom) (a) horizontal current in the along slope (topostrophic) direction; (b) horizontal current in the cross slope direction; (c) $\overline{w}$; (d) $\overline{T}$; (e) $log_{10}|\kappa_{eff}|$; (f) $-\vec{u} \cdot \vec{\nabla} b$; (g) $-\vec{\nabla} \cdot \vec{B}$ and (h) the cumulative vertical integral of $-\vec{\nabla} \cdot \vec{B}$, plotted against height above the bottom Hab. The averaging area is the whole of the ALB domain east of Gibraltar Strait and with steep slopes, applying different minimum thresholds for the slope steepness, $log_{10}|\nabla h|$ (as indicated by the line thicknesses). The green line is the mean BBL height.</figcaption></figure>

Figure 12a shows the decomposition of the advective tendency term, $-\vec{u} \cdot \vec{\nabla} b$, into the contributions of the time-mean flow, $-\vec{u} \cdot \vec{\nabla} \overline{b}$, and the eddies, $-\overline{\vec{u}' \cdot \vec{\nabla} b'}$. It demonstrates that the contribution
of the mean current (solid gray line) is always positive near the bottom, increasing the local buoyancy by advecting less dense water downhill, while the contribution of eddy advection (dashed gray line) is always negative. The net buoyancy advection is thus a competition between these two contributions. In the bulk of the boundary layer, the positive mean advection dominates. However, the negative eddy advection is the cause of the reversal of the mean buoyancy balance signs in the thin layer at the slope (Fig. 11e).

The source of this negative eddy buoyancy advection is indicated in Fig. 12b, which shows a vertically extensive profile of negative eddy buoyancy flux, $\overline{w'b'} < 0$, that is also the rate of conversion of eddy kinetic energy to eddy potential energy. This implicates mean-flow shear instability as the source of eddy energy followed by the conversion to potential energy. Furthermore, because the principal gradients of the mean currents on the broader scale of the abyssal slopes are vertical (Fig. 11a,b), it seems likely that it is a vertical shear instability phenomenon, i.e., an ageostrophic process on a scale larger than boundary-layer turbulence. The negative sign of $\overline{w'b'}$ is inconsistent with BBL baroclinic instability and thus is contrary to the circumstances analyzed in Callies (2018) and Drake et al. (2022).

The 1D profiles in Fig. 12 plotted vs. Hab show the area-averaged values for the Cartagena domain (solid) and for the whole Alboran Sea domain in regions with the sea-floor deeper than 800 m.

## 7. 1D topostrophy-mixing-slope BBL

To provide a framework for interpreting the flow structure and dynamical balances in the preceding sections, we redo the 1D BBL analysis following Garrett (1990) but insisting on a topostrophic interior flow and also including eddy buoyancy transport. The existence of the latter effects are demonstrated by the buoyancy and momentum analyses in Secs. 5-6, and they are due to dynamical processes that are not directly associated with the boundary layer turbulence.

The 1D model is for a steady flow and buoyancy profiles adjacent to a uniform slope. In the interior, the solution asymptotes to a uniform along-isobath horizontal velocity $U_\infty$ with uniform stratification $b_\infty(z) = N^2z$, both in geostrophic, hydrostatic balance with the geopotential field $\Phi_\infty(y, z) = -fU_\infty y + N^2z^2/2$.<fnref n="1" />

<footnote n="1">¹The background stratification is not necessarily a constant $N$ if $\kappa(z)$ is variable. Rather than develop a more general background balance model, we implicitly assume a background right-side maintenance term in the buoyancy balance, $B_0(z) = -\partial_z(\kappa N^2)$, to cancel any $N$ evolutionary</footnote>
<figure><img src="image_13.png" /><figcaption>FIG. 12. (a) Area- and time-averaged near-bottom vertical profiles for $-\overline{\vec{u} \cdot \vec{\nabla} b}$ (as shown in Fig. 11 for $log_{10}|\nabla h| > -1.3$; black); the corresponding input from the time-mean flow $-\overline{\vec{u} \cdot \vec{\nabla} b}$ (solid gray); and the eddy contribution $-\overline{\vec{u}' \cdot \vec{\nabla} b'}$ (dashed gray), plotted against height above the bottom Hab. The averaging area is the whole of the ALB domain, east of Gibraltar Strait, for $h \ge 800$ m. The green line is the mean BBL height. (b) Time- and along-isobath-averaged cross-sections on the Cartagena slope for the baroclinic term $\langle w' b' \rangle$, which entails a conversion from eddy available potential energy to eddy kinetic energy through baroclinic instabilities if the term is positive. The profiles in the inset box correspond to the area-averaged quantities plotted against Hab and for $h \ge 800$ m in the same region (Cartagena slope; solid) and for the whole ALB domain east of Gibraltar Strait (dashed).</figcaption></figure>

The model is expressed in rotated (X, Y, Z) coordinates to align with the sloping bottom at Z = 0. The coordinate transformation is

$$X = x, \quad Y = y \cos[\theta] + z \sin[\theta], \quad Z = -y \sin[\theta] + z \cos[\theta], \quad (15)$$

<footnote>tendency and to focus on the boundary layer dynamics. This is consistent with most previous 1D BBL analyses, and it might well implicate interior buoyancy advection and circulation, as in Callies (2018); Drake et al. (2022).</footnote>
where $\theta > 0$ is the angle of the sloping bottom upward to the north. In these coordinates,

$$ 
\begin{aligned}
b_{\infty} &= N^2(\sin[\theta]Y + \cos[\theta]Z) \\
\Phi_{\infty} &= -fU_{\infty}(\cos[\theta]Y - \sin[\theta]Z) + \frac{N^2}{2}(\sin[\theta]Y + \cos[\theta]Z)^2.
\end{aligned} 
\quad (16) 
$$

The velocity $\vec{u} = (U, V, W)$ is similarly rotated into this slope-aligned coordinate system. For $f > 0$ a topostrophic flow has $U_{\infty} < 0$. Notice that this configuration matches the sketch in Fig. 1; i.e., $\hat{\mathbf{s}} = -\hat{\mathbf{X}}$ for $f > 0$.

The 1D BBL problem matches this interior solution (with only along-isobath flow, i.e., $V_{\infty} = W_{\infty} = 0$) with boundary conditions at the slope with a no-slip velocity, $\vec{u} = 0$, and no buoyancy flux, $\kappa \partial_Z b = 0$, at $Z = 0$ for general, slope-normal eddy viscosity and diffusivity profiles, $v(Z)$ and $\kappa(Z)$. It seeks solutions that depend only on $Z$ apart from a uniform $Y$ dependency in $b$ such that $\partial_Y b = N^2 \sin[\theta]$, a constant, the same as for $\partial_Y b_{\infty}$. Thus, the ordinary differential equation (ODE) system for deviations from the interior fields (denoted by $b' = b - b_{\infty}$, etc.) is

$$ 
\begin{aligned}
-f \cos[\theta] V' &= \partial_Z (v \partial_Z U') + F' \\
f \cos[\theta] U' &= \sin[\theta] b' + \partial_Z (v \partial_Z V') \\
f \sin[\theta] U' &= -\partial_Z \Phi' + \cos[\theta] b' \\
N^2 \sin[\theta] V' &= \partial_Z (\kappa \partial_Z b') + B'_e.
\end{aligned} 
\quad (17) 
$$

The boundary conditions are

$$ 
\begin{aligned}
&U' = -U_{\infty}, \quad V' = 0, \quad \kappa \partial_Z b' = -N^2 \kappa \cos[\theta] \quad \text{at } Z = 0; \\
&U', V', b' \to 0 \quad \text{as } Z \to \infty.
\end{aligned} 
\quad (18) 
$$

These equations are linear in the dependent variables; the only advective term retained is the cross-isobath buoyancy advection, i.e., $V'$ times the constant $\partial_Y b$. $W' = 0$ everywhere. This system has forcing terms proportional to $F'(Z)$ and $B'_e(Z)$ in the interior of the ODE system and $U_{\infty}$ and $N^2 \cos[\theta]$ in the boundary conditions, but not all of these are independent as explained below.
$F'(Z)$ is the boundary layer deviation from the $ADV + PGF + COR = 0$ interior balance in the along-isobath mean momentum balance (Sec. 5). It represents the $ADV' + PGF'$ terms associated with the eddy advection or other process that causes topostrophy but is otherwise absent from (17); its role is to convey this effect into the BBL dynamics. $B'_e(Z)$ represents the eddy buoyancy advection, $-\langle \vec{u}' \cdot \vec{\nabla} b' \rangle$ in Fig 12a, also otherwise absent in (17). Both $F'$ and $B'_e \to 0$ as $Z \to \infty$.

Garrett (1990) shows how to solve this system when $F' = B'_e = 0$, partly by recasting it in terms of the cross-isobath streamfunction $\psi'(Z)$ with $V' = \partial_Z\psi'$. By elimination of other variables, the ODE system for $\psi'$ in a domain $[0, Z_{max}]$ is the following:

$$ 
\begin{aligned}
\partial_Z^2(\nu \partial_z^2 \psi') + 4\nu_v q^4 \psi' &= 4\nu q^4 \left( \kappa_0 \operatorname{ctn}[\theta] + \frac{1}{N^2 \sin[\theta]} \int_0^\infty B'_e dZ \right) \\
&\quad + \frac{f \cos[\theta]}{\nu} \int_Z^\infty F' dZ - \frac{\sin[\theta]}{\kappa} \int_Z^\infty B'_e dZ \\
\psi' &= \partial_Z \psi' = 0 \quad \text{at } Z = 0 \\
\psi' &\to \kappa_0 \operatorname{ctn}[\theta] + \frac{1}{N^2 \sin[\theta]} \int_0^\infty B'_e dZ, \quad \partial_Z \psi' \to 0 \quad \text{as } Z_{max} \to \infty, \quad (19)
\end{aligned}
$$

where $\kappa_0 = \kappa(0)$ and

$$ q^4 = \frac{1}{4} \left( \frac{f^2 \cos^2[\theta]}{\nu^2} + \frac{N^2 \sin^2[\theta]}{\nu \kappa} \right). \quad (20) $$

Given a solution for $\psi'(Z)$, we can evaluate the other profiles by

$$ 
\begin{aligned}
b'(Z) &= -N^2 \sin[\theta] \int_Z^\infty \left( \psi' - \kappa_0 \operatorname{ctn}[\theta] \right) - \frac{1}{N^2 \sin[\theta]} \int_0^\infty B'_e dZ'' \frac{dZ'}{\kappa(Z')} \\
&\quad - \int_Z^\infty \left( \int_{Z'}^\infty B'_e dZ'' \right) \frac{dZ'}{\kappa(Z')} \\
U'(Z) &= \frac{f \cos[\theta]}{\nu} \int_Z^\infty \left( \psi' - \kappa_0 \operatorname{ctn}[\theta] - \frac{1}{N^2 \sin[\theta]} \int_0^\infty B'_e dZ'' \right) \frac{dZ'}{\nu(Z')} \\
&\quad - \int_Z^\infty \left( \int_{Z'}^\infty F' dZ'' \right) \frac{dZ'}{\nu(Z')} \\
V'(Z) &= \partial_Z \psi',
\end{aligned}
\quad (21) $$

together with the bottom stress, $\tau = \nu \partial_Z(U', V')(0)$, and lateral transport anomaly, $T' = \int_0^\infty (U', V') dZ$. Notice that $U_\infty$ does not appear explicitly in (19)-(21); thus, the only forcing terms are $\kappa_0 \operatorname{ctn}[\theta]$, $F'$, and $B'_e$, while $U_\infty = -U'(0)$ is implicit.
Our purpose with this model is as a paradigm for explaining the topostrophic boundary layer structure evident in the ROMS simulations, rather than as a deliberate fit at any particular location. In this spirit, we choose simple boundary layer profile shapes for the forcing functions, e.g.,

$$F'(Z) = F_0 e^{-\ell Z}, \quad B'_e(Z) = B_{e0} e^{-\mu Z}, \qquad (22)$$

and “typical” values for the parameters to illustrate the solution behavior. With both analytic solutions (Appendix) and a numerical solver for (19)-(21), we have extensively explored different parameter values and different profile shapes for $v$, $\kappa$, $F'$, and $B'_e$. While these choices have meaningful quantitative consequences for the 1D solutions, their qualitative shapes are quite robust for positive diffusivities and mostly negative values for $F'$ and $B'_e$, consistent with the near-boundary mixing, momentum advection, and eddy buoyancy advection shown in Sec. 5-6. Therefore, we restrict our attention here to a representative set of solutions with constant $v$ and $\kappa$; this choice of parameter values is not delicate with respect to demonstrating the relevant behaviors. In our view, the value of the 1D model is as an idealized demonstration of the influence of topostrophy on the slope BBL — and as making a connection with the conventional 1D model of Garrett (1990) — rather than as a quantitative explanation of the simulation results.

The solution procedure is to decompose the solution into superimposable pieces,

$$\psi' = \psi'_{\kappa_0, B'} + \psi'_F, \qquad (23)$$

etc., where $\psi'_{\kappa_0, B'}$ satisfies the boundary-value problem (BVP) (19) for $F' = 0$ and $\psi'_F$ satisfies the BVP (19) for all the terms with $\kappa_0$ and $B'_e$ set to zero and $F' \neq 0$. The latter is initially specified with an arbitrary multiplicative coefficient that is then renormalized such that

$$U'_F(0) = -(U_\infty + U'_{\kappa_0, B'}(0)), \qquad (24)$$

where $U_\infty$ is specified. This renormalization factor is applied to all the $\cdot_F$ functions. The total solution is then evaluated by summation as in (23).
The 1D BBL solutions with representative parameter values are in Fig. 13 for three cases: (1) A mixing BBL where $F' = B'_e = 0$ and

$$U_{\infty \kappa_0} = \frac{f_{\kappa} \cos[\theta] \operatorname{ctn}[\theta]}{\nu q} > 0 \qquad (25)$$

is fully determined and directed opposite to topostrophic interior flow (cf., Garrett 1990); (2) A topostrophic BBL with momentum advection where $B'_e = 0$ and $U_\infty$ has the specified value by the choice of $F_0$; and (3) A topostrophic BBL with both momentum advection and eddy buoyancy advection where both $F'$ and $B'_e$ are nonzero.

Figure 13 shows the primary profiles for the three cases. As a broad summary statement, the topostrophic boundary layers are very different from the mixing boundary layer. $b'(z)$ is positive in all cases as a reflection of the requirement for the interior stratification to transition to a no buoyancy flux boundary condition. The along-isobath flow $U(z)$ monotonically increases from its no-slip zero value at the slope and has the opposite sign, by construction when $F' \neq 0$, due to the requirement that it match with the topostrophic interior value ($U_\infty < 0$ such that $\overline{T} > 0$). The cross-isobath flow $V'(Z)$ also is in the opposite direction near the boundary with topostrophy; i.e., it is downhill there, consistent with $\overline{w}_b < 0$. In all cases $V'$ reverses sign with distance above the boundary; its net lateral transport $T^{y'} = \psi'(\infty)$ is positive in the (1) mixing and (2) only-topostrophic cases but becomes net negative in case (3) with $B'_e \neq 0$. These topostrophic BBL solutions in cases (2) or (3) are the basis for the sketch in Fig. 1.

As $U_\infty$ becomes more negative (i.e., more topostrophic), the strength of the near-boundary downhill flow in $V'$ increases; e.g., for the parameters otherwise as in Fig. 13 for case (2), this $V' < 0$ regime occurs for all $U_\infty < -0.02 \text{ m s}^{-1}$. As $B_{e0}$ becomes more negative, the cross-slope transport $T^{y'}$ becomes more downhill; e.g., for the parameters otherwise as in Fig. 13 for case (3), the marginal value for this regime is $B_{e0} = -2.5 \times 10^{-10} \text{ m s}^{-3}$. I.e., the essential behaviors of downwhll flow, $V'(0) < 0$, and near-slope diapycnal downwelling are robust across wide variations in the 1D model with topostrophy, while the value and sign of the total cross-slope transport, $\psi'(\infty) = T^y$, is more variable depending on the relative importance of $B'_e$.

The signs of these near-slope $V'(Z)$ are consistent with Ekman-layer thinking in the sense that they are directed to the left of the along-isobath interior flows $U_\infty$. However, the sign-reversals in $V'(Z)$ are inconsistent with the relation between lateral Ekman transport $T'$ and the bottom stress $\tau$.
<figure><img src="image_14.png" /><figcaption>FIG. 13. Vertical profiles in the 1D BBL with constant diffusivities: from upper-left clockwise, cross-hill streamfunction $\psi'$, buoyancy anomaly $b'$, total along-isobath flow $U$, and uphill flow $V'$. The blue dashed lines are for the mixing BBL ($F' = B'_e = 0$); the blue solid lines are for topostrophic BBL with momentum advection ($F' \neq 0, B'_e = 0$); and the red solid lines are for the topostrophic BBL with both momentum advection and eddy buoyancy advection ($F', B'_e \neq 0$). For $b'$ the blue and red solid curves nearly coincide, and for $U$ they are also close. When $F' \neq 0$, $\ell = 0.09 \text{ m}^{-1}$, and $F_0 = -5.82 \times 10^{-6} \text{ or } -4.92 \times 10^{-6} \text{ m s}^{-2}$ for the cases without or with $B_{e0} \neq 0$, respectively, to assure the value of $U_\infty = -0.1 \text{ m s}^{-1}$. When $B'_e \neq 0$, $\mu = 0.05 \text{ m}^{-1}$ and $B_{e0} = -3 \times 10^{10}$. Other parameters are $f = 10^{-4} \text{ s}^{-1}$, $\kappa = \nu = 0.5 \times 10^{-2} \text{ m}^2 \text{ s}^{-1}$, $\theta = 0.1$, $N^2 = 10^{-6} \text{ s}^{-2}$, and $q = 0.12 \text{ m}^{-1}$.</figcaption></figure>

(both defined following (21)). In none of these three cases are either of the lateral Ekman transport relations satisfied. In particular, after vertically integrating the X-momentum equation in (17) and
using the $\partial_Z U'(0)$ boundary condition in (18), the cross-isobath transport relation is

$$ T^{y'} - \frac{\tau^x}{f \cos[\theta]} = - \frac{1}{f \cos[\theta]} \int_0^\infty F' dz, \qquad (26) $$

and the nonzero right-side momentum advection breaks the Ekman relation. (The along-isobath transport relation has a more complicated formula, not shown here.) In the limit of a flat bottom ($\theta \to 0$) and no topostrophic or eddy buoyancy advection forcing ($F' = B'_e = 0$), these relations revert to Ekman-layer transport.

<figure><img src="image_15.png" /><figcaption>FIG. 14. The 1D BBL buoyancy balance terms: mean plus eddy advection (solid) and slope-normal diffusion (dashed). From left to right the cases are the mixing BBL ($F' = B'_e = 0$), the topostrophic BBL with momentum advection ($F' \neq 0, B'_e = 0$), and the topostrophic BBL with both momentum advection and eddy buoyancy advection ($F', B'_e \neq 0$). The parameters are the same as in Fig. 13.</figcaption></figure>

The eddy buoyancy advection $B'_e$ has a moderate quantitative influence on the topostrophic BBL profiles, slightly so for $b'$ and $U$ and more so for $\psi'$ and $V'$, but not a qualitatively large one (Fig. 13), except perhaps for the sign of $T^{y'}$ that is determined by a competition between buoyancy mixing and eddy advection. Where the role of $B'_e$ becomes evident is in comparing cases (2) and (3) for their mean buoyancy balance in (17) (Fig. 14) between net buoyancy advection, $-N^2 \sin[\theta]V' + B'_e$, and turbulent buoyancy mixing, $\partial_Z(\kappa \partial_Z b')$. In the mixing-only case (1), advection tends to decrease the buoyancy as an upslope advection, while mixing tends to increase it (i.e., diapycnal upwelling).
In contrast, in the topostrophic cases (2) and (3) the lower part of the boundary layer has advective buoyancy increase and diapycnal downwelling. However, in a very thin layer near the slope, case (3) shows a small amount of diapycnal upwelling due to $B'_e \neq 0$ (cf., Fig. 11e).

The 1D topostrophic BBL model is therefore a useful paradigm for abyssal slope currents. Here we have only shown representative solutions. The parameter values and profile shapes for $\kappa, \nu, F'$ and $B'_e$ could be varied to explore their sensitivities, but it seems to us much more important to look to further measurements and simulation models like ours to better determine what the true distributions are for diapycnal mixing, topostrophic momentum mixing, and eddy buoyancy fluxes near the abyssal slopes.

The mixing-only 1D BBL model has been very influential in ideas about abyssal upwelling near slopes (Sec. 1), and the widespread prevalence of topostrophy is an important corrective to such arguments, as further discussed in Sec. 8. Topostrophy is also an important correction to the “arrested Ekman layer” idea that initial interior flows on slopes evolve toward the mixing-only steady state with a weak $U_\infty$ (e.g., Garrett et al. 1993; Brink and Lentz 2010; Ruan et al. 2019, 2021) when there are no confounding interior dynamics to counter this evolution. In 3D flows with sustained interior flows, the tendency toward arrest is weaker (Jagannathan et al. 2023), in particular in topostrophic situations.

We can contrast this 1D steady, stratified, slope BBL model with the conventional 1D steady Ekman layer model for unstratified flow above a flat bottom, both with constant eddy mixing coefficients. The connection between these simple models is not direct; i.e., in the limit of $\theta \to 0$ with other parameters held constant, the former has $(U_{\kappa\infty}, V') \to \infty$ because (17-18) imply non-zero stratification and lateral buoyancy advection (see Appendix formulas). In other words, the conventional Ekman BBL model does not encompass stratification, which is an essential part of the slope BBL models. To properly span these two regimes, a more general BBL model is required.

# 8. Discussion

This is an era of increasing appreciation of the dynamical importance of the ocean's abyssal topographic slopes: boundary-layer turbulence, inertia-gravity wave and vorticity generation, adjacent interior diapycnal mixing due to breaking waves and wakes, along-isobath (topostrophic) flows, and eddy fluxes.
The primary contributions in this paper are demonstrating widespread topostrophy and how this feature and momentum and buoyancy advection and turbulent mixing in the near-slope region lead to a down-slope mean flow and extensive diapycnal downwelling, taking a coarse-grained view beyond small-scale canyons and ridges. A similar result is reported in Schubert et al. (2024) for an analogous simulation in the Atlantic Ocean. Relative to the scientific literature, in our view, this amplifies the importance of eddy-driven topostrophy, and it provides an important modification to the recent advocacy of slope-adjacent diapycnal upwelling (Ferrari et al. 2016; de Lavergne et al. 2016; McDougall and Ferrari 2017; Callies and Ferrari 2018; de Lavergne et al. 2016; Holmes et al. 2018, 2019; Drake et al. 2020; Holmes and McDougall 2020; Drake et al. 2022; Peterson and Callies 2022; Baker et al. 2023; Peterson and Callies 2023; Wynne-Cattanach et al. 2024).

A strong rationale for the slope-upwelling view has been the 1D slope BBL model (e.g., Garrett 1990) and its elaborations with various depictions of the vertical profile of bottom-enhanced mixing. We show how this is modified by a topostrophic flow in Sec. 7 to reverse the slope-adjacent flow, but still have diapycnal upwelling in the outer part of the BBL and nearby interior. The global thermohaline circulation requires net diabatic and volumetric upwelling in the abyss, of course, and because of the locally enhanced diabatic mixing it seems likely that much if not most of it occurs in topographic neighborhoods. However, regional models, such as the ones analyzed here and in Schubert et al. (2024), are not the right tool for assessing the net global upwelling because they are strongly influenced by imposed abyssal inflows through their open lateral boundaries.

# 9. Summary and Prospects

Abyssal slope currents are significantly topostrophic, due, in a general way, to turbulent “mesoscale” oceanic momentum dynamics. From this phenomenon follow both the downhill mean flow and the diapycnal slope downwelling and nearby upwelling that are also widespread phenomena, supported by both boundary-layer turbulence and adjacent diapycnal mixing.

This paper and its companion (Schubert et al. 2024) comprise a synthesis and partial revision of several prevailing ideas about oceanic circulation and mixing: geostrophic topographic turbulence, slope bottom boundary layers, near-slope mean-flow eddy instability, and abyssal upwelling due to bottom-intensified diapycnal mixing. However, current-topography interaction in the abyss remains a significant scientific frontier. Apart from simple idealized configurations, it involves
respectively. For each of these, $U_\infty$ is fully determined. For $U'_F$ we will take the inverse view that $F'(Z)$ is chosen such that $U_\infty$ has the desired topographic value.

The solutions specified below are plotted in Fig. 13 for the parameter values listed in the caption.

### a. Mixing BBL

First consider the case with $F' = 0$ (referred to as case (1) in Sec. 7). Its solution with constant $v$ and $\kappa$ is presented in Garrett (1990). The formulas are

$$ \psi'_{\kappa}(Z) = \kappa \operatorname{ctn}[\theta] (1 - e^{-qz}(\sin[qz] + \cos[qz])) . \quad (A2) $$

$$ \begin{align*} b'_{\kappa}(Z) &= \frac{N^2 \cos[\theta]}{q} e^{-qz} \cos[qz] \\ U'_{\kappa}(Z) &= -U_* e^{-qz} \cos[qz], \\ V'_{\kappa}(Z) &= 2q\kappa \operatorname{ctn}[\theta] e^{-qz} \sin[qz], \end{align*} \quad (A3) $$

where the surface value of the along-isobath velocity anomaly is $U'_{\kappa}(0) = -U_*$ with

$$ U_* = \frac{f\kappa \cos[\theta] \operatorname{ctn}[\theta]}{\nu q} > 0. \quad (A4) $$

As anticipated, this solution is consistent only with a retrograde interior along-isobath velocity, i.e., the opposite of topostrophy.

### b. Topostrophic BBL

Now consider the solution with $F' \neq 0$ with no terms explicitly proportional to $\kappa$ in (19)-(21), because they are accounted for in the mixing BBL solution, but still with a $\kappa$ dependency in both $q$ and the denominator of the coefficient for $b'(Z)$. For the 1D BBL problem, $F'$ is not known a priori because it is considered to be a product of the 3D interior topostrophic dynamics. This product is conveyed into the BBL problem by requiring that the total $U(Z) \to U_\infty$ as $Z \to \infty$ with
a specified $U_\infty$ value. For the $U'_F$ surface boundary condition from (18), this implies that

$$U'_F(0) = -(U_\infty - U_*). \qquad (A5)$$

Then $F'(Z)$ is chosen such that this condition is satisfied. Thus, we can specify a simple representative shape,

$$F'(Z) = F_0 e^{-\ell Z}, \qquad (A6)$$

whose amplitude $F_0$ is chosen to satisfy (A5). The analytic solution for $\ell \neq q$ is

$$
\begin{align*}
\psi'_{F}(Z) &= \psi_{0} \left( e^{-\ell Z} - \cos[qZ] e^{-qZ} + \left(\frac{\ell}{q} - 1\right) \sin[qZ] e^{-qZ} \right) \\
b'_{F}(Z) &= -\frac{N^2 \sin[\theta] \psi_0}{\ell \kappa} \left( e^{-\ell Z} - \frac{\ell}{q} \cos[qZ] e^{-qZ} + \frac{\ell^2}{2q^2} (\sin[qZ] + \cos[qZ]) e^{-qZ} \right) \\
U'_{F}(Z) &= \left( \frac{f \cos[\theta] \psi_0}{\ell \nu} - \frac{F_0}{\ell^2 \nu} \right) e^{-\ell Z} - \frac{f \cos[\theta] \psi_0}{q \nu} \cos[qZ] e^{-qZ} \\
&\quad + \frac{f \ell \cos[\theta] \psi_0}{2q^2 \nu} (\sin[qZ] + \cos[qZ]) e^{-qZ} \\
V'_{F}(Z) &= \psi_0 \left( -\ell e^{-\ell Z} + 2q \sin[qZ] e^{-qZ} + \ell (\cos[qZ] - \sin[qZ]) e^{-qZ} \right), \tag{A7}
\end{align*}
$$

where

$$ \psi_0 = \frac{f \cos[\theta]}{\ell \nu^2 (\ell^4 + 4q^4)} F_0 \quad \text{and} \\ F_0 = (U_\infty - U_*) / \left( \frac{1}{\ell^2 \nu} - \frac{f^2 \cos^2[\theta]}{\ell^2 \nu^3} \left( \frac{1 - \ell/q + \ell^2/(2q^2)}{\ell^4 + 4q^4} \right) \right). \qquad (A8) $$

With the parameter values listed in the caption to Fig. 13, $F_0 = -5.8 \times 10^{-6} \text{ m s}^{-2}$ and $\psi_0 = -0.80 \text{ m}^2 \text{ s}^{-1}$ for the case referred to as (2) in Sec. 7.

### References

[779] Allen, J., D. Smeed, J. Tintoré, and S. Ruiz, 2001: Mesoscale subduction at the Almeria–Oran front: Part 1: Ageostrophic flow. *J. Mar. Syst.*, **30**, 263–285.

[781] Alvarez, A., J. Tintoré, G. Holloway, M. Eby, and J. M. Beckers, 1994: Effect of topographic stress on circulation in the western Mediterranean. *J. Geophys. Res. Oceans*, **99**, 16053–16064.
783 Baker, L., A. Mashayek, and A. Naveira Garabato, 2023: Boundary upwelling of Antarctic Bottom Water by topographic turbulence. *AGU Advances*, **4**, e2022AV000858.

784 Becker, J. J., and Coauthors, 2009: Global bathymetry and elevation data at 30 arc seconds resolution: SRTM30_PLUS. *Marine Geodesy*, **32**, 355–371.

787 Bretherton, F., and D. Haidvogel, 1976: Two-dimensional turbulence above topography. *J. Fluid Mech.*, **78**, 129–154.

789 Brink, K. H., and S. J. Lentz, 2010: Buoyancy arrest and bottom Ekman layer transport. Part I: Steady flow. *J. Phys. Ocean.*, **40**, 621–635.

791 Callies, J., 2018: Restratification of abyssal mixing layers by submesoscale baroclinic eddies. *J. Phys. Ocean.*, **48**, 1995–2010.

793 Callies, J., and R. Ferrari, 2018: Dynamics of an abyssal circulation driven by bottom-intensified mixing on slope. *J. Phys. Ocean.*, **48**, 1257–1282.

795 Capó, E., 2020: Submesoscale dynamics in the western Mediterranean Sea. Ph.D. thesis, University of the Balearic Islands, http://hdl.handle.net/11201/151819.

797 Capó, E., 2023: ALB500 ROMS simulation [Dataset]. 10.5281/zenodo.10215065, URL https://doi.org/10.5281/zenodo.10215065.

799 Capó, E., and J. McWilliams, 2022: Coherent Lagrangian pathways near an East Alboran front. *J. Geophys. Res. Oceans*, **127**, e2021JC018022.

801 Capó, E., J. McWilliams, and A. Jagannathan, 2023: Topographic vorticity generation along the Spanish coast in the Alboran Sea. *J. Geophys. Res. Oceans*, **128**, e2022JC019480.

803 Capó, E., J. McWilliams, E. Mason, and A. Orfila, 2021: Intermittent frontogenesis in the Alboran Sea. *J. Phys. Ocean.*, **51**, 1417–1439.

805 Chen, R., J. McWilliams, and L. Renault, 2021: Momentum governors of the California Undercurrent. *J. Phys. Ocean.*, **51**, 2915–2932.

807 de Lavergne, C., G. Madec, J. L. Sommer, A. Nurser, and A. N. Garabato, 2016: The impact of a variable mixing efficiency on the abyssal overturning. *J. Phys. Ocean.*, **46**, 663–681.
809

Drake, H., R. Ferrari, and J. Callies, 2020: Abyssal circulation driven by near-boundary mixing:

810

Water mass transformations and interior stratification. *J. Phys. Ocean.*, **50**, 2203–2226.

811

Drake, H., X. Ruan, J. Callies, K. Ogden, A. Thurnherr, and R. Ferrari, 2022: Dynamics of eddying

812

abyssal mixing layers over sloping rough topography. *J. Phys. Ocean.*, **52**, 3199–3219.

813

Ferrari, R., T. McDougall, M. Nikurashin, and J. Campin, 2016: Turning ocean mixing upside

814

down. *J. Phys. Ocean.*, **46**, 2239–2261.

815

Flexas, M., D. Gomis, S. Ruiz, A. Pascual, and P. León, 2006: In situ and satellite observations of

816

the eastward migration of the Western Alboran Sea Gyre. *Progr. Oceanogr.*, **70**, 486–509.

817

Garabato, A. N., K. Polzin, B. King, K. Heywood, and M. Visbeck, 2004: Widespread intense

818

turbulent mixing in the Southern Ocean. *Science*, **303**, 210–213.

819

Garabato, A. N., and Coauthors, 2019: Rapid mixing and exchange of deep-ocean waters in an

820

abyssal boundary current. *Proc. Nat. Acad. Sci.*, **116**, 13 233–13 238.

821

Garrett, C., 1990: The role of secondary circulation in boundary mixing. *J. Geophys. Res.*, **95**,

822

3181–3188.

823

Garrett, C., P. MacCready, and P. Rhines, 1993: Boundary mixing and arrested Ekman layers:

824

Rotating stratified flow near a sloping boundary. *Ann. Rev. Fluid Mech.*, **25**, 291–323.

825

Gula, J., M. Molemaker, P. Damien, E. Capó, R. Schubert, F. Lemarié, N. Schifano, and

826

J.C.McWilliams, 2024: Quantifying diapycnal mixing in a sigma coordinate model. In prepara-

827

tion.

828

Herring, J., 1977: On the statistical theory of two-dimensional topographic turbulence. *J. Atmos.*

829

Sci., **34**, 1731–1750.

830

Hersbach, H., and Coauthors, 2020: The ERA5 global reanalysis. *Q. J. R. Meteorol. Soc.*, **146**,

831

1999–2049.

832

Holloway, G., 1987: Systematic forcing of large-scale geophysical flows by eddy-topography

833

interaction. *J. Fluid Mech.*, **184**, 463–478.
834 Holloway, G., 1996: Neptune effect: statistical mechanical forcing of ocean circulation. *Stochastic Modelling in Physical Oceanography*, R. Adler, P. Muller, and B. Rozovski, Eds., Birkhauser, 207–219.

837 Holloway, G., 2008: Observing global ocean topostrophy. *J. Geophys. Res. Oceans*, **113**, C07054.

838 Holmes, R., C. de Lavergne, and T. McDougall, 2018: Ridges, seamounts, troughs, and bowls: Topographic control of the dian Neutral circulation in the abyssal ocean. *J. Phys. Ocean.*, **48**, 861–882.

841 Holmes, R., C. de Lavergne, and T. McDougall, 2019: Tracer transport within abyssal mixing layers. *J. Phys. Ocean.*, **49**, 2669–2695.

843 Holmes, R., and T. McDougall, 2020: Diapycnal transport near a sloping bottom boundary. *J. Phys. Ocean.*, **50**, 3253–3266.

845 Jackett, D. R., and T. J. McDougall, 1995: Minimal adjustment of hydrostatic profiles to achieve static stability. *J. Atmos. Ocean. Tech.*, **12**, 381–389.

847 Jagannathan, A., K. Srinivasan, J. McWilliams, M. Molemaker, and A. Stewart, 2023: Evolution of bottom boundary layers on three-dimensional topography: Ekman arrest and instabilities. *J. Geophys. Res. Oceans*, **128**, e2023JC019705.

850 Klinck, J., 1996: Circulation near submarine canyons: A modeling study. *J. Geophys. Res. Oceans*, **101**, 1211–1223.

852 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.

854 Legg, S., 2021: Mixing by oceanic lee waves. *Annu. Rev. Fluid Mech.*, **53**, 173–201.

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

858 Lemarié, F., J. Kurian, A. F. Shchepetkin, M. J. Molemaker, F. Colas, and J. C. McWilliams, 2012: Are there inescapable issues prohibiting the use of terrain-following coordinates in climate models? *Ocean Modelling*, **42**, 57–79.
861

Mahadevan, A., A. Pascual, D. L. Rudnick, S. Ruiz, J. Tintoré, and E. D'Asaro, 2020: Coherent pathways for vertical transport from the surface ocean to interior. *Bull. Am. Meteorol. Soc.*, **101**, E1996–E2004.

862

863

864

Marchesiello, P., L. Debreu, and X. Couvelard, 2009: Spurious diapycnal mixing in terrain- following coordinate models: The problem and a solution. *Ocean Modelling*, **26**, 156–169.

865

866

Mashayek, A., J. Gula, L. E. Baker, A. C. Naveira Garabato, L. Cimoli, J. J. Riley, and C. de Lavergne, 2024: On the role of seamounts in upwelling deep-ocean waters through turbulent mixing. *Proc. Natl. Acad. Sci.*, **121**, e2322163 121.

867

868

869

Mason, E., M. J. Molemaker, A. Shchepetkin, F. Colas, J. C. McWilliams, and P. Sangrà, 2010: Procedures for offline grid nesting in regional ocean models. *Ocean Modelling*, **35**, 1–15.

870

871

Mason, E., S. Ruiz, R. Bourdalle-Badie, G. Reffray, M. García-Sotillo, and A. Pascual, 2019: New insight into 3-D mesoscale eddy properties from CMEMS operational models in the western Mediterranean. *Ocean Sci.*, **15**, 1111–1131.

872

873

874

McDougall, T., and R. Ferrari, 2017: Abyssal upwelling and downwelling driven by near-boundary mixing. *J. Phys. Ocean.*, **47**, 261–283.

875

876

Merryfield, W., and G. Holloway, 1999: Eddy fluxes and topography in stratified quasi-geostrophic models. *J. Fluid Mech.*, **380**, 59–80.

877

878

Molemaker, M. J., J. C. McWilliams, and W. K. Dewar, 2015: Submesoscale instability and generation of mesoscale anticyclones near a separation of the California Undercurrent. *J. Phys. Ocean.*, **45**, 613–629.

879

880

881

Munk, W., 1966: Abyssal recipes. *Deep Sea Research*, **13**, 707–730.

882

Peterson, H. G., and J. Callies, 2022: Rapid spinup and spindown of flow along slopes. *J. Phys. Ocean.*, **52**, 579–596.

883

884

Peterson, H. G., and J. Callies, 2023: Coupling between abyssal boundary layers and the interior ocean in the absence of along-slope variations. *J. Phys. Ocean.*, **53**, 307–322.

885

886

Phillips, O., 1970: On flows induced by diffusion in a stably stratified fluid. *Deep-Sea Res.*, **17**, 435–443.

887
888

Rhines, P., and W. Holland, 1979: A theoretical discussion of eddy-driven mean flows. *Dyn. Atmos. Oceans*, **3**, 283–325.

889

890

Ruan, X., A. Thompson, and J. Taylor, 2019: The evolution and arrest of a turbulent stratified oceanic bottom boundary layer over a slope: Downslope regime. *J. Phys. Ocean.*, **49**, 469–487.

891

892

Ruan, X., A. Thompson, and J. Taylor, 2021: The evolution and arrest of a turbulent stratified oceanic bottom boundary layer over a slope: Upslope regime and PV dynamics. *J. Phys. Ocean.*, **51**, 1077–1089.

893

894

895

Salmon, R., G. Holloway, and M. Hendershott, 1976: The equilibrium statistical mechanics of simple quasi-geostrophic models. *J. Fluid Mech.*, **75**, 691–703.

896

897

Sánchez-Garrido, J. C., and I. Nadal, 2022: The Alboran Sea circulation and its biological response: A review. *Front. Mar. Sci.*, **9**, 393–390.

898

899

Saunders, P., and B. King, 1995: Bottom currents derived from a shipborne ADCP on WOCE Cruise A11 in the South Atlantic. *J. Phys. Ocean.*, **25**, 329–347.

900

901

Schubert, R., J. Gula, E. Capó, P. Damien, M. Molemaker, C. Vic, and J. McWilliams, 2024: The ocean moves downhill near the seafloor and recirculates uphill above. *Nature*, under review.

902

903

Shchepetkin, A. F., 2015: An adaptive, Courant-number-dependent implicit scheme for vertical advection in oceanic modeling. *Ocean Modelling*, **91**, 38–69.

904

905

Shchepetkin, A. F., and J. C. McWilliams, 2005: The Regional Oceanic Modeling System (ROMS): A split-explicit, free-surface, topography-following-coordinate oceanic model. *Ocean Modelling*, **9**, 347–404.

906

907

908

Shchepetkin, A. F., and J. C. McWilliams, 2009: Computational kernel algorithms for fine-scale, multi-process, long-term oceanic simulations. *Handbook of Numerical Analysis: Computational Methods for the Ocean and the Atmosphere*, R. Temam, and J. Tribia, Eds., Vol. 14, Elsevier, 121–183.

909

910

911

912

Shchepetkin, A. F., and J. C. McWilliams, 2011: Accurate Boussinesq oceanic modeling with a practical, “stiffened” equation of state. *Ocean Modelling*, **38**, 41–70.

913
Spingys, C. P., A. C. N. Garabato, S. Legg, K. L. Polzin, E. P. Abrahamsen, C. E. Buckingham, A. Forryan, and E. E. Frajka-Williams, 2021: Mixing and transformation in a deep western boundary current: A case study. *J. Phys. Ocean.*, **51**, 1205 – 1222, URL https://journals.ametsoc.org/view/journals/phoc/51/4/JPO-D-20-0132.1.xml, 10.1175/JPO-D-20-0132.1.

Spurgin, J., and S. Allen, 2014: Flow dynamics around downwelling submarine canyons. *Ocean. Sci.*, **10**, 799–819.

St. Laurent, L., J. Toole, and R. Schmitt, 2001: Buoyancy forcing by turbulence above rough topography in the abyssal Brazil Basin. *J. Phys. Ocean.*, **31**, 3476–3495.

Stewart, A., Y. Wang, A. Solodoch, R. Chen, and J. McWilliams, 2024: Formation of eastern boundary undercurrents via mesoscale eddy rectification. *J. Phys. Ocean.*, submitted.

Stommel, H., 1958: The abyssal circulation. *Deep-Sea Res.*, **5**, 80–82.

Thompson, L., and G. C. Johnson, 1996: Abyssal currents generated by diffusion and geothermal heating over rises. *Deep-Sea Res. I*, **43**, 193–211.

Tintore, J., P. L. Violette, I. Blade, and A. Cruzado, 1988: A study of an intense density front in the eastern alboran sea: The almeria-oran front. *J. Phys. Ocean.*, **18**, 1384–1397.

Vargas-Yáñez, M., F. Plaza, J. García-Lafuente, T. Sarhan, J. Vargas, and P. Vélez-Belchí, 2002: About the seasonal variability of the Alboran Sea circulation. *J. Mar. Syst.*, **35**, 229–248.

Venaile, A., 2021: Bottom-trapped currents as statistical equilibrium states above topographic anomalies. *J. Fluid Mech.*, **699**, 500–510.

Viúdez, Á., J. Tintoré, and R. L. Haney, 1996: Circulation in the Alboran Sea as determined by quasi-synoptic hydrographic observations. Part I: Three-dimensional structure of the two anticyclonic gyres. *J. Phys. Ocean.*, **26**, 684–705.

Waterhouse, A., and Coauthors, 2014: Global patterns of diapycnal mixing from measurements of the turbulent dissipation rate. *J. Phys. Ocean.*, **44**, 1854–1872.

Weijer, W., A. Barthel, M. Veneziani, and H. Steiner, 2020: The Zapiola Anticyclone: A Lagrangian study of its kinematics in an eddy-permitting ocean model. *Deep-Sea Res. I*, **164**, 103 308.
Wenegrat, J. O., and L. N. Thomas, 2020: Centrifugal and symmetric instability during Ekman adjustment of the bottom boundary layer. *J. Phys. Ocean.*, **50**, 1793–1812.

Wunsch, C., 1970: On oceanic boundary mixing. *Deep-Sea Res.*, **17**, 293–301.

Wynne-Cattanach, B. L., and Coauthors, 2024: Observations of diapycnal upwelling within a sloping submarine canyon. *Nature*, **630**, 884–890.

Xie, X., and Coauthors, 2022: Enhanced near-bottom circulation and mixing driven by the surface eddies over abyssal seamounts. *Prog. Ocean.*, **208**, 102–896.