## Abstract

The southward-flowing deep limb of the Atlantic meridional overturning circulation is composed of both the deep western boundary current (DWBC) and interior pathways. The latter are fed by “leakiness” from the DWBC in the Newfoundland Basin. However, the cause of this leakiness has not yet been explored mechanistically. Here the statistics and dynamics of the DWBC leakiness in the Newfoundland Basin are explored using two float datasets and a high-resolution numerical model. The float leakiness around Flemish Cap is found to be concentrated in several areas (hot spots) that are collocated with bathymetric curvature and steepening. Numerical particle advection experiments reveal that the Lagrangian mean velocity is offshore at these hot spots, while Lagrangian variability is minimal locally. Furthermore, model Eulerian mean streamlines separate from the DWBC to the interior at the leakiness hot spots. This suggests that the leakiness of Lagrangian particles is primarily accomplished by an Eulerian mean flow across isobaths, though eddies serve to transfer around 50% of the Lagrangian particles to the leakiness hot spots via chaotic advection, and rectified eddy transport accounts for around 50% of the offshore flow along the southern face of Flemish Cap. Analysis of the model’s energy and potential vorticity budgets suggests that the flow is baroclinically unstable after separation, but that the resulting eddies induce modest modifications of the mean potential vorticity along streamlines. These results suggest that mean uncompensated leakiness occurs mostly through inertial separation, for which a scaling analysis is presented. Implications for leakiness of other major boundary current systems are discussed.

## 1. Introduction

The Atlantic meridional overturning circulation (AMOC)^{1} connects disparate water masses, depths, and geographical locations (Buckley and Marshall 2016; Lozier 2012) and plays major roles in the broader climate system (Srokosz et al. 2012; Bullister et al. 2013). These include driving a significant fraction of the global atmosphere–ocean meridional heat flux, e.g., an estimated ≈15% at 40°N (virtually all of the oceanic component; Trenberth and Fasullo 2017, their Fig. 3), and influencing the CO_{2} sink in the North Atlantic (Takahashi et al. 2009). Despite its importance, the characterization of three-dimensional AMOC pathways remains incomplete, as does the understanding of their driving mechanisms (Lozier 2012).

A significant portion of the deep (southward) AMOC branch occurs within the deep western boundary current (DWBC). The occurrence and role of the DWBC was predicted by Stommel and Arons (1959), albeit on the basis of assumptions now partially outdated (Ferrari et al. 2016). The DWBC has nonetheless been observed from the subpolar North Atlantic southward to the southern Atlantic, forming an intensified boundary current that carries North Atlantic Deep Water (NADW) along the western Atlantic continental slope (Hogg and Johns 1995; Talley 2011).

However, in recent decades it has become clearer that the DWBC is not the only southward transport branch of the AMOC. A series of float experiments (Lavender et al. 2000; Fischer and Schott 2002; Bower et al. 2009) and tracer analyses (Rhein et al. 2002; Gary et al. 2012; Le Bras et al. 2017) have identified significant loss (leakiness) of material from the DWBC in the Newfoundland (Nfl) Basin. This leakiness was specifically targeted and quantified in the “Export Pathways” experiment (ExPath; Bower et al. 2011). The majority (≈90%) of floats seeded upstream within the DWBC at Labrador Seawater (LSW) depths^{2} leaked to the interior within the Nfl Basin. Much of the leakiness occurred between two large underwater capes (Fig. 1) in the DWBC’s path: Flemish Cap (FC) and the Grand Banks of Newfoundland (GB).

Within the 2-yr lifespan of the floats, ~20% of the floats that leaked out of the DWBC continued southward in the basin interior away from the boundary. Hence, these additional pathways are referred to as interior pathways. These findings of DWBC leakiness and interior pathways represent a significant revision of the classical picture of deep southward AMOC transport being confined to the DWBC. Furthermore, Argo observations (Biló and Johns 2018) and numerical simulations (Gary et al. 2011, 2012; Lozier et al. 2013) suggest that interior pathways continue south farther than the 2-yr ExPath observations demonstrate. Gary et al. (2012) show that 75% of simulated floats initialized within the DWBC and traveling from 44° to 30°N did so in the interior rather than within the DWBC.

Two contrasting views on the dynamical causes of interior pathways were examined hitherto: Gary et al. (2011) have shown that within realistic numerical models and in hydrography, interior pathways were largely collocated with Eulerian recirculation gyres, elevated eddy kinetic energy, and decreased potential vorticity gradients (see also Lozier 1997), all qualitatively consistent with the theory of eddy-driven gyres (Rhines and Young 1982). Furthermore, in the eddy-resolving model examined in Gary et al. (2011), eddy fluxes explained a large fraction of the potential vorticity balance. In contrast, Pedlosky (2018) has shown, in the context of an idealized, steady, flat-basin model, that interior pathways are necessary somewhere in the domain to provide westward flow into the boundary current at all latitudes to its south. That is, since inertial boundary currents need inflow from the east to avoid Rossby wave energy radiation away from the boundary.

Previous studies have thus addressed the locations of DWBC leakiness, interior pathways trajectories, as well as interior pathways dynamics. In contrast, the mechanism underlying the leakiness itself remains unclear. In the following paragraphs, we review four hypotheses that have been posited in the literature.

### a. DWBC–NAC interactions

The DWBC and the more energetic, surface-intensified North Atlantic Current (NAC, an extension of the Gulf Stream), pass quite close to each other in the GB–FC area. The currents come especially close together at the southern tip of the GB, and at the southeast corner of FC, where a large fraction of the floats leaked out of the DWBC. Therefore, interaction between these currents could plausibly cause material to leak from the DWBC (Fischer and Schott 2002; Lavender et al. 2005; Bower et al. 2009, 2011). The high eddy kinetic energy (EKE) values measured (e.g., Carr and Rossby 2001) near the GB region and east of FC imply that the loss of floats from the DWBC may be eddy driven. Additionally, the surface intensification of EKE in the region suggests that the eddies result from instabilities of the surface-intensified NAC.

### b. Inertial separation

Current systems throughout the Nfl region are strongly steered by topography, including the surface-intensified NAC and the DWBC (Rossby 1996; Kearns and Paldor 2000; Fischer and Schott 2002; Lavender et al. 2005). Boundary currents approaching coastal bends may separate from the coast if they have sufficient inertia (e.g., Ou and De Ruijter 1986; Klinger 1994). Pickart and Huang (1995) examined the inertial downstream adjustment of a DWBC-like current to changes in bathymetry in a steady, semigeostrophic, 1.5-layer model. They found that a substantial fraction of the current volume flux was lost to offshore or to a recirculating component, although these solutions lay outside the formal regime of applicability of the semigeostrophic model.

### c. SCVs

Previous studies have found that material may leak from boundary currents via shedding of submesoscale coherent vortices (SCVs) (McWilliams 1985; D’Asaro 1988; Bower et al. 1997). Bottom-reaching prograde boundary currents (propagating left of inshore in the Northern Hemisphere) can generally be expected to develop negative vorticity near the bottom boundary layer due to bottom drag (Molemaker et al. 2015). If the prograde boundary current then separates from the slope, e.g., at a bathymetric cape, the negative vorticity in the bottom boundary layer can cause a roll up into an anticyclonic SCV. Of the ExPath float dataset, Bower et al. (2013) indeed found that three floats became trapped within anticyclonic SCVs formed at the southern tip of the GB.

### d. Instabilities of the DWBC

Oceanic boundary currents may be unstable, and therefore intrinsically favor leakiness (e.g., Cimoli et al. 2017). Motivated by the observed leakiness around FC and GB, the effect of horizontal curvature of bathymetry (and streamlines) upon baroclinic instability was examined by Solodoch et al. (2016), in a two-layer quasigeostrophic model. They found that uniform parallel flow over curved bathymetry has similar baroclinic modal instability growth rates to the case of rectilinear bathymetry (i.e., the extended Phillips model; Mechoso 1980), if the mean flow has a weak barotropic component. The growth rate generally diminishes with increasing mean barotropic flow, an example of the barotropic governor effect (James 1987) in the presence of mean strain.

Based on Eulerian transport measurements at southeast FC and at southeast GB, Mertens et al. (2014, hereafter M14) estimated that out of ≈30 Sv (1 Sv ≡ 10^{6} m^{3} s^{−1}) of southward flowing NADW at southeast FC, 15 Sv are lost offshore before the southern tip of the GB. Biló and Johns (2018) analyzed interior pathways of LSW based on Argo data, and found that of the water leaked from the DWBC within the Nfl Basin, 9.3 ± 3.5 Sv recirculates within the subpolar basin, while 3.2 ± 0.4 Sv continues eastward. These studies therefore show that DWBC leakiness has a significant (on the order of multiple Sverdrups) uncompensated component, i.e., that there is a net loss of mass from the DWBC, rather than simply an exchange of mass with the ambient ocean. This defines a distinction between compensated and uncompensated leakiness, which we shall use in what follows.

In this paper we focus on DWBC leakiness in the Nfl Basin, rather than on the interior pathways which follow leakiness. We combine a new regional model of the northwest Atlantic with historical observations to characterize the leakiness process in detail, and to investigate the mechanisms via which it occurs. In section 2 we describe the regional model, a particle advection code, and the observational datasets used in this study. In section 3 we diagnose the leakiness of the DWBC around FC, using both Lagrangian trajectories (section 3a) and Eulerian mean flow patterns (section 3b). We then quantify the variability in the patterns of leakiness (section 3c) and use budgets of potential vorticity (PV) (section 3d) and energy (section 3e) to investigate the relative roles of mean flows and variability in driving the leakiness. In section 4 we relate our results to the mechanism of leakiness (1–4) summarized in this section, and we put forward a hypothesis for the dependence of leakiness on the geometry of the continental slope. In section 5 we summarize our findings and conclude.

## 2. Methods

### a. Numerical model

We use the Regional Oceanic Modeling System (ROMS), which solves the Boussinesq primitive equations with a free surface (Shchepetkin and McWilliams 2005). ROMS is appealing for use in modeling areas of varying bathymetry, such as the path of the DWBC in the Nfl Basin, due to the combination of terrain-following coordinates that allow fine resolution of the bottom boundary layer and accurate pressure gradient calculation (Shchepetkin and McWilliams 2011) to minimize spurious along-slope flows. The specific ROMS branch we use is the Coastal and Regional Ocean Community (CROCO) branch (Debreu et al. 2012).

We designed a North Atlantic domain ROMS configuration (hereafter GB_{B}), with the Nfl Basin close to the domain center. The model domain is shown in Fig. 1, along with the barotropic (depth-averaged) velocity magnitude averaged over model year 16. Several important topographic features discussed below are annotated in the figure. The domain extends to and beyond the Mid-Atlantic Ridge on the east, and to the Labrador and Irminger Seas on the north. The Gulf Stream enters from the western boundary, following its separation from Cape Hatteras within the parent grid (discussed below).

The GB_{B} horizontal resolution is approximately 2.5 km, which is small compared to the first baroclinic Rossby radius of deformation (*R*_{d} ≈ 10−20 km) in the Nfl Basin (Chelton et al. 1998). Therefore, the model configuration resolves the mesoscale, and possibly a portion of the submesoscale. Fifty (terrain-following) vertical levels are used. At middepths, the typical resolution is then ≈100 m in the deep ocean, and finer in shallower areas, e.g., the DWBC path along the continental slope. Top and bottom coordinate stretching (with stretching factors *θ*_{s} = 6 and *θ*_{b} = 4, respectively) further increases vertical resolution near the top and bottom boundaries. Vertical resolution is approximately 5 m near the surface. At continental slope to continental rise seabed depths (1000–4000 m), vertical resolution near the bottom is ≈15–50 m, respectively. The model bathymetry is derived from the 30-arc-second-resolution Shuttle Radar Topography Mission global product, SRTM30_PLUS (Becker et al. 2009), processed for use in ROMS as described by Renault et al. (2016).

Boundary conditions at open boundary segments are prescribed using an offline nesting approach (Mason et al. 2010). Model variables at the open boundaries are relaxed to values from a coarser parent domain, using radiation-like boundary conditions. These are as described in Marchesiello et al. (2001), except for the barotropic momentum and surface elevation boundary conditions, which are described in Mason et al. (2010). The parent (ROMS) solution is described in Renault et al. (2016). Its domain covers the entire North Atlantic Ocean, with ≈5-km horizontal resolution in the GB_{B} region, and 50 vertical levels as well. The parent configuration was spun up for 14 years using climatological forcing, and subsequently solved for five additional years with time-dependent forcing, corresponding to calendar years 2000–04. For boundary data used in the nesting procedure, in the first four GB_{B} years we use the last four parent solution years, since they were conducted with time-dependent forcing. For each following 4-yr GB_{B} period (years 5–8, 9–12, 13–16), we recycle the same four years of boundary data from the parent solution. Thus interannual variability is statistically limited in the model (see discussion in appendix B). To minimize shock-like numerical artifacts when the forcing cycle is restarted, the last 10 samples of the boundary data cycle (last 10 days of December 2004) are linearly interpolated toward its first sample (1 January 2001). Because radiation boundary conditions are generally not completely free of artifacts, such as boundary reflections, sponge layers are applied near the open boundaries, with a maximum viscosity of 300 m^{2} s^{−1} at the boundary, and a decrease as a cosine quarter cycle to zero over a distance of 25 km from the boundary. Air–sea fluxes are accounted for using bulk formulas (e.g., Fairall et al. 1996), with the atmospheric state interpolated from 6-h interspersed CFSR reanalysis data (Saha et al. 2010).

Vertical subgrid-scale mixing is parameterized via the *K*-profile parameterization (Large et al. 1994). For the tracer advection scheme we initially used the split-rotated scheme “RSUP3,” with the diffusive component aligned with the local neutral plane (Lemarié et al. 2012). However, we found severe numerical issues in our configuration (see online supplemental material). Therefore, we reverted to isopotential alignment of the diffusive part of RSUP3 (Marchesiello et al. 2009). We integrate the model for 16 ocean years, and save 2-day averages of output variables, on which all presented analysis are performed offline. Domain-integrated kinetic energy and available potential energy (Vallis 2017) are examined (not shown), to probe the degree to which the model has spun up. Both quantities have pronounced seasonal cycles, with no clear interannual drift, i.e., the solution appears close to a statistical steady state. Further model validation is presented in appendix B. Given that statistics of domain-integrated energy, water mass properties, and circulation pattern exhibit little variation after year 8 ( appendix B), the presented results (e.g., mean quantities) are based on model years 9–16, unless stated otherwise.

### b. Float datasets

Two observational datasets of subsurface Lagrangian floats are used here. One is Export Pathways from the Subpolar North Atlantic Experiment (ExPath) dataset (Furey and Bower 2009; Bower et al. 2011). In ExPath, RAFOS floats were seeded within the DWBC region west of Orphan Knoll (OKL; Fig. 1). These are isobaric (i.e., approximately depth maintaining) floats that are tracked by acoustic sound sources and hence do not need to surface during their trajectory (unlike Argo floats, see below).

Relative to float datasets used in prior analyses of DWBC leakiness in Newfoundland (Lavender et al. 2000; Fischer and Schott 2002), the ExPath dataset has the advantages that the floats used are not profiling (eliminating contamination of velocity from surfacing), and that the floats were all seeded within the DWBC and just north of the leakiness area, whereas previous floats were seeded farther upstream in the Labrador sea. In numerical simulations the isobaric nature of simulated ExPath-like floats did not appreciably change the interior pathways statistics compared with 3D simulated floats (Bower et al. 2011).

Approximately equal fractions of floats were ballasted for 700 and 1500 dbar (1 dbar ≈ 1 m) depths. Each float drifted for 2 years before resurfacing. We analyze the trajectories of the 55 floats deemed usable in Furey and Bower (2009). Floats positions are generally available with daily resolution. Exceptions include the positions of floats within Flemish Pass (the channel running between FC and GB), which was shielded from sound sources. Due to failure of sound sources during part of the experiment, position triangulation for some trajectories in the continental slope area south of FC were also not possible (Furey and Bower 2009).

The second dataset consisted of a subset of Argo floats (Riser et al. 2016). Argo floats drift at a set “parking depth.” After a typical period of 9 days, the float first descends to 2 km depth, and then ascends to the sea surface, while taking hydrographic measurements. At the surface the float transmits collected data via satellite communication. Then the float descends back to its parking depth, restarting the cycle. We compiled a dataset of all Argo floats that have ever crossed the DWBC cross section along which the ExPath floats were deployed. Specifically, the chosen area is west of OKL, between latitudes 49.5° and 50.5°N and longitudes 49.6° and 47.7°W. We find 67 floats that meet this criterion, with parking depths between 800 and 2000 m, between the years 1998 and 2017. Specifically, the number of floats that have parking depths of 800, 1000, 1500, and 2000 is 3, 43, 18, and 3, respectively. Unfortunately, not all floats in the assembled dataset have actual pressure readings stored from their drift periods, in which case we rely on the programmed parking depth. Despite this caveat, we find this dataset to be a useful complement to the ExPath observations.

### c. Particle advection

The phenomenon in question, Lagrangian leakiness, is most directly addressed in a Lagrangian framework. For that purpose, and for comparison with the float observations, we seed and track passive particles in the velocity fields obtained from the numerical model described in section 2a. We developed a FORTRAN code (named “TrajInt,” for trajectory integration) that allows three-dimensional (3D) integration of particle trajectories given their initial positions at a particular time. Particle advection experiments were performed offline, i.e., after running GB_{B} (section 2a). The main features of the code are described here. The particles are passively advected by solving the advection ordinary differential equation $\u2202tx\xaf=u\xaf$, where $x\xaf=x\xaf\u2061(t)$ is the particle position at time *t*, and $u\xaf$ is the ROMS velocity field interpolated to time *t* and position $x\xaf$. The temporal interpolation is done using cubic splines, the spatial interpolation is trilinear, and time stepping is done using the classical fourth-order Runge–Kutta method. The chosen advection time step is half an hour, which is 1/96 of the GB_{B} saved output rate (2-day averages). At characteristic speeds within the Nfl Basin at middepth of up to 0.3 m s^{−1}, the maximal displacement within a TrajInt time step is ≈0.5 km, or one-fifth of a gridcell side length. Therefore, the time step is likely sufficient to resolve the model output space and time scales. We confirmed this via sensitivity experiments in which trajectories were recomputed with refinement (repeated halving) of the time step size, and found that the differences in the trajectories were smaller than one grid point for at least 10 model days after initialization, when the larger time step is $\u2272$1 h. This time period is comparable to the observed velocity autocorrelation time in subthermocline depths in the northwest North Atlantic (Böning 1988; Lumpkin et al. 2002), hence we consider the convergence satisfactory.

We conduct several particle advection experiments. Floats are initialized along a line (Fig. 1c) within the DWBC, west of OKL, close to the seeding locations of ExPath floats (Furey and Bower 2009), between the 1 and 2.8 km isobaths, and with an initial depth at least 300 m above the bottom. The specified line is chosen since it is located upstream of Flemish Cap, where much of the leakiness occurs, and for comparison with the ExPath dataset (section 2b). The mean model velocity at the seeded depths along OKL is everywhere downstream (approximately southward) within the DWBC, except for a clockwise recirculation at Orphan Knoll.

Experiment 1 (Exp3d) employs a large number of deployments (≈550 000 particles) to get statistically robust estimates of leakiness metrics. We deploy up to 1000 particles at depths of 700 and 1500 m each, uniformly distributed along the entire OKL section, every 10 days between years 9 and 16. At each seeding date, particles were only seeded along the OKL in locations where the meridional component of the 2-day averaged velocity was directed southward.

These particles are advected for 200 days each. Because the velocity autocorrelation (integral) time scale in this region is generally between 5 and 10 days (Böning 1988; Lumpkin et al. 2002), seeding more often than 10 days would not likely have been effective in terms of relative contribution of additional effective degrees of freedom. Experiment 2 (Exp3dMean) uses the time mean (years 9–16) velocity field in place of the 2-day-averaged output velocity field. Its purpose is to delineate the mean offshore flow pathways and compare time mean with variable leakiness. We deploy 1000 particles each at depths of 700 and 1500 m, uniformly distributed along the OKL section. Only floats that drift southward past Orphan Knoll are considered here to delineate the mean DWBC trajectory, and thus are used in the analysis of Exp3dMean.

## 3. Results

### a. Lagrangian leakiness pathways

An estimated 73%–84% of all ExPath floats were lost (leaked) from the DWBC to the interior before circumnavigating FC (see supplemental material), demonstrating its relative importance in DWBC leakiness within the Nfl Basin. We therefore focus mostly on the FC area in this section. Figure 2a shows the trajectories of the ExPath floats around the time that each float makes its first (offshore) crossing of the 4 km isobath at FC, which is approximately the offshore limit of the DWBC. This diagnostic parameter is useful as floats which crossed offshore around FC^{3} do not appear to have reentered the DWBC (see individual trajectories in Furey and Bower 2009). While some do cross the 4 km isobath back near FC, these either recirculate immediately offshore again, or travel close to the same isobath upstream, apparently entrained in the NAC.

The distribution of trajectories leaving the DWBC around FC (Fig. 2a) suggests that leakiness of ExPath floats occurs in three main FC subregions (leakiness hot spots): at the northeast (NEC) and southeast (SEC) corners of FC, and in the southern face (SF) just following SEC. Their approximate locations are marked in Fig. 1b. The concentration of ExPath floats leakiness near SEC was previously reported by Bower et al. (2011).

At the NEC and SEC, ExPath floats leave the DWBC via trajectories that are oriented almost directly offshore. In the SF area, the offshore velocity component is weaker relative to the alongshore component, but some of the floats abruptly turn back upstream (approximately northeastward) midway through the SF. The hot spots are approximately collocated with local maxima in topographic changes: convex curvature at NEC and SEC, and a two- to threefold increase in bottom steepness in the SF area (section 4). Figure 2 also shows that as floats travel offshore, they tend to turn cyclonically, consistent with vortex stretching assuming conservation of potential vorticity of the layer below the pycnocline. The question of what sets the locations of the leakiness hot spots is discussed in sections 3d, 3e, and 4.

We further examine Lagrangian pathways in observations, by performing a similar analysis (Fig. 2b) on Argo floats traveling south within the DWBC.^{4} The subset of Argo floats is described in section 2b. While the temporal resolution of Argo floats locations is an order of magnitude lower, the clustering of the Argo floats’ crossings of the 4 km isobath is qualitatively similar to that of the ExPath floats.

Next, we examine Lagrangian pathways of particles seeded within the numerical model (section 2b), beginning with a small subset of the seeded particles for a qualitative visual comparison with the floats. Figure 2c is identical to Figs. 2a and 2b, but displaying the trajectories of a random batch of 60 model particles from Exp3d (section 2c)—30 from each seeding depth (700 and 1500 m). The leakiness hot spots and other related properties described above for the ExPath floats are largely reproduced in this case. These results are consistent in other random samples of the floats from Exp3d (not shown).

To examine leakiness within the full set of (~550 000) model particles, we first plot the distribution of the locations at which each particle in Exp3d first crossed from the DWBC to offshore of the 4 km isobath (Fig. 3a). We find the same clustering as suggested in Fig. 2, i.e., the offshore crossing density is highest at the NE corner, SE corner, and SF. The pattern appears qualitatively consistent with the ExPath observations (circles superimposed in the panel). For a quantitative comparison, we apply a two-sample Kolmogorov–Smirnov (KS) test. The two sample sets are the ExPath and Exp3d offshore crossing locations. We partition the 4 km isobath into consecutive 50-km-long sections, and count (bin) the number of floats or particles crossing each section. The cumulative distribution function (CDF) of the ExPath floats (*F*_{f}) offshore crossings is then compared with the same CDF for the Exp3d particles (*F*_{p}). The KS test statistic, defined by *D* = max_{n}|*F*_{f}(*n*) − *F*_{p}(*n*)|, is then compared with the theoretical KS distribution. The result from the comparison is that the two distributions are statistically indistinguishable (*p* value = 0.96). This indicates that the observed and modeled float trajectories are consistent with one another, to the extent that differences between them could be distinguished statistically.

Finally, we calculate the Lagrangian-mean velocity, based again on the full number of Exp3d model particles. The Lagrangian-mean velocity is defined for our purpose as the average velocity within a grid cell of all particles that have crossed it. Note that this is a conditional average, in that it includes solely particles that were released within the DWBC, and in that we apply a further restriction by including particles only before their first crossing of the 4.2 km isobath offshore. This differs from an Eulerian average because, for example, the velocities of parcels carried by intrusions into the DWBC from offshore will not directly contribute to the calculation. In Figs. 3b and 3c, the Lagrangian mean velocity is displayed and decomposed into along and cross-bathymetry components (hereafter *υ*_{a} and *υ*_{c}, respectively).^{5} Only statistically significant (see supplemental material) *υ*_{a} and *υ*_{c} values are displayed. We similarly calculate the average Lagrangian eddy kinetic energy, $EKE=\u2061(1/2)\u2061(\upsilon \xaf\u2212\upsilon \xaf\xaf)2\xaf$, where $\upsilon \xaf$ is the velocity of an individual particle sampled within a grid cell, and an overbar again denotes an average over all such samples within a single grid cell. In Fig. 3 only the results based on 1500-m-deep model particles are shown. The same diagnostics for the 700-m-deep model particles are quantitatively similar (supplemental material).

The along-bathymetry velocity component (Fig. 3b) exhibits a maximum along the path of the DWBC on the continental slope. The cross-bathymetry component (Fig. 3c) shows that offshore Lagrangian mean velocities occur in patches stretching across the DWBC to its offshore edge (≈4 km isobath) at the identified Lagrangian leakiness hot spots (NEC, SEC, and SF). The Lagrangian mean velocity follows pathways from the DWBC core to NEC, SEC and SF, which is most easily seen via the Lagrangian mean velocity vectors overlaid on Fig. 3d.

The EKE (Fig. 3d) is considerably lower (by roughly 50%) at leakiness hot spots NEC and SEC compared to adjacent patches along the same isobaths, suggesting that the cross-isobath Lagrangian transport at these hot spots is primarily due to an Eulerian mean flow. In contrast, if the Lagrangian mean offshore flow were an eddy-forced or eddy-rectified flow, one would expect it to be associated with elevated EKE values. That may be the case at SF, where EKE is indeed locally elevated (see sections 3c and 3e as well).

In summary, the analysis presented in this section shows that in observations Lagrangian leakiness trajectories are clustered in a few key locations (NEC, SEC, SF). Additionally, the numerical model compares well with the observations, and using a much larger number of (numerical) particles demonstrates that these leakiness hot spots are associated with high Lagrangian mean offshore velocities, offshore deflections of the peak alongshore velocity (*υ*_{a}) upstream, reductions in the magnitude of *υ*_{a}, and (except at SF) low variability (EKE).

### b. Eulerian characterization of leakiness

In section 3a we quantified the Lagrangian leakiness via the Lagrangian mean offshore flow. The Lagrangian mean flow may be locally represented as the sum of the Eulerian mean flow, and the rectified eddy flow. In the present section we analyze the Eulerian mean flow over the same time period (years 9–16), and thereby deduce the contribution of rectified eddy transport to the Lagrangian mean offshore flow.

We begin by examining cumulative (Eulerian) offshore transport on the 3 km isobath (Fig. 4),^{6} in comparison with the observational estimate of M14 (section 1). As in M14, we decompose the transport into densities greater or smaller than *σ*_{θ} = 27.68 kg m^{−3} (Fig. 4a), approximately the upper boundary of LSW. This partitioning is also useful because the bias in model isopycnal depths significantly decreases for *σ*_{θ} < 27.0 kg m^{−3} ( appendix B and Fig. B2).

Although there is a substantial offshore transport (~4 Sv) at OK, it is compensated by shoreward flow immediately downstream, resulting in negligible net offshore transport around OK (Fig. 4a). In contrast, around FC there is an offshore transport of 13–16 Sv, which is uncompensated in the *σ*_{θ} ≥ 27.68 kg m^{−3} density range. The offshore transport rate (slope of the curve) greatly increases around the SE corner and downstream from it (around SF), where much of the Lagrangian leakiness is clustered. Additionally, 3–5 Sv are lost around the southern tip of the GB. The cumulative loss from FC to GB is consistent with the M14 estimate, and our analysis further constrains (within the numerical model) the along-slope distribution of the offshore transport. Results for *σ*_{θ} ≤ 27.68 kg m^{−3} show a similar pattern, with ~3 Sv lost around FC, and ~1 Sv lost around GB. An examination of the cumulative offshore transport in depth layers (Fig. 4b) reveals that the transport is largely depth independent down to 2.5 km depth (and slightly surface intensified).

We now examine the Eulerian mean circulation patterns around the DWBC leakiness hot spots. Figure 5 shows the Eulerian mean velocity streamfunction on two representative isopycnal surfaces: *σ*_{1} = 32.43 and *σ*_{2} = 37.014 kg m^{−3}, averaged over years 9–16. The upper surface lies between depths of 800 and 1750 m in the DWBC (Fig. 5a), similar to ExPath floats and to the Lagrangian analysis in the previous section. It also corresponds to typical LSW depths (Bullister et al. 2013; M14). The lower surface lies between depths of 1800 and 2700 m in the DWBC, corresponding to lower LSW or upper Overflow Water. These surfaces are hereafter referred to as LSW and lLSW for clarity. However, we do not suggest they correspond accurately to observed water mass properties ( appendix B). The streamfunction is calculated by an adaptation of a flood-fill algorithm (supplemental material). Closed streamlines with inner minima (maxima) are cyclonic (anticyclonic) recirculations, and streamfunction values are only meaningful up to an addition of a global constant.

Figure 5 shows that the Lagrangian leakiness hot spots (NEC, SEC, SF; section 3a) coincide with mean streamlines exiting the DWBC. This indicates that the leakiness is at least partially attributable to Eulerian mean offshore flows at the NEC, SEC and SF hot spots. At NEC separating streamlines are apparent only in the deeper density surface (*σ*_{2} = 37.014 kg m^{−3}) plotted, although they appear if more streamlines are plotted in the shallower density surface (*σ*_{1} = 32.43 kg m^{−3}) as well. This is consistent with the larger offshore flux near SEC (Fig. 4).

Figure 5 also reveals the existence of three closed cyclonic recirculations with radii of *O*(100) km immediately offshore of the DWBC around FC. These recirculations stand between the DWBC and the NAC, complicating the potential NAC influence on DWBC leakiness (mechanism 1, section 1). Similar cyclonic recirculations around FC were reported in circulation estimates based on profiling floats (Lavender et al. 2005), and in numerical simulations by Xu et al. (2015), which noted that the recirculations are consistent with the distribution of Tritium [see also Fig. 2a in Biló and Johns (2018) and Fig. 1a in Getzlaff et al. (2006)]. The separating streamlines at NEC and SEC do not return to the DWBC, but rather turn (around the offshore recirculations) cyclonically east and northward post separation, and appear to join or travel adjacent to the NAC. The cyclonic turning of separated streamlines is visually similar to the cyclonic trajectories of the Lagrangian particles after they have left the DWBC (Fig. 2). These circulation patterns (including separation and recirculation) are similar on both density surfaces shown in Fig. 5, which are separated by around 1 km vertically. Similar results are found when the streamfunction is computed for other, intermediate density surfaces, or for the depth-integrated flow (not shown).

We investigated the role of eddies transport by comparing the thickness-weighted average (TWA) velocity streamfunction (Young 2012) to the simple time-averaged velocity streamfunction discussed above in this subsection. The patterns (not shown) and speeds are nearly indistinguishable between the two different averages. The mean speed difference in the area shown in Fig. 5 is 0.002 m s^{−1}. The maximal difference (≈0.01 m s^{−1}) occurs around the SF hot spots and in the confluence zone offshore of SF. This is consistent with the greater EKE diagnosed at SF relative to NEC or SEC from model particle motions (section 3a). Thus, the eddy-rectified circulation is generally negligible in comparison with the mean Eulerian circulation on these isopycnals. This is consistent with the qualitative similarity between the Lagrangian mean (section 3a) and Eulerian mean (Fig. 5) offshore flow velocity distributions. Furthermore, it suggests that the Eulerian mean flow accounts for the offshore transport of Lagrangian particles at the leakiness hot spots.

### c. Robustness of spatial patterns of separation

Diagnostics presented in the previous two subsections suggest that leakiness occurs, at least partially, as a spatially localized and temporally steady (time-mean) offshore flow pattern. In this section we examine the following question: how representative is the diagnosed time-mean circulation pattern of the time-varying circulation patterns? The answer permits dynamical interpretation of the mean circulation; for example, if the mean offshore flow is locally the result of infrequent but intense offshore flow events, while most of the time the velocity is inshore, the mean flow state itself would be atypical. Such a scenario may be consistent with rare but intense external events, e.g., NAC-derived eddies propagating inshore, causing the mean offshore flow. We will see, however, that the mean circulation patterns are in fact statistically quite representative of instantaneous patterns.

In Fig. 6 we present statistics of the cross-bathymetry flow as a function of distance along the 4 km isobath. The velocity is averaged between depths of 700 and 1500 m, but the findings are representative of velocity statistics in other layers between 500 m depth and the sea floor (not shown). Figure 6a shows the Eulerian mean, median, and mode of the cross-bathymetry flow, as well as the Exp3d Lagrangian mean cross-isobath velocity *υ*_{c}. Here the mode was defined relative to 1 cm s^{−1} resolution binning of all samples. Figure 6b shows a histogram of the cross-isobath velocity at SEC, the location of which is marked in Fig. 6a and in Fig. 1b. In constructing the histogram, all time samples from locations up to two grid cells distant from the indicated point along the isobath were used. The error in estimation of the Eulerian mean, $std/Ne$, is everywhere <0.01 m s^{−1}, where std is the standard deviation over the *N* = 1460 time samples (years 9–16, 2-day intervals), and *N*_{e} = *N*/(10/2) is the number of effective degrees of freedom, assuming an integral time scale of 10 days (section 2c). Hence, the mean offshore velocity at NEC, SEC, and SF is statistically significant (*p* < 0.05).

The Eulerian mean and median are very close to one other along this section, and particularly so at the mean leakiness hot spots (Fig. 6a). The mode fluctuates strongly, but generally follows the mean values well over length scales $\u227350\u2212100km$. The mode is very close to the mean at the SEC and (slightly less so at) NEC. Figure 6b shows that offshore flow is indeed the typical occurrence, and the distribution is quite symmetric around the mean. The distributions at NEC (not shown) and SEC are both center-heavy (an excess kurtosis magnitude of |*κ*| ≈ 0.25), and symmetric (a skewness magnitude of |*γ*| ≤ 0.1). At SF the distribution (not shown) remains center-heavy, although to a lesser degree: (*κ* = 0.75 and *γ* = 0.5). In summary, the Eulerian mean offshore flow is statistically representative, i.e., typical values are close to the mean. In the supplemental material we present a cluster analysis that demonstrates that the spatial Eulerian pattern of mean separation (including separating streamlines) is statistically representative as well, in a similar sense to that described above.

The Lagrangian and Eulerian means are very similar around most of FC, confirming that the eddy-induced rectified offshore flow is relatively low in this area (Fig. 6a). Along eastern FC outside of the hot spots the Lagrangian mean is generally slightly higher, increasing mean leakiness there (cf. with Fig. 4). Another exception is that at SF the Eulerian mean only accounts for around 50% of the Lagrangian mean offshore velocity, suggesting that the remainder of the transport is due to the rectified eddy mean flow. This is consistent with the elevated Lagrangian mean EKE at SF (section 3a). However, the Lagrangian mean offshore flow is also generally weaker at SF than it is at NEC or SEC.

We emphasize that although Lagrangian and Eulerian mean velocities are almost identical at the leakiness hot spots, particularly NEC and SEC, time variability nonetheless has a nonnegligible influence on the Lagrangian leakiness. In Exp3d, more than 90% of the particles are exported across the 4 km isobath before they can reach the GB (longitude ≈−55°). In contrast, in Exp3dMean, in which particles are advected by the time-mean velocity fields, only ~59% (35%) of the particles initialized at a depth of 700 m (1500 m) are exported across the 4 km isobath before they can pass GB (longitude ≈−55°). Therefore flow variability contributes substantially to the leakiness. This contrasts with the high quantitative similarity demonstrated between Eulerian and Lagrangian mean offshore flow, and the relatively low magnitude of eddy-rectified offshore flow.

Although mean streamlines leave the DWBC offshore at the three identified leakiness hot spots, fewer particles in Exp3dMean reach those hot spots. This is to be expected given approximate planetary vorticity (*f*/*h*) conservation. Indeed, the DWBC at and upstream of the particle seeding locations (Fischer et al. 2004; Bower et al. 2011) is confined inshore of the 3 km isobath, but around FC (Fig. B2) and GB (Schott et al. 2004) it extends to the 4 km isobath. With temporal variability (i.e., in Exp3d), particles cross *f*/*h* lines and populate the mean leakiness hot spots offshore of the 3 km isobath, where mean velocity can propel them further offshore. This may be a manifestation of the phenomenon known as chaotic advection (e.g., Shepherd et al. 2000; Rypina et al. 2010). Chaotic advection (Aref 1984) refers to complex Lagrangian trajectories which often result even from simple Eulerian fields by the kinematics of superimposed eddies and nonuniform mean circulation.

### d. PV distribution and balance

Given the separation of the mean flow from the DWBC into the interior, one might ask: how does the mean flow cross the dynamical barrier presented by the cross-bathymetry PV gradient? To address this, we now examine the TWA PV budget (Smith 1999; Young 2012). The TWA of a variable *a*, and its deviation from TWA, are defined by $a^=\u2061(ha\xaf/h\xaf)$, and $a\u2032=a\u2212a^$, respectively. Here an overbar denotes a time average and *h* = −*ρ*_{0}(∂*z*/∂*ρ*) is the isopycnal “thickness density,” where *ρ*_{0} = 1027.4 kg m^{−3} is a constant reference density. All averages are performed on a selected isopycnal. The TWA PV $q^$, and its balance, are then respectively defined by Smith (1999),

From left to right, the terms in (1b) are time tendency, advection of the mean PV by the mean velocity, the eddy PV flux divergence, and all nonconservative terms (nct) lumped together. Figure 7 shows the TWA PV, its budget, and the TWA eddy enstrophy, all calculated on *σ*_{1} = 32.43 kg m^{−3}, which is the same (LSW) isopycnal as in Figs. 5a and 5b. The analysis was also repeated (not shown) on *σ*_{2} = 37.014 kg m^{−3} (as in Figs. 5c and 5d), and we find that the patterns described below are similar on this deeper isopycnal as well.

The PV (Fig. 7a) is generally lower near the western boundary, due to the low stratification imparted to LSW in its formation via deep convection (Talley and McCartney 1982; Rhein et al. 2002). We observe that in addition to the large-scale offshore gradient, low-PV pockets extend away from the DWBC along the mean flow streamlines at the NEC, SEC, and SF areas, and into the adjacent recirculations. Thus, separation occurs across (up) the mean PV gradient, and mean PV dilution or modification by the eddy and nct terms is sufficiently weak that low PV contours protrude offshore. Figure 7b displays TWA potential eddy enstrophy, $Z=\u2061(1/2)\u2061(q\u2033)2^$, which peaks inshore within the DWBC, and upstream of FC. Values are lower further offshore, including the areas offshore of the leakiness hot spots.

The conservative terms of the PV equation (1b) are displayed in Figs. 7c and 7d. At the leakiness hot spots the offshore mean PV advection (Fig. 7c) results in a negative contribution to the local PV tendency, because PV increases along the path of the TWA flow. Elsewhere the pattern at the offshore edge of the DWBC is generally less coherent. The eddy PV flux divergence (Fig. 7d) approximately matches the pattern and amplitude of the mean PV advection term, but is opposite in sign (with pattern correlation = −0.85). Therefore, at separation streamlines, mean PV advection is upgradient and balanced by eddy PV flux convergence, with only a secondary role for nonconservative processes.^{7} The magnitudes of the conservative PV terms are generally largest just downstream of the leakiness hot spots and within the cyclonic recirculations. Because eddy PV flux divergence is order one in the PV budget, we gauge its influence on the mean PV distribution by evaluating the change in mean PV along a mean streamline. We specifically pick a streamline that separates from the DWBC (at SEC), marked in Fig. 7a. The PV values at three points along this streamline, and an additional point along a mean NAC streamline, are given in the caption. The total growth in mean PV along the DWBC streamline after its separation (occurring between points 1 and 2 in the plot) is ≈10% of the DWBC–NAC mean PV difference (between points 1 and 4). The maximal cumulative growth along the streamline, ≈20%, occurs at point 3.

To summarize, eddy PV flux divergence is a first-order term in the PV budget, largely balancing the offshore PV advection. However, cumulative mean PV change along mean separating streamlines (which is dominated by eddy stirring), is relatively modest. In contrast, if leakiness occurred mainly via eddies derived from the NAC (mechanism 1, section 1), then along a separating streamline eddy stirring should result in significant [*O*(1)] changes in PV relative to the NAC–DWBC mean PV difference. Furthermore, under mechanism 1, we would expect that variability would either peak offshore at the eddy source (NAC) or be more homogeneous in between the NAC and DWBC. That does not appear to be the case, based on our diagnostics of the eddy potential enstrophy and the Lagrangian EKE (section 3a).

### e. Energy conversions

To more completely characterize the role of eddies in leakiness of the DWBC, we examine in this section the energy balance around FC. Given that the PV budget is primarily a balance between the inviscid (mean advection and eddy flux divergence) terms, we focus on the conversion terms between the mean and eddy energy reservoirs. We define the mean kinetic energy (MKE), mean potential energy (MPE), and eddy kinetic energy (EKE) as

Here the time mean (between model years 9 and 16) and deviation from the mean are denoted by the overbar and prime symbols, respectively. The MKE [Eq. (3)] and EKE [Eq. (4)] budgets are given by Harrison and Robinson (1978):

Here we define

where double indices imply summation. The left-hand-side terms of Eqs. (3) and (4) are the time tendency terms, which we find are negligible compared with all other terms in each equation. The *T* symbols represent nonlocal transport and pressure work terms, whereas nct symbols denote nonconservative terms. The local energy conversion terms in both the MKE and EKE budgets are the Reynolds stress work terms (RSW), and potential energy conversion terms (PEC). In the MKE equation (3), RSW_{e2m} is the eddy-to-mean Reynolds stress work, which corresponds to conversion from EKE to MKE. The PEC_{m} term corresponds to conversion of potential energy to MKE. In the EKE equation (4), the RSW_{m2e} term is the mean to eddy Reynolds stress work, and PEC_{e} corresponds to conversion of potential energy to EKE. We calculate the local energy conversion terms in model (terrain-following) coordinates, and later sample them on the time mean *σ*_{1} = 32.43 kg m^{−3} (LSW) isopycnal surface, for comparison with the previous diagnostics on the same surface. However, very similar conversion patterns are obtained at other LSW depths, and in a full depth integral (not shown). The results are also robust in that they vary little when alternative averaging periods are used in place of years 9–16, e.g., when averaging over individual years.

The EKE to MKE conversion term, RSW_{e2m}, is displayed in Fig. 8a. If the mean flow is driven by eddy fluxes, that should be reflected by positive values of Reynolds stress work by the eddies. Upstream of the leakiness hot spots, the mean flow is accelerated by positive RSW_{e2m}. However, RSW_{e2m} is low and close to a sign change at the leakiness hot spots, indicating that mean separation is not forced energetically by eddies. In particular, RSW_{e2m} becomes negative upstream of the SEC separating streamlines. The MPE to MKE conversion term, PEC_{m} (not shown), is positive at and following the mean separation areas, due to column stretching and downwelling which occurs in the offshore crossing of isobaths.

The conversion term of MPE into EKE (by slumping of sloping isopycnals, for example in the form of baroclinic instability) is shown in Fig. 8b.^{8} This term is positive over most of the extent of the DWBC in the figure, including upstream of FC. At and downstream of separation points, as well as within the adjacent recirculating streamlines, PEC_{e} values are elevated. It is plausible that the separation of the mean flow from the continental slope contributes to the growth in PEC_{e}: a parallel flow over sloping bathymetry exhibits lower linear growth rates compared to flows over a flat bottom or a free jet crossing isobaths (Mechoso 1980; Gula and Zeitlin 2014; Solodoch et al. 2016).

We note that the release of MPE to EKE may locally contribute to the following conversion of EKE to MKE. Indeed, in most areas within the DWBC where RSW_{e2m} is positive, PEC_{e} is also positive. This path to MKE increase does not appear to facilitate the leakiness itself, but rather appears to be a consequence of the mean flow departing the continental slope, as noted above. However, positive PEC_{e} upstream of NEC may be sufficient to locally increase EKE, and so may contribute to the “diffusion” of particle trajectories across the DWBC and toward the leakiness hot spots (section 3c).

## 4. Discussion

We now relate our results to the potential mechanisms of DWBC leakiness identified in section 1. We first argue that our results are consistent with inertial separation of the DWBC at FC, and we briefly review previous theoretical works on the conditions for inertial separation. Given that existing theories are not applicable to the DWBC, in section 4b we present a scaling argument for inertial separation. Finally, possible dependence of separation on model resolution and physics is briefly discussed (section 4c).

### a. Mechanism of DWBC leakiness

Taken together, the results presented in section 3 are consistent with inertial separation of the DWBC at FC (mechanism 2 in section 1). The evidence in support of this claim is as follows:

The Lagrangian leakiness hot spots coincide with relatively sharp bathymetric variations, namely, convex turns and steepening of the continental slope (section 3a).

The offshore Lagrangian mean flow coincides with the Eulerian mean flow, which is a typical (rather than intermittent) offshore flow pattern at the leakiness hot spots.

Mean DWBC PV contours are deformed in the offshore flow direction at leakiness hot spots NEC and SEC, indicating advection of PV from the continental slope into the open ocean by the separating mean flow. The PV exhibits relatively modest changes, mainly due to eddy stirring, along mean separating streamlines (section 3d).

Separating streamlines (Fig. 5) and floats leaving the DWBC (Fig. 2) tend to turn anticlockwise, consistent with potential vorticity conservation and thus vertical stretching.

The main hypothesis put forward in previous studies is that high NAC-generated EKE is responsible for DWBC leakiness (mechanism 1 in section 1), which may be expected based on the spatial proximity between the currents at separation areas. While rectified offshore eddy transport indeed accounts for ≈50% of the Lagrangian mean offshore velocity at SF, it is negligible at NEC and SEC. Eddying effects also play a significant role in shifting particles from the upper continental slope toward the leakiness hot spots at NEC and SEC via chaotic advection, as revealed in a comparison of Exp3d with Exp3dMean. However, the majority of the uncompensated, cumulative, leakiness occurs as an Eulerian time-mean offshore flow (sections 3b and 3c). Additionally, the mean offshore flow does not appear to be directly forced by either internally or externally generated eddies (mechanisms 1 and 4 in section 1); baroclinic eddy production is relatively weak within the DWBC, and Reynolds stress work by the eddies on the mean flow is negative close to the separation of mean streamlines (section 3e). Finally, we do not directly address the possible role of SCV formation in DWBC leakiness (Bower et al. 2013) here (mechanism 2 in section 1), a topic reserved for future study.

We note that the Nfl Basin lies close to the latitude of zero wind stress curl. This marks the border between the subpolar and subtropical wind gyres in Sverdrup theory. Furthermore, the Sverdrup “streamfunction” predicts 10–20 Sv leaving the western boundary near FC (Talley 2011, Fig. S9.3), which is similar to the observed (uncompensated) leakiness of the DWBC (section 3b). However, previous studies have demonstrated that Sverdrup balance is significantly compromised in the subpolar gyre due to bottom pressure torque (Hughes and DeCuevas 2001; Spence et al. 2012) and eddy terms (Gary et al. 2011). In a high-resolution (2 km) numerical model, Le Corre et al. (2019) show that in the North Atlantic subpolar gyre, Sverdrup balance does not hold even to first order (including in the gyre interior). Rather, bottom pressure torque, nonlinearity, and other terms are dominant. Thus, it is not clear whether boundary current separation should be expected in the vicinity of the latitude of zero wind stress curl in the subpolar gyre.

The results of Le Corre et al. (2019) further show that nonlinear terms, representing fluxes from the slope region, are a dominant positive (cyclonic) term in the interior vorticity balance of the simulated subpolar gyre. Le Corre et al. (2019) show that this flux is high near FC and is mainly due to the mean rather than the eddying circulation. We interpret this result as supportive of locally determined mean inertial separation. Therefore, although large-scale gyre constraints may play a role in DWBC separation, we focus on the local constraints and dynamics here and leave the role of the gyre-scale circulation in the FC separation as a topic for future study.

We likewise do not analyze here nonlocal energy transfer terms (pressure work, and eddy transport of EKE). These terms may be important in interactions between the DWBC and the NAC, but remain outside the scope of this work. This is related to the concept of boundary currents collision (Cessi 1991; Agra and Nof 1993), which occurs when two western boundary currents converge, and can substantially modify their latitude of separation. In the present case, however, the DWBC and NAC occupy distinct ranges of isobaths, and the NAC separates farther north than the (partially separating) DWBC.

### b. A scaling analysis of inertial separation

Given the evidence for inertial separation as the primary mechanism of DWBC leakiness (section 4a), we now examine the distributions of offshore flow *υ*_{c} and bathymetric changes. Figure 9 shows that larger offshore values of *υ*_{c} tend to be collocated with sharp increases in curvature and steepness (see Fig. 4 as well). Everywhere around FC *υ*_{c} > 0, and it peaks around NEC and SEC. This occurs to a lesser degree around GB, where steepening and curvature are not as pronounced as around FC. Elsewhere, away from FC, *υ*_{c} is either lower in magnitude or oscillates in sign along other bathymetric features. More quantitatively, we plot cross correlation of offshore velocity with isobath curvature and steepening along the 3 km isobath in Fig. 9. The correlation between steepness gradient and offshore velocity reaches *r* = 0.47 at a downstream lag of 73 km, while the correlation between curvature and offshore velocity reaches *r* = −0.56 at a downstream lag of 45 km (see appendix C). These correlations are consistent with inertial separation initiated by sharp changes in the geometry of the continental slope. However, the correlations are at least partially due to meandering of the DWBC along the entire length of the isobath, rather than just the separation points around FC.

Inertial separation of currents flowing around capes or ridges was studied theoretically by Pickart and Huang (1995), Ou and De Ruijter (1986), Klinger (1994), and Jiang (1995).^{9} As reviewed in section 1, Pickart and Huang (1995) specifically studied a DWBC-like current traversing a ridge, and demonstrated that a significant flux is lost to offshore. However, these studies all made the semigeostrophic approximation, which is invalid when along-stream variations are of similar or shorter length scales than cross-stream variations. At FC, the radius of curvature at the SEC is around 10 km, and a few tens of km at the NEC. In comparison, the width of the DWBC (50–100 km) is considerably larger. Hence, at the convex corners, where much of the separation happens, the semigeostrophic approximation fails, and these models become inapplicable.

Furthermore, these works (except Pickart and Huang 1995) all find separation happens within their respective models due to surface outcropping of a density surface bounding a surface current from below. The DWBC is not a surface current, but rather has significant deep and depth-independent components. Indeed leakiness and separation at FC are to a large degree depth independent in our numerical model (e.g., Figs. 4, 5, and 9). Hence theories derived for separation via isopycnal outcropping in a buoyant boundary current are not applicable to the DWBC. Finally the cited works do not cover downstream changes in bottom slope, which the scaling analysis we employ (section 4b) suggests is a significant factor at FC. In fact, most of the cited works assumed a flat bottom and vertical sidewalls.

Greenberg and Petrie (1988) presented a barotropic numerical model over an Nfl-like bathymetry, where the only prescribed current (by boundary conditions) is DWBC-like. The solution indeed displayed significant offshore transport around FC (their Fig. 3a), consistent with inertial separation. A caveat is that the eastern domain boundary was very close to FC.

In the absence of a closed-form theory of inertial separation relevant for the present case, we present a scaling analysis that seeks to determine a simple condition for cape separation of a prograde deep boundary current. The analysis first assumes that mean streamlines continuously curve around a convex corner, while conserving PV (as noted above, mean PV changes moderately along mean streamlines around FC in GB_{B}, including at separation). Then, a condition is derived under which offshore excursions or recirculations form. Because we do not explicitly solve for or use constraints related to the global streamfunction, it is still conceivable that offshore excursions may be followed by meandering and reattachment downstream, rather than permanent separation. However, we are concerned with bathymetric turns of large angles (~90° for SEC), which are likely more favorable for permanent separation. From a kinematic standpoint, larger bathymetric turn angles require larger inshore displacements to compensate for a set offshore detachment distance. A contributing factor in that regard is that for a prograde slope current, vorticity stretching upon offshore excursions may enhance separation, as discussed below.

We use the fact that vorticity can generally be decomposed into shear vorticity *ζ*_{s} = ∂_{n}*U* and curvature vorticity, *ζ*_{c} = *U*/*r*_{c}. Here *t* and *n* are the tangential and normal components of the “natural” coordinate system (Holton 1973), where *t* is locally downstream, and *n* is to its left. The subscript *n* denotes differentiation in the *n* direction, *U* is the magnitude of the DWBC velocity, which is directed in the *t* direction, and *r*_{c} is the streamline radius of curvature, negative for clockwise turns as for, e.g., the DWBC around FC NEC or SEC. The expression for PV in isopycnal coordinates is then

where *h* is the thickness (distance between two chosen density values), and *f* is the Coriolis parameter. The downstream change in PV along a mean streamline under the above assumptions is

The first term is the change in planetary vorticity, which is neglected over the scales relevant in the present analysis, i.e., *O*(100) km. The last is the vorticity stretching term, which is linearized under the assumption that vertical excursions are a modest fraction of the total thickness for mesoscale motions.

At the turn itself, $d\zeta c\u2248U^/Rc$. Here $U^$ is a velocity scale, while *R*_{c} is the bathymetric radius of curvature. The scale of the downstream change in shear vorticity is written as *dζ*_{s} ≈ Δ*U*/*W*, where *W* is the current width upstream, and Δ*U* is the downstream change in cross-current shear integrated in the positive *n* direction. Note that if width decreases (increases) downstream, then Δ*U* below is an overestimate (underestimate). We assume that the change in current width will be a modest fraction, if the current does not partially separate. Therefore, we have a scaling equation relating cross-stream shear changes to bathymetric curvature and deepening (the latter related to steepening across the current):

At FC, the turn of the DWBC is clockwise, hence *R*_{c} is negative, especially at the leakiness hot spots NEC and SEC. Additionally, a steepening occurs at and prior to NEC, SEC, and SF. The deepening should be accompanied by vertical stretching (*dh* > 0), given that the current fills a significant part of the water column (Fig. B2). By (8) both clockwise curvature and vertical stretching each add to a drop in velocity per unit distance offshore (Δ*U*), tending to reduce current flux downstream along the isobaths. In the next two paragraphs we examine the contribution of each of these terms in turn.

It follows from Eq. (8) that if the radius of curvature *R*_{c} is similar in magnitude to or shorter than the current width *W*, then its contribution to Δ*U* is of the same magnitude as the mean current speed. At the outer rim of the current the added shear is then of sufficient magnitude to reverse its direction,^{10} with speed comparable to $U^$. The associated large relative reduction in downstream flux in the steady circulation is a manifestation of inertial separation. As noted above, the radius of curvature at the SEC (NEC) is around (a few times) 10 km, while the width of the DWBC is 50–100 km. Thus, *R*_{c} is in fact significantly lower than *W*.

The stretching term in Eq. (8) (second right-hand-side term) has a similar effect in reducing downstream along-isobaths flux as does the curvature term, and is of a similar magnitude. The deepening of streamlines originating on, say, the *h* = 3 km isobath upstream of the SEC, is greater than *dh* = 500 m, resulting in cumulative vorticity stretching as great as that from curvature vorticity, $(Rcfdh/hU^)~1$, assuming $U^=0.15\u2009m\u2009s\u22121$, and *R*_{c} = 10 km. A similar but slighter steepening occurs around the NEC. Bathymetric steepening also limits streamline shoaling around bathymetric turns, which adds confidence to the scaling analysis since shoaling kinematically reduces streamline curvature.

Furthermore, if the (prograde) flow does “begin” to separate rather than turn around the cape, as a parcel travels offshore additional vortex stretching occurs. The increased vorticity may be expressed as increased cyclonic path curvature, steering the parcel further away from the upstream isobath. That can result in a positive feedback, e.g., by creating more positive curvature vorticity (by vortex stretching), and further angular separation from the continental slope. Note that floats trajectories separating around FC do tend to turn cyclonically offshore (Fig. 2), as do mean separating streamlines (Fig. 5).

To summarize, the scaling analysis suggests that partial separation (loss of outer streamlines to offshore) of a prograde current is a plausible outcome where significant downstream bathymetric steepening occurs, and especially where it is accompanied by anticyclonic bathymetric turning. The dependence on curvature radius is particularly simple to express—an order one flux reduction may result for curvature radius *R*_{c} < *W*. Both conditions are met at the FC leakiness hot spots (section 3a).

As an additional but still preliminary consistency check, we compare these conditions (downstream steepening and *R*_{c} < *W*) with the conditions at several locations of separation of other prograde currents: the western boundary current flowing around the southern tip of Greenland (Holliday et al. 2009), the Mediterranean Overflow current propagating around the Iberian Peninsula (McDowell and Rossby 1978; McWilliams 1985; Bower et al. 1997),^{11} and the California Undercurrent at the mouth of Monterrey Bay (Molemaker et al. 2015). The width of these currents is *O*(150, 50, 20) km, respectively, while the capes they traverse have *R*_{c} = *O*(10) km. Furthermore, steepening occurs as well on the upstream side of these capes. It is difficult to determine the relative contribution of steepening to vorticity stretching without knowledge of trajectories or mean streamlines, but the relative contraction of cross-isobaths distance at these capes is a large fraction, as in the FC separation locations.

Several assumptions and idealizations were made in deriving this scaling that remain to be tested:

The assumption of mean PV conservation along mean streamlines is only qualitatively motivated by the modest cumulative effect of eddy terms in the GB

_{B}PV budget.If the current width decreases downstream (to a value

*W*_{2}), then the magnitude of Δ*U*is overestimated in (8) by a factor ~*W*/*W*_{2}. However, the magnitude of Δ*U*estimated from our scaling at FC is such that even, e.g., a factor of two width decrease is relatively minor.^{12}If separation does occur, reattachment cannot be excluded within the analysis. It could only be suggested that lack of reattachment is likely for a prograde current separating from a large-angle bathymetric bend, due to cyclonic turning past separation resulting from additional vortex stretching.

In light of these assumptions and simplifications, the analysis needs to be further refined and tested in dedicated and controlled (e.g., numerical) experiments.

### c. Sensitivity to model circulation and resolution

The scaling analysis also suggests that inertial separation at a bathymetric turn should depend mostly on the local conditions: radius of curvature, bottom steepness changes, current width, and speed. Two implications may be that 1) the leakiness at FC should be largely insensitive to external variations in the Nfl circulation pattern and 2) as long as numerical model resolution is fine enough that, e.g., bathymetric curvature radii are similar to or smaller than model DWBC width, separation should still occur to some degree.

The observed leakiness patterns are reproduced well in the FLAME model employed by Bower et al. (2011), which has a coarser resolution, ≈6.5 km, despite water mass biases ( appendix B) in the Nfl Basin. Thus, the leakiness may not be strongly dependent on the detailed structure of the DWBC and surrounding currents, or on good resolution of baroclinic instabilities at the DWBC boundary (where the Rossby radius is ≈10 km). Indeed, with further decrease in model resolution, at least up to 0.5°, leakiness around FC and interior pathways still appear, but seem to gradually change and eventually severely deteriorate relative to observations (Gary et al. 2011; Spence et al. 2012).

## 5. Summary and conclusions

### a. Phenomenology

Using two observational float datasets and a realistic, high-resolution, numerical model, we demonstrate that within the Newfoundland (Nfl) Basin, the DWBC has a few well-defined geographical hot spots of maximal Lagrangian leakiness (Figs. 2 and 3). At the leakiness hot spots, local maxima of time-mean Lagrangian and Eulerian offshore velocities occur in the numerical model, while Lagrangian EKE is minimal (Figs. 3–6). These hot spots are further characterized by convex curvature and/or downstream steepening isobaths (Fig. 9). The localized, and time-mean nature of the leakiness, and its apparent correlation with bathymetric variations, suggests that it occurs largely via an inertial separation mechanism (mechanism 2 in section 1). This contrasts with previous hypotheses that suggested that the DWBC leakiness was due to interaction with NAC eddies (mechanism 1 in section 1).

The Eulerian mean circulation is examined within potential density layers, revealing that mean DWBC streamlines separate offshore at the identified Lagrangian leakiness hot spots (Fig. 5). Following separation, the streamlines revolve around deep cyclonic recirculations that reside between the DWBC and NAC. The Eulerian mean and Lagrangian mean DWBC velocities are very similar in the region. Consistently, the thickness-weighted average (TWA) flow is almost identical to the Eulerian mean flow, which means that the rectified eddy mean flow is negligible. Cluster analysis (supplemental material) supports these conclusions as well.

The mean offshore flow is associated with cumulative downstream reduction in DWBC mass flux (Fig. 4). Thus, we distinguish uncompensated leakiness from compensated leakiness: the former (latter) is associated with a net (zero) loss of material flux to offshore. The time-mean flow only contributes uncompensated leakiness, since no mean streamlines appear to join the DWBC from the interior (Fig. 5). The eddy component (with respect to the time mean) can contribute to uncompensated leakiness only via rectified eddy transport, which is found to be negligible compared to the Eulerian mean circulation around most of FC. The eddies may contribute substantially to compensated leakiness, but this has not been examined in this study. The model DWBC volume flux decreases by 13–16 Sv around FC, within LSW and deeper waters. The reduction primarily takes place at the Lagrangian leakiness hot spots identified in this study. The result is generally consistent with observational estimates (M14) showing ≈15 Sv loss between FC and east GB.

### b. Dynamics

The dynamics of separation are addressed from a (TWA) PV perspective (Fig. 7) as well as in terms of energetic transformations (Fig. 8). We find that mean separation deforms the PV contours offshore at the leakiness hot spots, which is consistent with inertial separation. Indeed, the cumulative change in mean PV along mean separating streamlines is modest (≈10%–20%) relative to the DWBC–NAC contrast. However, the mean PV advection is found to be balanced to first order by eddy PV flux divergence, indicating that eddies play a role in guiding the separated mean flow offshore. We therefore examined energy conversion processes in the region (Fig. 8). We find that the separation of mean DWBC streamlines is not directly forced by conversion of EKE to MKE (RSW_{e2m}). In fact RSW_{e2m} decreases and becomes negative prior to separation. This is consistent with the low magnitude of the eddy-rectified flow relative to the mean flow; as well as with cluster analysis (supplemental material), which shows that the separation of streamlines is statistically typical.

Outside of the separation areas, patches of positive PEC_{e} are collocated with positive RSW_{e2m}, which may be interpreted as a forcing of the mean flow by eddies spawned locally from baroclinic instability. This seems to occur in the recirculations, as well as in the DWBC itself (except at the separation areas). Rectified eddy mean flow is indeed toward offshore and significant at SF, but not elsewhere around FC. Our Lagrangian experiments using the time-mean model flow field (Exp3dMean) highlight another role of eddies in DWBC leakiness: this experiment exhibits ≈50% less leakiness compared with Lagrangian experiments using the time-dependent velocity (Exp3d). We attribute this difference to the eddies chaotically advecting Lagrangian particles from the upper continental slope toward the leakiness hot spots.

In contrast with previous hypotheses, our findings are in line with the main fraction of *mean uncompensated* leakiness occurring by inertial separation. Leakiness hot spots and mean streamline separation are localized at areas of convex and/or steepening bathymetry, where inertial separation may be expected. Furthermore, cumulative leakiness is demonstrated to be a persistent and typical occurrence, rather than eddying or intermittent. Along these mean separating streamlines, eddy PV flux divergence does not induce a dramatic change in mean PV, in support of inertial separation. The separation process is likely inviscid, since nonconservative terms have a small role in the TWA PV balance. Finally, past separation, Lagrangian trajectories as well as mean streamlines tend to turn cyclonically, which is consistent with vortex stretching in inertial motion into deeper water.

Previous theoretical frameworks determining conditions for inertial separation are not suitable for treating the DWBC conditions near FC. This is partially due their focus on buoyant rather than deep boundary currents. Additionally, the semigeostrophic approximation (made in these studies) is violated in areas of high curvature of the slope (section 4). Instead, a scaling analysis is presented (section 4b) for the downstream evolution of a boundary current due to bathymetric variations. The result suggests that a steady and continuous DWBC flow around the convex corners of FC is unlikely. A significant reduction in flux (e.g., partial separation) is a plausible outcome, due to influence of bathymetric curvature and steepening. Several assumptions made in the scaling analysis cannot be validated in the present study, and they require detailed examination in dedicated numerical experiments.

### c. Outlook

We note several caveats of the present investigation (also see appendix B):

The numerical model configuration developed and presented here suffers from water mass biases that make detailed comparisons with observations delicate at times, although mean circulation features and their variability appear to agree favorably with observations ( appendix B). Similar water mass biases plague numerical models of the area, and have been partially resolved in some studies using relaxation of water properties toward climatology (e.g., Tréguier et al. 2005; Rattan et al. 2010), a method not without drawbacks for dynamical analysis.

Likewise, total model DWBC transport east of FC is anomalously high in comparison with observations ( appendix B). This may have an impact on the leakiness process. For example a faster current may be more likely to inertially separate. However, in the supplemental material we show that the model DWBC transport is in good agreement with observations elsewhere in several other locations in Newfoundland including along FC, and that the anomaly is likely related to the cyclonic recirculations east of FC rather than to the DWBC itself.

While the model output frequency of two days is likely sufficient to resolve mesoscale processes, it may not provide sufficient representation of the submesoscale. For that reason leakiness by SCV formation and escape (mechanism 3 in section 1) is not addressed. Indeed, Bower et al. (2013) found that several ExPath floats were trapped in SCVs at or near their leakage from the DWBC around the GB southern tip.

On a related note, while vertical resolution is at a relatively high present standard (section 2a), it is not sufficient to resolve bottom boundary layer processes in deeper regions.

We note that the interpretation used here of binned and conditionally sampled Lagrangian velocities as the Lagrangian mean velocity (section 3a) is only approximately representative of true Lagrangian mean velocity. However, the low amplitude of rectified eddy flow as calculated independently of the defined Lagrangian mean (section 3b), corroborates that the deviation of Lagrangian mean from the Eulerian mean flow is small.

The results of this study suggest that the leakiness and separation mechanism depend strongly on the bathymetric environment of the current. Therefore, future work should examine the circulation in idealized scenarios where a DWBC-like current traverses a region of bathymetry resembling FC and GB. Within a simplified setting the dynamical mechanisms can be better isolated in experiments where factors such as bathymetry and the presence of a NAC-like countercurrent can be varied. Additionally, the geographical distribution of cumulative (uncompensated) leakiness was evaluated in our model, inspired by observational estimates (M14; Biló and Johns 2018). While they are consistent in terms of total flux, the observational record is not yet extensive enough to test the model distribution in detail.

Comparison of model particle trajectories transported by time-averaged versus unaveraged currents (section 3c) suggests that chaotic advection significantly increases the offshore leakiness of particles, including at the mean leakiness hot spots. We do not distinguish quantitatively the roles of pure eddy variability and of eddy interaction with spatial gradients in mean flow (i.e., chaotic advection). While several metrics were previously suggested to evaluate the relevance of chaotic advection in particular scenarios (e.g., Shepherd et al. 2000; Brett et al. 2019), it remains challenging to do so locally in a realistic flow such as examined here. Hence, we do not attempt in the present study to determine quantitatively the enhancement of leakiness by chaotic advection.

This study has concentrated on the mechanisms of DWBC leakiness in the Nfl Basin. Previous studies had a greater focus on characterization of the interior pathways that follow from the subpolar to the subtropical region. It has also been shown previously that most leaked particles recirculate in the Nfl Basin for years (Bower et al. 2009; Gary et al. 2011; Lozier et al. 2013). In this regard, the robustness of cyclonic mesoscale recirculations demonstrated in the present model also merits further study. Their relation to the larger-scale interior pathways and recirculation is also of interest. Furthermore, it remains to be determined if and how diapycnal mixing and water mass transformations are associated with the leakiness process or with the long recirculation period water parcels spend within the Nfl Basin.

## Acknowledgments

This material is based in part upon work supported by the National Science Foundation with Grants OCE-1538702, OCE-1751386, and OCE-1419450. The authors thank Dr. Amy Bower for providing the Export Pathways dataset, and Dr. Christian Mertens for providing the processed shipboard ADCP observations at 47°N, east of Flemish Cap. The SLA altimeter products were produced by Ssalto/Duacs and distributed by AVISO, with support from CNES (http://www.aviso.altimetry.fr/duacs/). The Argo data were collected and made freely available by the International Argo Program and the national programs that contribute to it (http://www.argo.ucsd.edu, http://argo.jcommops.org). The Argo Program is part of the Global Ocean Observing System. Argo float data can be downloaded from http://doi.org/10.17882/42182.

### APPENDIX A

#### Terms and Acronyms

Terms, acronyms, and symbols used often in the text are contained in Table A1.

### APPENDIX B

#### Numerical Model Validation

In this section we describe model validation against observations and discuss possible caveats in model setup. We begin with examining sea surface height (SSH), because it determines surface geostrophic velocity. We compare model SSH to the measurements of absolute dynamic topography from satellite altimetry. Model SSH is averaged over model years 9–16. The observational product we use is the DUACS L4 merged reprocessed product (Pujol et al. 2016), with 1/4° grid resolution and product samples spaced 1 day apart, with data from 1993 to 2017. The absolute dynamic topography to model SSH comparison is shown in Figs. B1a and B1b. Some differences in mean SSH and EKE are to be expected due to differences in averaging periods. There is general agreement in SSH patterns and amplitudes of the main circulation features, including the mean paths of the Gulf Stream and the Labrador Current, and the standing meanders of the NAC, including the Mann eddy.

Geostrophic surface $EKE=\u2061(1/2)u\u2032g2$ is compared between the model and the altimetric observations over the same periods as for the mean SSH. The model geostrophic component of surface eddy velocity $u\u2032g$ is calculated from eddy SSH (Vallis 2017). The observed eddy velocity is an available variable within the DUACS product. The model (observed) eddy component is defined as the instantaneous deviation from the time-mean SSH. Within the area shown in Figs. B1c and B1d, the model EKE is higher on average by a factor of ≈5. Higher EKE is generally to be expected in the model because its grid resolution is about 10 times higher compared with the altimetric product grid resolution, and because *R*_{d} in this region is close to or lower than the altimetric product grid resolution. Low-pass filtering of the model output shows that the unresolved scales likely account for the majority of the EKE difference (supplemental material, where effective resolution is taken into account). In addition, the spatial patterns of model and observed EKE are generally in good agreement. Both peak along the trajectories of the Gulf Stream and NAC. The model EKE also has a local peak of EKE along the 1 km isobath in the Labrador Sea and Nfl Basin. The peak is related to the Labrador Current, the inshore and upper ocean component of the western boundary current in the Subpolar North Atlantic, the deep component being the DWBC. The absence of Labrador Current signature in the observed EKE is again likely due to the coarser resolution.

We compare the depth and cross-stream structure and amplitude of the DWBC east of FC at 47°N (Fig. B2) with the observations of M14 reproduced in Fig. B2b. The observational estimate was obtained by averaging over six individual vessel ADCP cross-DWBC sections, taken at various dates between April and August at six different years (M14). The ROMS data presented in Fig. B2a are an average over model years 9–16. However, model averages over single years are generally quite similar (e.g., for DWBC width variation). Considering the very different averaging details, the spatial patterns are visually similar between the model and observational estimates. Above the continental slope there is an intensified, quasi-barotropic DWBC core, while over the continental rise a bottom-intensified DWBC core is present. The northward flow to the east is related to the NAC and is similar in its structure between the model and observational estimates as well. The multiple surface-intensified cores present in the observational northward flow may be smeared out in the longer model time average.

The maximal DWBC velocity magnitude is just under 0.3 m s^{−1} in both model and observational averages within both current cores, except within very limited regions in the observational estimate where the magnitude exceeds 0.3 m s^{−1}. The total width of the DWBC compares well with the observations. Here an operational definition of the current width is taken as the distance between the 0.05 m s^{−1} velocity contours near the bottom, west of the western (continental slope) core, and east of the eastern (continental rise) core. With this definition, the model (observed) width is 156 ± 13 (153) km. The model width error estimate quoted here is the standard deviation in annual-mean widths between years 9–16. We did not obtain the results for the individual (six) observational cruises on which Fig. B2b is based. However, based on Fig. 4 in M14, we estimate the observed width std at *O*(50) km, across the six cruises.

Despite the agreement with observations in patterns, widths, and maximal velocities along this section, the model flow is seen to be more barotropic than the observational estimate at the 47°N section, and therefore carries a higher total volume transport. In what follows, GB_{B} transport uncertainties are calculated from interannual variations of annual mean flow, unless otherwise stated. Mean transport is calculated as the total southward transport west of 41°W of the averaged velocity across the section. Note that this straightforward Eulerian mean transport definition is different than that of M14 (supplemental material). From the mean section of M14 observations, we calculate a depth-integrated transport estimate of 30.8 Sv. The DWBC transport in the model along this section is 58.5 ± 29.8 Sv, where the standard deviation is over all model (2-day) output samples, while interannual standard deviation in annual mean transports is 4.5 Sv. The difference between the model and observed transport sample means is statistically significant (supplemental material). The model transport estimate, as well as the stronger barotropic tendency relative to the observations along the section, are similar to the results of the VIKING20 numerical model employed by M14, 60.3 ± 23.6 Sv. However, further validation deferred to the supplemental material shows that the GB_{B} model DWBC transport in other sections is very similar to observations, and a cause for this difference is suggested.

We also compared model Eulerian EKE with observations. Fischer et al. (2018b) have gridded velocity data from Argo floats (Lebedev et al. 2007) at 1500 m depth around FC as well as farther north. We use their Gaussian-interpolated product (Fischer et al. 2018a), with gridcell size (generally not equivalent to resolution) of 1/4° (1/2°) latitude (longitude). Note Fig. 5b in Fischer et al. (2018b) is somewhat saturated in some areas around FC. We find that around east and south FC, within the DWBC and NAC, the model EKE is of similar magnitude or higher (by up to a factor ~3) than the Fischer et al. (2018b) gridded EKE. As in the altimetric observations (see above) this is likely related to the coarser observational product not fully resolving smaller-scale fluctuations. Fischer et al. (2018b) also provide EKE values at two moorings (K18 and B227) within the DWBC around FC. Mooring B227 is near the M14 section. We find that model EKE at these mooring locations is only ~20% higher than observed, and within the uncertainty range (the difference being equal to about one standard deviation of model EKE values).

The model suffers from a bias in the middepth density field. As seen in Fig. B2, the *σ*_{θ} = 27.8 kg m^{−3} isopycnal is 400–700 m too deep in the model. The model density bias is mostly related to (not shown) a salinity bias. A salinity-related density bias, especially at middepth, is very common in subpolar North Atlantic numerical models. See, for example, Figs. 3 and 6 in Bower et al. (2011) in comparison with Fig. B2, as well as Fig. 2 in Handmann et al. (2018). This common problem was previously attributed (Tréguier et al. 2005; Rattan et al. 2010) largely to salt transport biases appearing in model boundary currents. Typically, nudging model salinity to climatological values is required, although not always sufficient, to reduce or eliminate the bias in present models. A disadvantage associated with a nudging procedure may be reduction in frontal features and sharpness of boundary currents in high-resolution models, because the resolution of climatological datasets is generally coarser. We therefore did not apply such a nudging procedure. In our model the bias gradually appears during spinup and appears to be fully developed by year 9, without further increase in the bias amplitude in the following years. It is difficult to determine with certainty to what degree our key results are affected by the water mass bias. We expect, however, that such effects should manifest mainly indirectly, through the effects on the mean circulation and on EKE. The good agreement of leakiness and recirculation patterns with other models, and with observations (here, and in section 3a), is encouraging in this regard, as is the comparison with EKE observations (above). It appears, however, that DWBC flow east of Flemish Cap has a stronger barotropic component than observations suggest. The implications of this possible bias are discussed in section 5c.

There are additional caveats concerning the temporal extent of the surface and horizontal boundary fields used to determine the model boundary conditions. These fields only have a 4-yr length, corresponding to a 2001–05 atmospheric state, and are recycled after the first four model years (section 2a). Since the domain (open) boundaries are very far (over 500 km) from the analyzed area, the transient effects of the recycling method are likely very limited. Indeed, we do not observe any significant changes at the 4-yr period (e.g., in mean kinetic or potential energy), other than the seasonal cycle similar to that observed in other years. Furthermore, the analysis presented in section 3c confirms that rare events are not important for either the mean or eddy components of offshore flow. However, years 2001–05 cover only negative to moderate North Atlantic Oscillation index values, and therefore the model boundary forcing is likely not representative of the full range of DWBC variability. Interannual and decadal variability in atmospheric forcing, including that due to the North Atlantic Oscillation, influences the depth of deep convection in the Labrador Sea, and hence the variability in LSW thermohaline properties (Yashayaev and Loder 2016) as well as DWBC transport (Zantopp et al. 2017).

Finally, we qualitatively compare in Fig. B3 pathways of (3D) Lagrangian floats in the model (Exp3d, section 2c) to the (isobaric) ExPath floats (Bower et al. 2011). At each deployment depth (700 and 1500 m), a batch of 30 particles are randomly selected and their trajectories extended to a 2 year duration. The full trajectories of these floats and ExPath floats are displayed in Fig. B3. The transport patterns are generally similar to those sampled by the ExPath floats: the majority of particles were caught in recirculations within the Newfoundland Basin. A smaller fraction traveled south in the interior of the ocean. Some particles crossed the Mid-Atlantic Ridge eastward at the Charlie–Gibbs Fracture Zone. Only a few particles traveled within the DWBC continuously past the GB, although more model particles did so compared with ExPath floats. This likely due to fewer model particles traveling through Flemish Pass, which we speculate happens either due to model velocity output frequency not being high enough, or due to water mass biases. A second bias appears in the model in that more particles appear to cross north to the subpolar area compared to the number for ExPath floats. Some of these differences from the ExPath floats manifest similarly in the 3D float trajectories of Bower et al. (2011, their Fig. 7a).

### APPENDIX C

#### Correlations Between Offshore Velocity and Bathymetric Variation

We present correlations between offshore velocity and variables related to bathymetric variation along the 3 km isobath. Offshore velocity is averaged over depths greater than 500 m. The bathymetric variables examined are curvature, steepness, and steepening. The latter is defined as the change in steepness with along-isobath distance. Offshore velocity and the former two bathymetric variables are displayed in Fig. 9. The correlation data are summarized in Table C1.

## REFERENCES

*Deep-Sea Res. I*

*J. Fluid Mech.*

*Mar. Geod.*

*Geophys. Res. Lett.*

*Deep-Sea Res.*

*J. Geophys. Res.*

*J. Phys. Oceanogr.*

*Nature*

*Deep-Sea Res. II*

*Deep-Sea Res. II*

*Nonlinear Processes Geophys.*

*Rev. Geophys.*

*Ocean Circulation and Climate: A 21st Century Perspective*

*J. Geophys. Res.*

*J. Mar. Res.*

*J. Phys. Oceanogr.*

*Ocean Sci.*

*J. Geophys. Res.*

*Ocean Modell.*

**50**,

*J. Geophys. Res.*

*J. Phys. Oceanogr.*

*J. Phys. Oceanogr.*

*J. Phys. Oceanogr.*

*Ocean Sci.*

*Deep-Sea Res. II*

*Geophys. Res. Lett.*

*Geophys. Res. Lett.*

*J. Geophys. Res.*

*Modeling Atmospheric and Oceanic Flows: Insights from Laboratory Experiments and Numerical Simulations, Geophys. Monogr.*

*J. Geophys. Res. Oceans*

*Dyn. Atmos. Oceans*

*Rev. Geophys.*

*J. Phys. Oceanogr.*

*Amer. J. Phys.*

*J. Phys. Oceanogr.*

*J. Atmos. Sci.*

*Geophys. Res. Lett.*

*J. Phys. Oceanogr.*

*Rev. Geophys.*

*Nature*

*Deep-Sea Res. I*

*J. Geophys. Res. Oceans*

*Ocean Sci.*

*Ocean Modell.*

*Science*

*Annu. Rev. Mar. Sci.*

*Deep-Sea Res. II*

*J. Phys. Oceanogr.*

*Ocean Modell.*

*Ocean Modell.*

*Ocean Modell.*

*Science*

*Rev. Geophys.*

*J. Atmos. Sci.*

*J. Geophys. Res. Oceans*

*J. Phys. Oceanogr.*

*J. Phys. Oceanogr.*

*J. Phys. Oceanogr.*

*J. Mar. Res.*

*Ocean Sci.*

*Ocean Modell.*

*J. Phys. Oceanogr.*

*J. Phys. Oceanogr.*

*J. Fluid Mech.*

*Nat. Climate Change*

*Rev. Geophys.*

*J. Phys. Oceanogr.*

*Bull. Amer. Meteor. Soc.*

*J. Phys. Oceanogr.*

*Ocean Modell.*

*Ocean Modell.*

*J. Geophys. Res.*

*J. Phys. Oceanogr.*

*J. Fluid Mech.*

*J. Phys. Oceanogr.*

*Bull. Amer. Meteor. Soc.*

*Deep-Sea Res.*

_{2}, and net sea–air CO

_{2}flux over the global oceans

*Deep-Sea Res. II*

*Descriptive Physical Oceanography: An Introduction*

*J. Phys. Oceanogr.*

*J. Phys. Oceanogr.*

*Geophys. Res. Lett.*

*Atmospheric and Oceanic Fluid Dynamics*. 2nd ed

*J. Geophys. Res.*

*Statistical Methods in the Atmospheric Sciences*

*J. Phys. Oceanogr.*

*J. Geophys. Res. Oceans*

*J. Phys. Oceanogr.*

*J. Geophys. Res. Oceans*

## Footnotes

^{1}

A table of acronyms and terms commonly used in the text appears in appendix A.

^{2}

LSW, formed mainly in Labrador Sea deep convection events, comprises the NADW upper component, typically ≈400–2000 m (Yashayaev and Loder 2016; Bullister et al. 2013). The lower component is Overflow Water (Talley 2011).

^{3}

Farther downstream, around GB, several floats did come back into the DWBC (Bower et al. 2009).

^{4}

Leakiness of profiling floats in this region was investigated by Lavender et al. (2000) and Fischer and Schott (2002) as well.

^{5}

The cross-bathymetry component *υ*_{c} points toward deeper water, and the along-bathymetry component *υ*_{a} is defined to point to the right of *υ*_{c}, i.e., generally downstream for the DWBC.

^{6}

The same calculation applied to the 4 km isobath yields very similar results (e.g., ~15 Sv offshore flux at FC). The 3 km isobath is used here since it extends farther north past Orphan Knoll.

^{7}

Nonconservative terms are almost certainly even lower in magnitude than indicated by the pattern correlation result. That is because diagnostics are based on 2-day averaged output and higher frequency variability is unresolved, i.e., aliased.

^{8}

The RSW_{m2e} term (not shown) of the EKE equation is mostly similar in pattern and opposite in sign to RSW_{e2m} in the area.

^{9}

Laboratory experiments related to the same parameter regimes as these theoretical works were conducted by, e.g., (Whitehead and Miller 1979; Bormans and Garrett 1989).

^{10}

Offshore of the downstream stagnation streamline, the present analysis cannot determine the circulation pattern, since offshore streamlines do not necessarily originate upstream. Rather than a reversal or recirculation, a split in the current may emerge for example. That does not affect the result inshore, however.

^{11}

Note that leakiness of the Mediterranean Overflow current is at least in some cases associated with SCV formation and interactions (Bower et al. 1997).

^{12}

Additionally, even following separation at SEC (within GB_{B}) width decreases by only a small fraction locally (Figs. 5b,d).