Exceptional freshening and cooling in the eastern subpolar North Atlantic caused by reduced Labrador Sea surface heat loss

Observations of the eastern subpolar North Atlantic in the 2010s show exceptional freshening and cooling of the upper ocean, peaking in 2016 with the lowest salinities recorded for 120 years. Published theories for the mechanisms driving the freshening include: reduced transport of saltier, warmer surface waters northwards from the subtropics associated with reduced meridional overturning; shifts in the pathways of fresher, cooler surface water from the Labrador Sea driven by changing patterns of wind stress; and the eastward expansion of the subpolar gyre. Using output from a high-resolution hindcast model simulation, we propose that the primary cause of the exceptional freshening and cooling is reduced surface heat loss in the Labrador Sea. Tracking virtual fluid particles in the model backwards from the eastern subpolar North Atlantic between 1990 and 2020 shows the major cause of the freshening and cooling to be an increased outflow of relatively fresh and cold surface waters from the Labrador Sea; with a minor contribution from reduced transport of warmer, saltier surface water northward from the subtropics. The cooling, but not the freshening, produced by these changing proportions of waters of subpolar and subtropical origin is mitigated by reduced along-track heat loss to the atmosphere in the North Atlantic Current. We analyse modelled boundary exchanges and water mass transformation in the Labrador Sea to show that since 2000, while inflows of lighter surface waters remain steady, the increasing output of these waters is due to reduced surface heat loss in the Labrador Sea beginning in the early 2000s. Tracking particles further upstream reveals that the primary source of the increased volume of lighter water transported out of the Labrador Sea is increased recirculation of water, and therefore longer residence times, in the upper 500–1000 m of the subpolar gyre.

1 Introduction

Upper ocean temperature and salinity and their variability in the eastern subpolar North Atlantic are governed by the convergence and divergence of heat and salt and local ocean–atmosphere exchanges. The convergence and divergence are associated with: advection by midlatitude ocean circulation (Piecuch et al., 2017) and the strength of the Atlantic meridional overturning circulation (AMOC) bringing subtropical-origin water northwards (Robson et al., 2016; Foukal and Lozier, 2018; Bryden et al., 2020); northward transport of heat to high latitudes (Keil et al., 2020); shoaling isotherms and isohalines bringing colder, fresher water closer to the surface (Josey et al., 2018); and eastward expansion and contraction of the subpolar gyre (SPG) (Bersch, 2002; Hátún et al., 2005; Bersch et al., 2007; Sarafanov, 2009; Koul et al., 2020), with the accompanying movement of subpolar fronts and changes in the North Atlantic Current (NAC) pathways. The ocean circulation variability is in turn associated with local and remote changes in atmospheric forcing, including surface heat fluxes (Josey et al., 2018), shortwave cloud feedbacks (Keil et al., 2020), the leading modes of atmospheric variability in the North Atlantic (Häkkinen et al., 2011; Koul et al., 2020), and wind stress curl patterns (Häkkinen et al., 2011; Chafik et al., 2019; Holliday et al., 2020); as well as upstream freshwater input (Peterson et al., 2006) and internal ocean modes (Johnson et al., 2019).

Multiple observations show exceptional freshening and cooling of the upper 500–1000 m of the eastern subpolar gyre (Josey et al., 2018; Holliday et al., 2020) starting around 2012 and running to the end of the 2010s (Figs. 1 and 2). At its peak, around 2016, this represented the strongest freshening event for 120 years (Holliday et al., 2020). Model- and observation-based analyses have shown decadal variability in heat and salt content to be primarily driven by variability in northward transport of warm, salty waters from the subtropics (Burkholder and Lozier, 2014; Desbruyères et al., 2015; Robson et al., 2016; Foukal and Lozier, 2018; Desbruyères et al., 2021). This provides decadal-scale preconditioning of the upper ocean by a reduced AMOC which, when coupled with intense local winter mixing events bringing colder fresher water to the surface (Josey et al., 2018; Bryden et al., 2020), could produce sufficient cooling and freshening to explain the observed recent exceptional event. Alternatively, increased volumes of cold fresh water from the Labrador Current due to changing pathways of this flow driven by winter wind stress have also been invoked to explain the exceptionally low upper ocean salinity observed (Holliday et al., 2020). This more closely follows ideas that redistribution of heat and salt within the subpolar North Atlantic (SPNA) and variability in the intensity and extent of the SPG drive interannual to decadal variability in the eastern subpolar North Atlantic (Koul et al., 2020; Chafik et al., 2019; Kenigson and Timmermans, 2020). More recent observations in years since 2016 reveal a cooling-to-warming transition induced by increasing transport of warm subtropical waters in the NAC (Desbruyères et al., 2021).

Our aims here are to quantify the relative contributions of annual to decadal variability in ocean circulation, upstream water properties, and local atmosphere–ocean fluxes to the recent exceptional freshening and cooling event of 2012 to 2016 in the eastern subpolar North Atlantic. Using the knowledge gained, we then look at possible mechanisms driving the changes. We use outputs from a hindcast with a high-resolution, eddy-rich ocean model (VIKING20X) (Biastoch et al., 2021), which reproduces the timing and spatial scales and patterns of the exceptional freshening. Using Lagrangian model analysis software (OceanParcels) (Lange and Sebille, 2017; Delandmeter and Van Sebille, 2019), particles, each representing a fixed volume transport, are tracked backwards in time from positions along the upper 500–1000 m of a vertical section spanning the eastern subpolar North Atlantic. The resulting time series of variability in the volume transports from the water sources, thermohaline properties at source, pathways, transit times, and along-track property transformation of the modelled water parcels are used to quantify the contributions to the observed freshening and cooling. We couple this Lagrangian analysis with an examination of Eulerian model time series and water mass transformation analysis in the model (Walin, 1982; Tziperman, 1986; Speer and Tziperman, 1992; Nurser et al., 1999) to produce a more complete theory of mechanisms driving the exceptional freshening and cooling.

Figure 1 The 2014–16 mean (a, c) salinity and (b, d) temperature anomalies, relative to the 2004–2019 mean, in the upper 500 m from the (a, b) EN4 dataset and (c, d) VIKING20X model. The OSNAP (Overturning in the Subpolar North Atlantic Programme) observational section (e, f and red lines in a–d), occupied since 2014 and used here to spawn tracked particles in VIKING20X, traverses the fresh patch in the Iceland and Irminger basins. In 2014–2016 there is a widespread region of fresher than average water in the upper 500 m of the eastern subpolar gyre, with a “compensating” region of raised salinity on the North West Atlantic Continental Shelf and Slope region (NWACSS, 40–50^∘ N, 50–70^∘ W, which includes offshore Slope Sea deep water north of the Gulf Stream). The blue polygon in the Labrador Sea outflow region (a–d) shows the region used for comparison with observations of isopycnal depths (Sect. 6.3). Vertical red lines in (e, f) show the positions of the time series in Fig. 6.

2 Methods

2.1 Hydrodynamic model

We make use of the eddy-rich, nested ocean–sea-ice model configuration VIKING20X (Biastoch et al., 2021), based on the Nucleus for European Modelling of the Ocean code (NEMO, version 3.6) (Madec et al., 2017) and the Louvain la Neuve Ice Model (LIM2) (Fichefet and Maqueda, 1997; Goosse and Fichefet, 1999). In the vertical, VIKING20X uses 46 geopotential z-levels with layer thicknesses gradually increasing from 6 m at the surface to 250 m in the deepest layers. The bottom topography is represented by partially filled cells, allowing for an improved representation of the bathymetry and its slopes (Barnier et al., 2006). In the horizontal, VIKING20X has a tripolar grid with 0.25^∘ global resolution, which is refined in the Atlantic Ocean to 0.05^∘, yielding an effective grid spacing of 3–4 km in the subpolar North Atlantic. The run used here, VIKING20X-JRA-short, is an experiment forced from 1980 to 2019 by the JRA55-do forcing (version 1.4) (Tsujino et al., 2018), which is branched from a previous hindcast experiment forced by the CORE dataset (version 2) (Griffies et al., 2009; Large and Yeager, 2009) covering the period 1958 to 2009. Three-dimensional model output fields are saved as 5 d means, which are used here for the offline Lagrangian particle tracking and water mass transformation analysis. Hindcasts of the past 50–60 years in this eddy-rich configuration realistically simulate the large-scale horizontal circulation, including the AMOC, the distribution of the mesoscale, overflow, and convective processes, and the representation of regional current systems in the North and South Atlantic (Biastoch et al., 2021; Rühs et al., 2021).

Figure 2 Model (plain solid lines) and observational (EN4, dashed lines; OSNAP, bold lines from 2014) time series of (a) salinity and (b) temperature in three zones – in the Irminger Sea, Iceland Basin, and Rockall Trough – along the eastern section of the OSNAP (Overturning in the Subpolar North Atlantic Programme) line, OSNAP_E (see Fig. 1 for the positions of the zones). All time series are monthly means and are averaged over the upper 500 m.

2.2 Lagrangian particle tracking

Particles are tracked in the VIKING20X model using the Parcels Lagrangian framework v2.2.2 (Lange and Sebille, 2017; Delandmeter and Van Sebille, 2019). About 70 000 particles are released every 5 d for 30 years, 1990–2019, throughout the top 1000 m of the OSNAP (Overturning in the Subpolar North Atlantic Programme) line, OSNAP_E (Fig. 3), east of 37^∘ W (Fig. 3). The western extreme was set at 37^∘ W and the eastern extreme at the Scottish coast to sample all the northward upper ocean flow in the eastern SPG while excluding the southward upper ocean flow east of Greenland. Each particle is tracked backwards in time for 10 years or until it leaves the VIKING20X domain or reaches 20^∘ N. This is a total of 1.5 billion particle years. The 1990 to 2019 period was chosen to keep the experiment within the post-1980 hindcast period of the VIKING20X-JRA-short simulation and to compare and contrast the early 1990s and late 2010s subpolar gyre index maxima (Biastoch et al., 2021).

Particles are advected by the VIKING20X 5 d mean, 3-dimensional velocity fields north of 20^∘ N using a fourth-order Runge–Kutta time step scheme with a 10 min time step. No additional diffusion is used. Particle position, velocity, temperature, salinity, and model mixed layer depth are recorded every 5 d along-track.

The initial particle distribution along the 2-dimensional OSNAP_E section is random, with local particle density scaled by the magnitude of model velocity normal to OSNAP_E at release time. Each particle therefore represents the same volume transport (here 0.001797 Sv), with positive transport northward and negative southward. Particles are assumed to maintain these along-track transport values throughout the track, as for streamtubes in steady flow (van Sebille et al., 2018). This assumption is only formally valid when using analytical advection methods for steady flows, but Schmidt et al. (2021) have shown that the numerical advection method in Parcels can be reliably used for volume transport estimations.

We attempt to estimate the random errors associated with the particle sampling strategy and the effect of unsteady flow on the streamtube assumption. As we increase particle numbers, the metrics we discuss will converge towards “final” values. Using random release positions allows us to split the particle set into a number of subsets (32 here), each of which is still random but with each particle representing a (32 times) larger transport. The spread of the results from these 32 subsets gives us a measure of errors due to the finite sample size. Particle release numbers were set high enough to keep sampling errors small compared to the variability being investigated, as is the current standard in following water masses using Lagrangian tools (Schmidt et al., 2021). Possible systematic biases present in all 32 subsets and associated with the model flow fields or the streamtube assumption as well as errors due to diffusive and eddying processes not being represented in the tracks will not be quantified by this method.

Figure 3 Example particle paths showing the identification of particle sources and pathways. The black line is the section of OSNAP_E where particles are released for backward tracking. The coloured straight lines delimit the source regions and pathways: orange – Labrador Sea, green – Gulf Stream, dark orange – the loop pathway, pink – Hudson Bay, brown – Davis Strait, and yellow – Greenland Sea. Panels (a–c) show the main source regions: (a) the Labrador Sea, (b) the Gulf Stream, (c) no source determined – “other” source. Panels (d, e) show further subsets of (a), highlighting different pathways from the Labrador Sea to OSNAP_E: (d) the direct path (tracks do not cross the dark orange line at 60^∘ W), (e) the loop path (tracks cross and recross the dark orange line). Panels (f–i) show an alternative subset of pathways of Labrador Sea origin (a), classifying particle tracks by source regions further upstream of the Labrador Sea: (f) Hudson Bay, (g) Davis Strait, (h) Greenland Sea, and (i) particles recirculating in the subpolar gyre.

Particle source regions – the Gulf Stream and Labrador Sea – were defined as the region last visited before particles arrived at OSNAP_E (Fig. 3). A small proportion of water (5 %–10 %, 2–3 Sv) has an undetermined source even after 10 years of tracking. These particles with an origin of “other” remained circulating within the eastern subpolar gyre for the full 10 years of tracking (apart from a very few which came from the Mediterranean Sea). The tracks of Labrador Sea origin were additionally classified by whether they passed through the western Slope Sea region between the Gulf Stream and the shelf slope (which we define as west of 60^∘ W, see New et al., 2021) – “loop” path – or took a more direct route – “direct” path – between leaving the Labrador Sea and arriving at OSNAP_E. We further subdivided Labrador Sea origin tracks by source further upstream (Fig. 3f–i): “Hudson Bay” describes particles which have entered the Labrador Sea from Hudson Bay; “Davis Strait” describes particles which entered the Labrador Sea southward through the Davis Strait; “Greenland Sea” describes particles which have crossed the Greenland–Scotland ridge southward, mostly through the Denmark Strait from the Greenland Sea, before travelling round the Labrador Sea; and “SPG” describes particles of Labrador Sea origin which have remained within the SPG for a complete circuit.

In the analysis we exclude any particles which pass north of OSNAP_E on the backward tracks between their release in OSNAP_E and leaving their source region (see, e.g. van Sebille et al., 2018). Here, we mainly present analyses based on a subset of particles which cross OSNAP_E in the surface 500 m where the freshening is concentrated. We refer to this section as OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. Two other sets of particles were explored – those crossing OSNAP_E in the top 1000 m or the upper limb of AMOC (σ₀<27.62 kg m⁻³). Results from these alternative release configurations are qualitatively the same as the main results; we present some quantitative comparisons between configurations in the text and tables.

Time series are primarily presented plotted against the time the particles cross OSNAP_E. It is also useful to present some time series plotted against the time the particles left the major source regions. By grouping particles by their source-leaving time, using the same 5 d intervals as for the particle releases, we obtain the 5 d mean transport that leaves the source regions of water and subsequently crosses OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. With a 10 year tracking duration, these transports measured at source will only be fully “saturated” between 1990 and 2009; outside this window, particle numbers drop off due to the transit times from OSNAP_E. To extend this window, we use the fractional distribution, g(t), and cumulative distribution, $h\left(t\right)={\int }_{\mathrm{0}}^{t}g\left(t^{\prime }\right)d{t}^{\prime }$ (where $\mathrm{0}<t<\mathrm{10}$ years and $\mathrm{0}<g\left(t\right),h\left(t\right)<\mathrm{1}$) of transit times to construct a time-dependent normalization factor, f(d), where d is the date in decimal years:

$$\begin{array}{}\text{(1)}& f\left(d\right)=\left\{\begin{array}{ll}\mathrm{1}-h(\mathrm{1990}-d)& \text{if}\mathrm{1980}\le d<\mathrm{1990},\\ \mathrm{1}& \text{if}\mathrm{1990}\le d<\mathrm{2010},\\ \mathrm{1}-h(d-\mathrm{2010})& \text{if}\mathrm{2010}\le d<\mathrm{2020}.\end{array}\right.\end{array}$$

Particle numbers leaving a source at any time are normalized by dividing by this factor. We cut off the extremes where fewer than 60 % of the transports would be sampled, giving a window between mid-1987 and mid-2017 (Fig. 4c).

Figure 4 Construction of the factor used to normalize the transports at source and select a window of maximum confidence. (a) Histogram of transit times between leaving the source region and crossing OSNAP_E. (b) Cumulative histogram showing the fraction of particles transiting in less than each number of years. (c) Particle number normalization factor. Numbers of particles leaving the source are saturated between 1990 and 2010; outside this range, numbers fall off at a rate determined by the transit time. Particle numbers leaving a source at any time are normalized by dividing by this factor. We cut off when fewer than 60 % of the particles are expected to be recorded (shaded zones).

Average temperature and salinity (TS) properties are calculated both as particles cross OSNAP_E and by averaging T and S from the time each particle leaves the source regions. These averages are performed for the combined sources and for separate sources and pathways. Note that the TS averages we present here are averaged over the volume being transported, rather than being area averages.

2.3 Water mass analysis

The rate of formation of a water mass between two isopycnals M(ρ) in some domain can be expressed as

$$\begin{array}{}\text{(2)}& M\left(\mathit{\rho }\right)=-{\displaystyle \frac{\partial F}{\partial \mathit{\rho }}}+{\displaystyle \frac{{\partial }^{\mathrm{2}}{D}_{\mathrm{diff}}}{\partial {\mathit{\rho }}^{\mathrm{2}}}},\end{array}$$

where F represents the the transformation driven by surface fluxes along the surface outcrop of the ρ isopycnal (Speer and Tziperman, 1992), and D_diff represents the total diapycnal diffusive density flux across the ρ isopycnal resulting from interior mixing (Walin, 1982; Speer and Tziperman, 1992; Nurser et al., 1999; Fox and Haines, 2003). Here F is given by

$$\begin{array}{}\text{(3)}& F=\underset{\mathrm{\Delta }\mathit{\rho }\to \mathrm{0}}{lim}{\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }\mathit{\rho }}}\underset{\mathrm{outcrop}}{\int }{\mathcal{D}}_{\mathrm{in}}dA,\end{array}$$

where 𝒟_in is the density influx per unit area and A is the surface area outcrop between isopycnals ρ and ρ+Δρ. The transformation driven by surface fluxes is divided into heat (F_H) and freshwater fluxes (F_F) (Tziperman, 1986):

$$\begin{array}{}\text{(4)}& F=F_{\text{H}}+F_{\text{F}},\end{array}$$

with

$$\begin{array}{}\text{(5a)}& {\displaystyle }{F}_{\text{H}}& {\displaystyle }=-\underset{\mathrm{\Delta }\mathit{\rho }\to \mathrm{0}}{lim}{\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }\mathit{\rho }}}\underset{\mathrm{outcrop}}{\int }{\displaystyle \frac{\mathit{\rho }\mathit{\alpha }\mathcal{H}}{C_{\text{p}}}}dA,\text{(5b)}& {\displaystyle }{F}_{\text{F}}& {\displaystyle }=\underset{\mathrm{\Delta }\mathit{\rho }\to \mathrm{0}}{lim}{\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }\mathit{\rho }}}\underset{\mathrm{outcrop}}{\int }\mathit{\rho }\mathit{\beta }S\mathcal{Q}dA+\underset{\mathrm{surface}{\mathit{\rho }}^{\prime }<\mathit{\rho }}{\int }\mathcal{Q}dA.\end{array}$$

Here ℋ is the downward surface heat flux, 𝒬 is the upward freshwater flux, α and β are the thermal expansion and saline contraction coefficients respectively, and C_p is the specific heat capacity.

Integrating Eq. (2) in density from the lightest ρ_min up to ρ, we have the volume transformation budget

$$\begin{array}{}\text{(6)}& \begin{array}{rl}& \underset{\mathrm{volume}\mathrm{tendency}}{\underbrace{{\displaystyle \frac{\partial V\left(\mathit{\rho }\right)}{\partial t}}}}+\underset{\mathrm{outflow}}{\underbrace{\mathrm{\Psi }\left(\mathit{\rho }\right)}}\\ & =-\underset{\mathrm{surface}\mathrm{heat}}{\underbrace{F_{\text{H}}}}-\underset{\mathrm{surface}\mathrm{freshwater}}{\underbrace{F_{\text{F}}}}+\underset{\mathrm{diffusive}\mathrm{outflux}}{\underbrace{{\displaystyle \frac{\partial D_{\mathrm{diff}}}{\partial \mathit{\rho }}}}},\end{array}\end{array}$$

where V is the total volume of water with density<ρ, and Ψ is the volume flux of water with density<ρ out of the domain.

The volume tendency, outflow, and surface influx terms in Eq. (6) can be calculated from the VIKING20X 5 d mean model outputs. The diffusive flux – which we cannot calculate directly from the standard outputs – can then be deduced as the sum of the other four terms.

2.4 Observational data

We use openly available EN4 gridded ocean analyses (EN.4.2.2, bias corrections .g10, downloaded 12 January 2021) (Good et al., 2013; Gouretski and Reseghetti, 2010; Gouretski and Cheng, 2020) for the period with good Argo float coverage (after 2004), OSNAP section gridded data (Lozier et al., 2019; Li et al., 2021) for 2014 to 2018, and RAPID AMOC time series (Frajka-Williams et al., 2021) from 2004 to 2020. These data are used with no notable further analysis, and details of the collection, quality control, and analysis can be found in the relevant publications listed.

For comparisons of model results to observed isopycnal depths in the Labrador Sea, we use data from ship, mooring, and profiling float observations collected and analysed as part of the Deep-Ocean Observation and Research Synthesis (DOORS), a follow-up to the World Ocean Circulation Experiment (WOCE). All profiles from both vessel surveys and Argo floats operating in the Labrador Sea during 1990–2020 have been quality controlled through semi-automated and visual inspections, calibrated (in the case of Argo floats, ship–platform and inter-platform optimizing calibration was developed and used), and validated. The steps of subsequent vertical interpolation and data processing, including the evaluation and removal of the seasonal cycle, are discussed in https://www.dfo-mpo.gc.ca/csas-sccs/Publications/ResDocs-DocRech/2022/2022_039-eng.html, last access: 14 July 2022. The isopycnal depths were computed for selected profiles covering the Labrador Sea outflow region (blue polygon in Fig. 1a–d). Analysis was performed separately for profiles over the Labrador Slope (in water depths of 750 to 3275 m) and in deeper water in the central Labrador Sea. The resulting time series of isopycnal depths had the regular seasonal cycle removed, and residuals were low-pass filtered.

3 The 2012–2016 eastern subpolar North Atlantic freshening event

First, we briefly establish the spatial and temporal characteristics of the exceptional freshening and cooling event in the eastern subpolar North Atlantic beginning around 2012 (Figs. 1 and 2) and show that the VIKING20X model captures enough detail to be a useful tool to evaluate the possible causes (for a more detailed description of the observations, see Holliday et al., 2020). Spatial patterns of salinification/freshening and warming/cooling in the observations and model are very similar (Fig. 1a–d show 2014 to 2016 anomalies from the 2004 to 2019 means). Both the EN4 observations and the VIKING20X model show an extensive area of upper ocean freshening and cooling in the eastern subpolar North Atlantic in 2014 to 2016, with strong warming and salinification in the NWACSS (North West Atlantic Continental Shelf and Slope, 40–50^∘ N, 50–70^∘ W, which includes offshore Slope Sea deep water north of the Gulf Stream) region. There is also warming and salinification in the subtropical gyre (south of 35^∘ N), the southern Labrador Sea, and north of the Greenland–Scotland ridge in both model and observations. While the large-scale spatial anomaly patterns are similar, VIKING20X anomalies show more small-scale structure than EN4.

Time series of mean temperature and salinity averaged over the upper 500 m in three areas along OSNAP_E (Rockall Trough, Iceland Basin, and Irminger Sea, Fig. 2) from VIKING20X and EN4 show that the fastest freshening and cooling at OSNAP latitudes occurred in the Iceland Basin through 2014 and 2015, with the lowest temperatures and salinities recorded in 2016. During the relatively short OSNAP observational period, 2014 to 2018, the model temperature and salinity signals closely track the OSNAP observational data on annual timescales, though they do not not match all the higher-frequency variability. Compared to EN4, in the upper 500 m, VIKING20X has a warm, salty bias in the Iceland Basin between 2004 and 2014, and a cool, fresh bias in the Irminger Sea from 2008 to the present day.

Vertical sections of model temperature and salinity anomalies (Fig. 1e, f) show the cooling and freshening to be confined to the surface 1000 m east of the central Irminger Sea, the region occupied by generally northward-flowing waters of the upper limb of the AMOC. Deeper waters show warming and salinification. This compares well with changes described in OSNAP section observational data (Fig. 6 of Holliday et al., 2018). While the OSNAP time series is short, it spans the period of fastest upper ocean freshening and cooling.

The evidence presented here shows that VIKING20X reproduces much of the characteristic spatial structure and timing of the exceptional freshening event. As a hindcast, VIKING20X is driven by reanalysis surface forcing but runs freely without data assimilation or relaxation to observations, so the similarity to observations is a product of the external forcing and modelled ocean dynamics. While we do not attempt to quantify the model's skill in reproducing the observations, the qualitative comparisons support the idea that we can use the VIKING20X model to help diagnose the causes of the exceptional freshening event. We discuss the possible implications of inaccuracies in the model representation of the freshening event on our conclusions, and caveat those conclusions appropriately, in Sect. 7.

There are fundamentally only three possible contributions to the exceptional freshening and cooling observed in the eastern subpolar North Atlantic. Firstly, changing proportions of transport coming from the various source regions (including the introduction of new source regions). Each different source region has different average temperature and salinity (TS) properties, so combining them in changing ratios will change TS properties in the eastern subpolar North Atlantic. Secondly, changing TS properties of transports out of the individual source regions. This could alter TS properties in the eastern subpolar North Atlantic even without any change in relative volume transports. And, finally, local processes – surface fluxes and internal mixing – including processes along-track between the source and the eastern subpolar region, causing changing fluxes of freshwater and heat between the upper eastern subpolar North Atlantic waters and the surrounding waters and the atmosphere.

4 Lagrangian particle tracking results

4.1 Changing volume transports by source and pathway

Annual mean total upper ocean volume transports in VIKING20X across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (Fig. 5) are high (20 Sv) prior to 1996, reduce by 10 %–20 % between 2000 and 2009, and return to 20 Sv after 2016. Isolating these transports by source and pathway shows gradually decreasing volumes of Gulf Stream source waters since 1996. For Labrador Sea source waters, the picture is mixed, with the loop path transport declining similarly to the Gulf Stream source transports, while the direct path transports increase rapidly between 2008 and 2016. As a consequence, prior to 1996, water crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ had 2–3 times as much Gulf Stream source water as Labrador Sea source water; by 2016, the ratio was 0.8–1.2:1 Labrador Sea source to Gulf Stream source. Notice the differences here with the modelling work of Koul et al. (2020), who – also using particle tracking, but in a coarser resolution model – find a maximum of 11 % of near-surface water in the eastern subpolar North Atlantic to be of subpolar (Labrador Sea) origin. In Asbjørnsen et al. (2021), 26 % of the Iceland–Scotland Ridge inflow has a subpolar or Arctic origin (via the Labrador Sea). Here we find never less than 22 % and up to 45 % of northward transport to be of subpolar origin. These proportions vary depending on the depths and longitude ranges considered. In the 1990s, both loop and direct paths from the Labrador Sea paths are equally important, with contributions of ∼2.5 Sv each; in the 2010s, the Labrador Sea contribution mainly consists of waters that followed the direct path (∼7.5 Sv), while the loop path contribution becomes small (∼1.0 Sv). The transport volumes are shown (Fig. 5) for the tracks crossing OSNAP_E in the upper 500 m, the maximum and minimum transport proportions are the same for tracks crossing in the upper 1000 m or in the upper limb (σ₀<27.62 kg m⁻³) but the total transports increase to 30 and 26 Sv respectively.

Figure 5 Time series of 1-year running mean volume transport across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ by source and pathway plotted against (a) the time at which water parcels cross OSNAP_E and (b) the time that water parcels leave the source region. The total northward transport (black line) is higher in the 1990s and 2010s and lower in the 2000s. Transport from the Gulf Stream (green) reduces steadily after 1996. From the Labrador Sea, loop path transport (orange) also reduces steadily after 1996, while direct path (blue) increases rapidly after about 2008. Transport with an unidentified source after 10 years of backward tracking is coloured grey. Error bars at ±1.96 standard deviations around the mean are included but are mostly too narrow to be visible here. The dashed blue line (and right-hand y axis) in (a) and (b) is the VIKING20X subpolar gyre index (SPGI) from Biastoch et al. (2021).

Due to the range of particle travel times (Fig. 4a), plotting volume transports against the time at which particles cross OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ temporally smooths changes happening at the sources. This introduces the risk of interpreting rapid changes at the sources as gradual changes at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$, so we also plot volume transports as particles leave the source region. Transports leaving the major sources (Fig. 5b) show similar patterns of variability compared to those recorded at OSNAP, but with source volume transport changes leading OSNAP changes by 3–5 years (corresponding to the most frequent transit times between the source sections and OSNAP, shown in Fig. 4a). In particular, transport changes at the sources still appear gradual and are maintained over several years, rather than being sudden step changes. There are signs of an increase in the Gulf Stream and loop path transports leaving the source in 2016, but these increases are yet to be reflected in corresponding increases across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ by the end of 2019.

The general picture over the 30 years of modelled transports across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ is therefore one of gradually reducing volume transports of Gulf Stream origin combined with a more rapid increase in total volume transports of Labrador Sea origin after about 2008. This later increase in volume transports of Labrador Sea origin occurs entirely via the direct path and is new transport rather than a switch in pathways from loop to direct; loop path transports continue to slowly decline in parallel with Gulf Stream origin transports. In addition, previously reported VIKING20X model results (Biastoch et al., 2021) show that the denser deep western boundary current (DWBC) transport out of the Labrador Sea decreases after 1996. Together, these point to a shift in export volume through the major export routes from the Labrador Sea.

4.2 Temperature and salinity variability, source regions, and along-track fluxes

Averaged TS properties at the source for all particles crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ at time t are calculated from the temperature and salinity of each particle as it leaves its source region. These are volume averages of the water as it leaves the source region. Neither Labrador Sea origin water nor Gulf Stream origin water show freshening or cooling of comparable magnitude to the total freshening and cooling at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (Fig. 6). In fact, at their Labrador Sea source, the waters which cross OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ between 1998 and 2012 appear to have become warmer and saltier, with only slight at-source freshening and cooling of waters which cross OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ after 2012. Gulf Stream origin water crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ after 2012 showed a small amount of freshening and cooling as it left the source, but as for the Labrador Sea origin, this signal is small compared to the total overall freshening and cooling at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. There is more variability (and more uncertainty) in the at-source (i.e. 10 years before crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$) properties of the “other” particles, (Fig. 6, grey lines), but these form less than 10 % of the overall transport and only show marked at-source cooling and freshening in particles crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ after 2016.

From the transports from the different sources and the mean TS properties at each source, we can calculate the “expected” mean properties of water transported across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ in the absence of along-track external fluxes (surface fluxes or mixing with surrounding waters) (dashed black lines in Fig. 6). The difference in properties between these “expected” values and the recorded properties (solid black lines in Fig. 6) gives the contribution from along-track external fluxes. For salinity, these along-track external fluxes act to slightly increase the overall mean salinity between source and OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. There is a small consistent trend in the magnitude of this contribution, resulting in weaker along-track salinification after about 1998, hence adding to the net freshening during this period. For temperature, the along-track external fluxes cool the water, as expected with the overall heat loss to the atmosphere in the subpolar North Atlantic. In contrast to salinity, the mean along-track external cooling reduces in time from 2.2 ^∘C in 1998 to 1.0 ^∘C in 2017. Note that this reducing trend in along-track externally driven cooling opposes the observed cooling at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$.

4.3 Quantifying the relative contributions to eastern subpolar North Atlantic freshening and cooling

Noting that, by construction, each individual particle transports the same volume, we can write down expressions for the mean temperature and salinity of the volume transported northwards across OSNAP_E at time t in terms of the properties and transports from each source region and the in-transit changes in mean temperature and salinity:

$$\begin{array}{}\text{(7a)}& {\displaystyle }\underset{\mathrm{at}{\mathrm{OSNAP}}_{\mathrm{E}}}{\underbrace{\stackrel{\mathrm{‾}}{T_{\text{o}}}}}& {\displaystyle }=\underset{\mathrm{at}\mathrm{source}}{\underbrace{\stackrel{\mathrm{‾}}{T_{\text{g}}}{\displaystyle \frac{V_{\text{g}}}{V_{\text{all}}}}+\stackrel{\mathrm{‾}}{T_{\text{l}}}{\displaystyle \frac{V_{\text{l}}}{V_{\text{all}}}}+\stackrel{\mathrm{‾}}{T_{\text{u}}}{\displaystyle \frac{V_{\text{u}}}{V_{\text{all}}}}}}+\underset{\mathrm{transit}}{\underbrace{\stackrel{\mathrm{‾}}{\mathrm{\Delta }{T}_t}}},\text{(7b)}& {\displaystyle }\underset{\mathrm{at}{\mathrm{OSNAP}}_{\mathrm{E}}}{\underbrace{\stackrel{\mathrm{‾}}{S_{\text{o}}}}}& {\displaystyle }=\underset{\mathrm{at}\mathrm{source}}{\underbrace{\stackrel{\mathrm{‾}}{S_{\text{g}}}{\displaystyle \frac{V_{\text{g}}}{V_{\text{all}}}}+\stackrel{\mathrm{‾}}{S_{\text{l}}}{\displaystyle \frac{V_{\text{l}}}{V_{\text{all}}}}+\stackrel{\mathrm{‾}}{S_{\text{u}}}{\displaystyle \frac{V_{\text{u}}}{V_{\text{all}}}}}}+\underset{\mathrm{transit}}{\underbrace{\stackrel{\mathrm{‾}}{\mathrm{\Delta }{S}_t}}},\end{array}$$

where all variables are functions of the particle release time, t. Here, $\stackrel{\mathrm{‾}}{T_{\text{o}}},\stackrel{\mathrm{‾}}{S_{\text{o}}}$ represent mean T, S properties at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$; V_all is the total volume transport across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ at time t; V_g, V_l, and V_u represent volume transports from the Gulf Stream, Labrador Sea, and unidentified “other” source respectively (and $V_{\text{g}}+V_{\text{l}}+V_{\text{u}}=V_{\text{all}}$); $\stackrel{\mathrm{‾}}{T_{\text{g}}},\stackrel{\mathrm{‾}}{S_{\text{g}}},\stackrel{\mathrm{‾}}{T_{\text{l}}}$, etc. represent mean properties of particles leaving the relevant source; and $\stackrel{\mathrm{‾}}{\mathrm{\Delta }{T}_t},\stackrel{\mathrm{‾}}{\mathrm{\Delta }{S}_t}$ denote in-transit changes to mean temperature and salinity.

The dependence of V_all on V_g, V_l, and V_u means we cannot completely separate the contributions to changing temperature and salinity at OSNAP_E from changes to individual properties at, and transports from, the different source regions. But we can estimate each of these contributions by varying the six independent variables one at a time while holding the others at their mean values. This will tend to underestimate (overestimate) the influence of changes at one source when transports from the other sources are below (above) their mean values.

Figure 6 Salinity (a) and temperature (b) time series at the source and at OSNAP_E. For salinity, (a), time series of the 1-year running mean salinity at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (solid black lines) and at source (dashed black lines) show the mean freshening at the source to be reflected in an overall freshening recorded at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. This contrasts with temperatures (b), where an overall at-source cooling (dashed black line) is mitigated by reduced along-track cooling (difference between the solid and dashed black lines), resulting in reduced cooling at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (black solid line). Dividing the properties at source into the three origins (Gulf Stream (green), Labrador Sea (orange), and other (grey)), none shows consistent freshening or cooling. The x axis is the time that water parcels cross OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. Shaded areas show ±1.96 standard deviations between the 32 subsets. Uncertainties are larger for the “other” source because the sample size is smaller and these particles originate in a wide range of conditions throughout the subpolar North Atlantic south of OSNAP_E rather than in tightly defined source regions.

Table 1 Contributions to the freshening and cooling of water crossing OSNAP_E in 2016–2017 (peak freshening and cooling) relative to 2011–2012 (the start of the more rapid freshening and cooling). Contributions are listed by the components described in Sect. 4.3 and Fig. 7 for water crossing OSNAP_E between 37^∘ W and 0^∘ in (1) the top 500 m, (2) the top 1000 m, and (3) the upper limb (σ₀<27.62 kg m⁻³). The largest contribution to both freshening and cooling is made by the changing volume of water of Labrador Sea origin. The single component contributions do not necessarily sum to the total, as the component contributions are not independent (see Eq. 7).

Following this process (Fig. 7), qualitatively we find that the dominant driver of variability for both salinity and temperature is the volume transported from the Labrador Sea to OSNAP_E, V_l. The variability in the volume transported from the Gulf Stream, V_g, plays an additional but smaller role. Changes in source temperature and salinity show weaker relationships to either overall mean properties at origin or properties at OSNAP_E. Only for temperature, we find a notable influence of in-transit, externally driven changes, $\stackrel{\mathrm{‾}}{\mathrm{\Delta }{T}_t}$.

Figure 7 Relative contributions of specific causes to the salinity (a, c, e, g) and temperature (b, d, f, h) changes observed at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. These are plots of 1-year running mean salinity and temperature anomalies. We see that, for salinity (a), the overall freshening at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (black line) closely follows the freshening of the combined sources (blue line), with a small additional contribution from changes during transit (pink). Temperature (b) shows in-transit changes (pink) opposing source changes (blue). In the panels below those, the contributions of six individual processes (orange and grey lines) to the overall source characteristics (blue lines) are estimated by holding the other five processes at their mean values (see Eq. 7). These six processes are: volume transport variability (orange lines, $V_{\text{l}},V_{\text{g}},V_{\text{u}}$) and salinity/temperature variability (grey lines, $\stackrel{\mathrm{‾}}{S_{\text{l}}},\stackrel{\mathrm{‾}}{T_{\text{l}}},\stackrel{\mathrm{‾}}{S_{\text{g}}},\stackrel{\mathrm{‾}}{T_{\text{g}}},\stackrel{\mathrm{‾}}{S_{\text{u}}},\stackrel{\mathrm{‾}}{T_{\text{u}}}$) for the Labrador Sea source (c, d), Gulf Stream source (e, f), and other source waters (g, h).

We now quantify the contributions of each of the variables on the RHS (right-hand side) of Eq. (7) to the changing properties of water crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ during the period of most rapid freshening and cooling, between 2011–2012 and 2016–2017 (Table 1). For salinity, over 60 % of the freshening is due to increasing volumes of water from the Labrador Sea, and 27 % due to decreasing volumes from the Gulf Stream. The small remainder is a combination of along-track external fluxes and changes in TS properties in the source regions, particularly slight freshening of the Gulf Stream source. For temperature, the picture is more complicated: half of the overall mean cooling of water leaving the combined source regions is mitigated along-track by reduced heat loss from the cooler water. The net all-source mean cooling of water leaving the source regions is, as for salinity, driven mostly by increasing volumes of water from the Labrador Sea (over 55 %) and decreasing volumes of water from the Gulf Stream (20 %). As for freshening, a small contribution to the cooling at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ stems from temperature changes at the sources; specifically, from the reducing mean temperature of waters leaving the Gulf Stream source. Table 1 also includes results for water crossing OSNAP_E in the upper 1000 m and in the AMOC upper limb (σ₀<27.62 kg m⁻³). These results are qualitatively the same as for OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ but with a slightly lower – but always the largest – proportion of the freshening and cooling being due to changing volume transport from the Labrador Sea source region, and with higher proportions due to changing volume transport and source TS from the Gulf Stream source region. These different proportions are due to water of Gulf Stream origin on average occupying deeper levels further east; so, as the averaging region goes deeper (to 1000 m) and focuses further east (AMOC upper limb), the proportion of Gulf Stream origin water in the mix increases, and with it the relative importance of changes in transports and TS properties of Gulf Stream origin water.

It is important to note that in the later years modelled, when volume transports of Gulf Stream origin waters are below their 30-year mean and those of Labrador Sea origin are above their mean, our analysis underestimates the influence of Labrador current origin transports on OSNAP_E properties and overestimates the influence of Gulf Stream origin transports. However, this only further stresses the importance of the described changes in the Labrador Sea for the development of the recent eastern cold and fresh anomalies.

4.4 Summary of particle tracking results

We find only small changes in the modelled TS properties of the Gulf Stream and the Labrador Sea origin waters as they leave their respective source regions and also no substantial freshening or cooling through along-track changes. Hence, the exceptional freshening and cooling in the eastern subpolar North Atlantic is driven by the increasing proportion of fresher Labrador Sea origin waters in the mix compared to saltier Gulf Stream origin waters. We have not identified any rapid – annual or intra-annual – changes in source properties, transports, or pathways driving the freshening and cooling, but rather pentadal to decadal, persistent changes beginning in the early 2000s. The increasing volume of Labrador Sea origin water accounts for over 60 % of the freshening, in combination with the reduced transport from the Gulf Stream source, accounting for about 30 %. For temperature, the shift in relative proportions of water from the two main source regions drives the cooling, but the average cooling of the combined-source waters is mitigated (reduced by about half) by reduced heat loss along-track.

In the following sections, we explore and discuss the possible mechanisms responsible for the modelled changes in the quantities and proportions of waters transported northwards through the eastern subpolar North Atlantic from the two major sources identified.

5 Mechanisms

5.1 North Atlantic Current variability and SPG eastward expansion

The distribution of transit times for particles to travel from their source region to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (Fig. 4a) show shortest times of under a year and a distribution mode of about 2 years. The distributions are skewed right, with a long tail out to longer transit times. The transit time distributions have a similar shape for both identified source regions and all years, with a steep early rise and long tail at longer transit times, so we can characterize transit time variability by examining the variability in the mean. Splitting the section at 21^∘ W (Fig. 8), we find contrasting changes over time of mean transit times, speeds, and track lengths for particles crossing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ in the east (Rockall–Hatton Bank and Rockall Trough) and west (Iceland Basin and Irminger Sea).

Figure 8 Evolution over time of particle transit times. Time series of mean volume transports (a, b), transit times (c, d), speeds (e, f), track lengths (g, h), and OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ crossing longitudes (i, j) of particles by source region. Results are divided into OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ crossings west (left column) and east (right column) of 21^∘ W, the western flank of the Rockall–Hatton Bank. Times on the x axis are crossing times of OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$, so values represent weighted averages of conditions along-track in the preceding years. Time series are smoothed with a 1-year rolling mean.

In the east, over Rockall–Hatton Bank and in the Rockall Trough, northward transports are dominated throughout by waters of Gulf Stream origin. Some increase in transports is seen after 2004 (Fig. 8b) as Gulf Stream source waters are increasingly confined to this eastern section in the later half of the period studied. Transit times, speeds, and track lengths to this eastern section of OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ from the individual sources show little consistent change over time (Fig. 8d, f, h) except for a slight slowing of the loop path (though volumes on this path are very small).

In the Iceland Basin and Irminger Sea, across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ west of 21^∘ W, the Gulf Stream influence reduces steadily between 2004 and 2016 and is largely replaced by direct path waters from the Labrador Sea (Fig. 8a). Mean transit times in the west, from all sources, increase by around 30 % between 1998 and 2016 (Fig. 8c). This is due to a combination of reduced current speeds and increasing track lengths (Fig. 8e, g). The increasing track lengths result from an eastward shift of the mean OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ crossing longitudes (Fig. 8i) and increased along-track eddying (results supported by Eulerian analysis of VIKING20X NAC, not shown). After 2016, transit times from all sources to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ west of 21^∘ W start to decrease, with faster speeds and shorter paths recorded, perhaps associated with a strengthening subpolar gyre from 2013.

These results are consistent with increasing transport in the upper layers of the SPG and eastward movement of the subpolar front between 2006 and 2016, with the reduced volumes of Gulf Stream origin surface waters increasingly confined to the east, passing OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ over Rockall–Hatton Bank and Rockall Trough. The Irminger Sea and Iceland Basin sections of OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ are correspondingly increasingly dominated by waters from the Labrador Sea. The OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ total is dominated by the larger transports west of 21^∘ W. Note that our methodology, calculating the mean salinity and temperature of volumes transported northwards across the whole upper eastern SPNA (rather than area-average salinities and temperatures), precludes pure subpolar gyre eastward expansion (with no new transport) as a cause of the overall cooling and freshening of transport across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. However, dividing the section zonally suggests that eastward movement of the subpolar front could contribute to redistribution of heat and salt transport along the OSNAP_E), causing enhanced local freshening and cooling in the Iceland Basin and reduced freshening and cooling over Rockall–Hatton Bank and in the Rockall Trough. The fractional evolution seen for the eastern section is similar to that seen for the Nordic Seas inflow (Asbjørnsen et al., 2021), consistent with the dominant pathways to the Nordic Seas being on the western flank of Rockall–Hatton Bank and in the Rockall Trough.

Changing current speeds and pathways in the SPG and NAC have been linked to eastern subpolar North Atlantic freshening (Bersch, 2002; Sarafanov, 2009; Bersch et al., 2007; Hátún et al., 2005; Koul et al., 2020). The proposed mechanism is that eastward expansion and strengthening of the SPG regulates the northward transport of warm salty water from the subtropics, modulating the proportion of subpolar and subtropical waters reaching the eastern North Atlantic (Häkkinen et al., 2011). In turn, the strengthening and eastward SPG expansion cause increased freshwater flow into the Nordic Seas (Kenigson and Timmermans, 2020; Asbjørnsen et al., 2021). Using a subpolar gyre index (SPGI) based on the second principal component of sea surface height variability (Koul et al., 2020), Biastoch et al. (2021) show a weak SPG in VIKING20X from 1996–2012, with increasing strength from 2013 (we reproduce this index in Fig. 5). Here we find that, while reduced salinity in the upper eastern subpolar North Atlantic is due to an increasing proportion of subpolar waters, this freshening, starting as early as 2008, appears to lead the increase in SPG strength. We can compare the late 2010s with the previous strong SPG of the early 1990s, during which VIKING20X shows an increased SPGI but a much smaller increase in the proportion of subpolar water at OSNAP_E (Fig. 5) and a correspondingly smaller cooling and freshening signal in the combined source TS (Fig. 7a, b).

SPG decadal-scale dynamics are complex and still poorly understood. Häkkinen et al. (2011) and Chafik et al. (2019) observe that a strong SPG during periods of eastern subpolar North Atlantic cooling is associated with a strengthened wind-stress curl, a horizontally expanded gyre, a southeastwards-shifted NAC pathway, and reduced advection of warm and saline subtropical waters into the northeast Atlantic. Kenigson and Timmermans (2020) show how the 2014–2016 period was characterized by strong positive NAO forcing, strengthened SPG, and eastward expansion of the SPG, resulting in freshening in the northeastern North Atlantic Ocean. Alternatively, Robson et al. (2016) propose that cooling of the upper eastern subpolar North Atlantic is due to reduced transport of warm water from the south associated with decreased AMOC, which in turn is linked to record low densities (lighter water) in the deep Labrador Sea.

The VIKING20X model results which we present here suggest a decoupling of the eastward expansion of the SPG from the strength (as measured by an SSH-based SPGI), with eastward expansion beginning during a period of weaker SPG. The cooling we observe in transports at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ beginning before 2012 is mainly linked to an increasing proportion of subpolar water (Fig. 7) comprising increasing volumes of water which left the Labrador Sea from around 2008 (Fig. 5). While we cannot discount that SPG variability contributes to the exceptional cooling, the trends in upper-layer volume transports, which drive the freshening and cooling in the eastern subpolar North Atlantic, begin before 2008 in a period of weak, and possibly still weakening (Biastoch et al., 2021), SPG.

5.2 Subtropical-origin water and meridional overturning

There is currently no consensus on whether the AMOC has declined since the 1990s (see, for example,  Jackson et al., 2022; Latif et al., 2022, for summaries of the debate). Models and proxies suggest a decline (e.g. Rahmstorf et al., 2015; Caesar et al., 2018, 2021), while observation-based analyses give a mixed picture, with some finding no significant AMOC decline (Fu et al., 2020; Worthington et al., 2021; Caínzos et al., 2022), the RAPID Meridional Overturning Circulation and Heat-flux Array (RAPID-MOCHA) (Cunningham et al., 2007) finding reduced AMOC after 2009 (Bryden et al., 2020), and the application of the Bernoulli inverse to hydrography (Fraser and Cunningham, 2021) showing weakening AMOC during our experiment period.

Here, we show results at 29^∘ N (Fig. 9), a latitude close to the RAPID observation line and to our Gulf Stream source definition latitude. We choose 29^∘ N, near the meridional centre of the subtropical gyre, because of the simple configuration of the large-scale flow at this latitude. The western boundary current (Gulf Stream) containing almost all the northward warm-water transport sits entirely on the 1000 m deep shelf; the balancing return southward flow is made up of the deep western boundary current – the lower limb of AMOC – along the slope at the shelf edge and the broad southward flow in the upper 1000 m of the interior and east of the subtropical gyre, which we designate subtropical gyre recirculation. At 29^∘ N in VIKING20X-JRA-short, the Gulf Stream flow weakens between 1990 and 2020. This is accompanied by reduced AMOC and therefore reduced lower limb return flow. The remaining southward return flow, the subtropical gyre recirculation, shows no overall trend. So the Gulf Stream weakening in VIKING20X is associated with reduced flow northwards from the subtropical towards the subpolar North Atlantic in the upper layers and corresponding reduced deep return flow – that is, reduced AMOC. This weakening of AMOC, by about 30 %, is associated with the reduction in modelled transport of Gulf Stream origin to OSNAP_E between the 1990s and 2016 (Fig. 5) diagnosed by the Lagrangian tracking.

Figure 9 Northward transports and AMOC strength at 29^∘ N in VIKING20X. Pale lines are monthly mean data; darker lines are low-pass filtered using a Butterworth filter and 18-month cutoff. The total transport in the Gulf Stream (black line) reduces throughout the experiment duration, while the volume recirculating in the subtropical gyre above 1000 m (light orange line) shows a smaller trend or remains largely constant. This implies weakened AMOC and reduced northward flow in the AMOC upper limb (blue, AMOC in z-space). At 29^∘ N, the AMOC in σ-space (green) is very similar to AMOC_z and also reduces consistently. The dark orange line is AMOC in z-space from RAPID observations, low-pass filtered at 18 months, for comparison. The reducing AMOC strength (blue and green) closely reflects the reducing contribution of Gulf Stream origin transport at OSNAP_E found from the particle tracking.

Previous longer-term, lower-resolution modelling work (Koul et al., 2020) linked periods of higher subpolar gyre index to weaker flow from the subtropical gyre to the NAC, suggesting a possible “blocking” mechanism. But the AMOC weakening in VIKING20X between 1990 and 2020 appears to be independent of the subpolar gyre strength. In this period, the AMOC weakening contrasts with a high (1990s) to low (2000s) and back to high (late 2010s) cycle in the subpolar gyre index.

Examination of the decadal-scale drivers of the AMOC in either the real ocean or in VIKING20X is beyond the scope of the current work, but we refer the reader to Biastoch et al. (2021) for a detailed analysis. They find outflow of deep water from the Labrador Sea at 53^∘ N to be a good indicator of the subpolar AMOC trend and also a strong dependence on the choice of freshwater forcing in the model. The VIKING20X-JRA-short run used here shows a stronger decline in the AMOC than either observations or other VIKING20X realizations since 1990, but all consistently suggest a decline for this period.

Comparing the VIKING20X AMOC time series since 2005 with the observational time series from RAPID (see e.g.  Bryden et al., 2020), in VIKING20X the AMOC appears 3–4 Sv too weak, but the model reproduces the observed overturning minimum in 2009 and the reduction in overturning by 2–3 Sv from before 2009 to after 2009. The 2009 minimum is suggested by Bryden et al. (2020) as the start of a period of reduced northward transport of heat and salt from the subtropics to the subpolar gyre. The similarities in AMOC at the RAPID array between model and observations since 2005 give us confidence that we could be, at least qualitatively, reproducing the contribution of changes in overturning circulation to reduced eastern North Atlantic salinity and temperature. In the first half of our time series, before the start of the RAPID observations, VIKING20X shows declining AMOC, with consistent weakening at latitudes throughout the North Atlantic (Biastoch et al., 2021). Compared to observation-based estimates, this AMOC decline between 1990 and 2005 in VIKING20X may be an overestimate, which could in turn lead to us overestimating the role of the AMOC weakening in the observed eastern subpolar gyre freshening and cooling.

In summary, in Sect. 4.3 we found that reduced transport from the Gulf Stream source accounted for about 30 % of the cooling and freshening observed in the model at OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. Examining this in more detail shows a longer-term gradual decline in subtropical water transported across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ between the mid-1990s and 2016, contrasting with the accelerating freshening and cooling from 2011–2016. This decline is associated with AMOC weakening in VIKING20X. The reduced transports of Gulf Stream origin will amplify TS variability in the eastern subpolar North Atlantic due to variability in subpolar gyre transport from the Labrador Sea origin (Eq. 7).

5.3 Subpolar and Arctic origin water

5.3.1 Pathways from the Labrador Sea

Upper water from the Labrador Sea reaches the eastern subpolar North Atlantic along one of two dominant pathways (Holliday et al., 2020; New et al., 2021): a direct pathway and a loop pathway which follows the coast southwestwards through the Slope Sea north of the Gulf Stream before reversing, joining the Gulf Stream origin waters travelling to the northeast (Fig. 3d, e). A change in the dominance of the direct pathway for cold fresh Arctic-origin water from the Labrador Current has been offered as a cause of the exceptional freshening event (Holliday et al., 2020).

Examining the changing volume transports along these two pathways between the Labrador Sea (Fig. 5) and OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ shows that the ratio of the direct to the loop source volume changes from close to 1:1 before 2000 to 8:1 after 2016. Prior to 2006, dominance switches repeatedly between the two pathways, with little change in total transport (Fig. 5b). This early period supports the ideas of changing pathway dominance of Holliday et al. (2020), but after 2006, the net increase in volumes leaving the Labrador Sea towards OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ in VIKING20X is entirely along the direct path.

It is this more recent increase in the total volume of Labrador Sea source water, rather than details of the pathways taken, which we showed in Sect. 4 is responsible for the major part of 2012–2016 freshening and cooling in the eastern subpolar North Atlantic. In the following sections we explore possible mechanisms driving this increased total transport from the Labrador Sea.

5.3.2 Source regions upstream of the Labrador Sea outflow

Increased transport from the Labrador Sea must be associated with changes in transport from source regions upstream of the Labrador Sea outflow. Tracking the water of Labrador Sea origin back in time further upstream (Fig. 10), we can distinguish between four upstream source regions: Hudson Bay; the Davis Strait; from the Greenland Sea, southward across the Greenland–Scotland ridge; and SPG water (particles circulating in the SPG that previously passed northwards across OSNAP_E).

The overall increase in volume of Labrador Sea origin water (4 Sv) primarily originates from an increase in circulating SPG water (3 Sv), with a smaller contribution from increases in waters of Greenland Sea and Davis Strait origin (1 Sv) (Fig. 10). While volumes taking the direct path increase and those taking the loop path decrease, we find little variability, and no trend, in the percentages of waters from the four upstream source regions forming either path. The direct path is about 80 % SPG origin, 10 % Davis Strait origin, and 10 % Greenland Sea origin water; the loop path is 40 % SPG origin, 40 % Davis Strait origin, 10 % Greenland Sea origin, and 10 % Hudson Bay origin water. All the flows we are examining here are in the light upper-layer waters which cross OSNAP_E in the surface 500 m.

Figure 10 Source regions upstream of the Labrador Sea outflow. Quantifying contributions from sources upstream of the Labrador Sea outflow for the direct path (a) and the loop path (b). Upper-layer inflow to the Labrador Sea comes from Hudson Bay (pink), Davis Strait (brown), and the West Greenland Current. The West Greenland Current is further divided here into water from the Greenland Sea (yellow) and recirculating SPG origin water (light blue). There is little change in the proportions of water from different upstream sources forming the different paths (loop or direct). Overall, the increased contribution of the direct path (a) is reflected in an increase in the total volume of recirculating SPG transport and a small increase in the total contribution from the Greenland Sea. Particles from Hudson Bay, Davis Strait, and the Greenland Sea have come directly from those regions; SPG origin waters are circulating in the subpolar gyre and have previously crossed northwards across OSNAP_E.

There is little change in total volume contribution from Hudson Bay or from Davis Strait upstream source regions, although the dominant pathway taken by the Davis Strait origin water switches from loop to direct. This is exactly the switch of dominance described by Holliday et al. (2020).

The small increase in the contribution of upper waters from the Greenland Sea is interesting and perhaps suggests increased Greenland meltwater flowing southward in the East Greenland Current. The fate of this fresh water in regions of deep water formation has led to speculations that the accelerating melting of the Greenland Ice Sheet could stratify the subpolar gyre, change the spatial deep convection patterns, and slow or stop the AMOC (see, for example, Böning et al., 2016; Foukal et al., 2020; Rühs et al., 2021). This increase accounts for <20 % of the increase in transport from the Labrador Sea to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$.

The largest part ∼75 %, 3 Sv, of the increased transport from the Labrador Sea is due to increased volumes recirculating within the upper North Atlantic subpolar gyre. This recirculation doubles, from 3 to 6 Sv, between the early 2000s and 2016, representing an increase from ∼15 % to ∼30 % of the northward transport across OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. It is this recirculation which is driving the increase in transport from the Labrador Sea to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$, and in turn the freshening and cooling.

There are various possible mechanisms behind this increased contribution from water recirculating within the near-surface SPG. One is an accelerating SPG, but in Sect. 5.1 we showed that the timing of the recirculation increase (from 2008) and the major SPG acceleration (after 2013) suggest this is an unlikely cause. To leading order, water transported northwards in the upper layers of the eastern subpolar North Atlantic generally does one of three things: it crosses the GSR northwards, is made denser in the Nordic and Arctic seas, and returns in deeper overflows; it is made denser within the SPG and sinks to intermediate and deep layers; or it continues to circulate within the SPG. To increase the proportion recirculating, one of the other pathways needs to weaken, suggesting reduced densification processes in either the subpolar gyre or Nordic and Arctic seas. We return to these ideas in Sect. 6).

5.3.3 Labrador Sea outflow

A VIKING20X mean velocity section across the shelf and shelf break at the western end of the OSNAP_W line spans the upper layer outflow, the Labrador Current, the source region of all the modelled transport from the Labrador Sea to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ (Fig. 11a). This shows the upper layer outflow to be concentrated in two main cores: one close to the coast – the Labrador Coastal Current containing Hudson outflow and water from the Davis Strait – and the other (larger and stronger) over the shelf break and slope. Most of the water exiting the Labrador Sea with a density (σ₀) lighter than 26.65–26.70 kg m⁻³ is transported round to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ in the surface 500 m. Pentadal mean velocity anomalies (Fig. 11b–g) across this Labrador Sea outflow section show no coherent overall increase in outflow velocities (negative, blue, in Fig. 11b–g).

The offshore current has previously been separated into distinct shallower inshore (upper 200 m inshore of the 600 m isobath) “LC-Arctic” and deeper offshore “LC-Atlantic” components (Florindo-López et al., 2020; Holliday et al., 2020; New et al., 2021). While the definitions of these outflow components do not map precisely onto our upstream origins (Sect. 5.3.2), the LC-Arctic component is water of Arctic origin and is primarily made up of waters we describe as being of Davis Strait and Greenland Sea origin, and the LC-Atlantic outflow corresponds closely to our recirculating SPG origin water. The LC-Arctic and LC-Atlantic components are not easily separable in the VIKING20X mean outflows at OSNAP_W (Fig. 11a), though some variability in the LC-Arctic component is visible near-surface at about 360 km offshore in the pentadal velocity anomalies (Fig. 11b–g). The VIKING20X anomalies in the LC-Arctic outflow are consistent with Florindo-López et al. (2020), showing higher outflow in the early 1990s and late 2010s and reduced outflow between. Note that the total LC-Arctic outflow is typically less than 2 Sv (New et al., 2021), so the increased LC-Arctic outflow after 2005 forms a minor component of the total outflow increase seen in VIKING20X, which is dominated by the offshore LC-Atlantic, SPG water.

Pentadal mean density fields show deepening of all isopycnals between 27.60 and 27.75 kg m⁻³ between the 1990–1994 and 2010–2014 pentads (Fig. 11b–g), followed by a slight recovery to shallower depths by 2015–2019. It is the increase in cross-sectional area occupied by these outflowing light waters (27.60 kg m${}^{-\mathrm{3}}<{\mathit{\sigma }}_{\mathrm{0}}<\mathrm{27.75}$ kg m⁻³), rather than increased flow speeds, which produces the increased transports of light upper waters out of the Labrador Sea.

Figure 11 Vertical sections of Labrador Sea upper-layer outflow 1990–2019 mean velocity (a) and pentadal velocity anomalies from the mean (b)–(g), with period-mean density contours superimposed, at the western end of the OSNAP line. Velocities, and velocity anomalies, are shown as positive into the Labrador Sea. Lighter waters (<27.65 kg m⁻³) leaving the Labrador Sea are mostly transported to OSNAP_E in the subpolar gyre; this transport is dominated by strong currents at the shelf break. Outward currents in the lighter density classes are weaker (positive – red – anomalies) than the 30-year average between 2000 and 2015 (d–f). The increasing contribution of Labrador Sea origin waters to eastern SPG transports in the 2010s is associated with the deepening of these lighter layers rather than increasing currents; notice the highlighted 27.65 kg m⁻³ isopycnal in (b)–(g).

Note that there is no obvious increase or decrease in the cross-sectional area occupied by waters with σ₀<27.30 kg m⁻³. These surface waters on the shelf form the Labrador Coastal Current, which, combined with some inshore near-surface water from the shelf break, mainly follows the loop path. Labrador Current outflow along the loop path has been found to be related to Gulf Stream transport and AMOC strength (New et al., 2021; Sanchez-Franks et al., 2016), possibly via wind stress in the subpolar gyre.

Qualitatively, we see the potential for increasing volume transport from the Labrador Sea to the upper eastern subpolar North Atlantic to be due to a deepening of the surface layer of lighter water in the Labrador Sea rather than faster currents. We now explore that idea further, examining the increasing outflow of lighter waters using a water mass transformation budget in the Labrador Sea (Walin, 1982; Speer and Tziperman, 1992).

6 Labrador Sea water mass analysis

6.1 Full water mass balance

To first order, the Labrador Sea is fed by water from the West Greenland Current, through the Davis Strait, from Hudson Bay together with river, meltwater, and precipitation inflows; this water is transformed to denser water within the Labrador Sea by winter cooling and ice formation before exiting southwards via the Labrador Current in the upper ocean and the DWBC below.

The mean water mass transformations (Walin, 1982; Speer and Tziperman, 1992) on the ocean model volume bounded by OSNAP_W, Davis Strait at 67^∘ N, and the exit of Hudson Bay at 68^∘ W (Fig. 12) confirm this balance. On average, net boundary inflow (inflow − outflow) acts to increase the volume of water lighter than σ₀=27.70 kg m⁻³ and reduce the volume of denser waters. This net inflow is primarily balanced by the surface heat loss, cooling the lighter inflow waters and hence transforming them to higher densities before they flow out of the Labrador Sea. The contribution to transformation by surface freshwater fluxes (including river inflow, ice melt/formation, and precipitation-evaporation) is much smaller than the cooling. There is also a small change in the relative volumes of water within the Labrador Sea, with increasing volumes of waters in most density classes balanced by decreasing volumes of the densest waters.

Figure 12 VIKING20X Labrador Sea mean water mass balance. Positive values of volume transformation are transformations across the isopycnal towards lighter densities. Net production of water of a particular density is therefore shown as negative gradients and net removal as positive gradients (i.e. surface heat flux acts to remove lighter waters by cooling – positive gradients, producing denser waters – negative gradients). The primary balance is between net inflow of lighter water and conversion to denser water by surface heat flux. Solid lines are average annual mean transformations, 1990–2019; the shaded area shows ±1 standard deviation in the annual mean.

6.2 Changing water mass balance in the upper Labrador Sea

From water mass transformation theory, we know that the net rate of inflow water ρ<ρ₀ is balanced by the rate of change of volume with ρ<ρ₀ and the rate of transformation of water across the ρ=ρ₀ isopycnal. This transformation rate is a function of the surface heat and freshwater flux at the ρ=ρ₀ isopycnal outcrop and the transformation by interior mixing across the ρ=ρ₀ surface. Figure 13 shows how this balance evolves over time for lighter waters of σ₀<27.65 kg m⁻³. We also further split the net inflow into “inflow” and “outflow” components by assuming the inflows to be through Davis Strait, Hudson Bay, and across OSNAP_W north of 57^∘ N (mostly the west Greenland Current), and all outflow to occur across OSNAP_W south of 57^∘ N. It is the variation over time of this “outflow” in the southern half of OSNAP_W we are interested in.

Figure 13 Labrador Sea water mass balance time series for σ₀<27.65 kg m⁻³. The inflow (a, light orange) shows no clear trend after 1998, whereas the outflow (a, purple) shows a clear increase. Summing the inflow and outflow shows a reducing net inflow of lighter waters to the Labrador Sea after 1998 (blue line, a, b). The transport leaving the Labrador Sea and heading to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$ from the particle tracking is included in (a) (dashed black line) for comparison. Note the slightly different scale on the right-hand y axis for this line to help comparison, since there is not a 1:1 correspondence between water with σ₀<27.65 kg m⁻³ leaving the Labrador Sea and water tracked to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$. In (b), this reduced net inflow of light water must be balanced by water mass volume transformations – heat flux out of the ocean (orange/red line), freshwater flux (light blue) – or show as a changing volume in the volume tendency term (green). Time series are generated from 5 d VIKING20X output that was low-pass filtered with a Butterworth filter and an 18-month cutoff.

The σ₀=27.65 kg m⁻³ density level is chosen because, in VIKING20X, it represents the approximate density lighter than which most water exiting the Labrador Sea is transported across to OSNAP_E in the upper layer. Coincidentally, it is close to the division between the upper and lower limbs of the meridional overturning at the OSNAP line. We will use “upper limb water” and “lower limb water” as shorthand for water lighter and denser than σ₀=27.65 kg m⁻³, respectively. The precise choice of density level between σ₀=27.65 kg m⁻³ and σ₀=27.70 kg m⁻³ does not affect the conclusions. The density σ₀=27.70 kg m⁻³ is often considered to be the lighter density limit of the deeper LSW intermediate water masses in the North Atlantic.

On interannual and longer timescales, after 2000, the modelled inflow of upper limb water to the Labrador Sea shows some variability but no trend. The outflow of these waters after 2000, while always being less than the inflow, has an increasing trend (Fig. 13a) closely resembling the timing of, but slightly larger in amplitude than, the increase of the Labrador Sea origin transport observed in the particle tracking. The amplitude difference is because there is not a 1:1 correspondence between Labrador Sea outflow with σ₀<27.65 kg m⁻³ and water parcels tracked to OSNAP${}_{\mathrm{E}-\mathrm{37}\mathrm{W}-\mathrm{500}\mathrm{m}}$.

Hence, the decline in net inflow (difference between inflow and outflow) of these upper limb waters after the year 2000 is the result of increasing outflow and must be balanced either by reduced transformations across the σ₀=27.65 kg m⁻³ isopycnal or a reduction in the volume of water with σ₀<27.65 kg m⁻³ (negative volume tendency) in the Labrador Sea (Fig. 13b). On these multiyear timescales, the declining net inflow (and therefore the increasing Labrador Sea outflow) of upper limb waters is found to be due to reduced transformation to denser layers by surface heat flux (i.e. reduced cooling) along the σ₀=27.65 kg m⁻³ isopycnal outcrop. The contributions from freshwater flux (Fig. 13b) and diapycnal mixing are smaller. The contribution from volume tendency shows variability around a mean near zero, but integrated over time it leads to important longer-timescale changes in isopycnal depths in the Labrador Sea and Labrador Current (see Sects. 5.3.3 and 6.3).

Prior to 2000, there is a period of increasing inflow of upper limb waters in the West Greenland Current (Fig. 13a). This increase is balanced by increasing transformation to lower limb densities within the Labrador Sea by surface heat flux, resulting in little change in outflow of upper limb waters.

In summary, in the VIKING20X model, reduced cooling in the Labrador Sea between 2000 and 2013 leads to reduced transformation to denser water and so to increased transport of lighter upper waters out of the Labrador Sea. This increased transport, when mixed with the warmer, saltier waters of Gulf Stream origin in the NAC, increases the relative contribution of fresher, cooler Labrador Sea origin waters and acts to freshen and cool the transports through the eastern subpolar North Atlantic. Note that the increasing lighter water outflow occurs at a time of relatively weak and constant (or further weakening) SPG strength in VIKING20X (Biastoch et al., 2021). So, rather than being associated with increased subpolar gyre circulation (which is mostly barotropic), the increased outflow of lighter waters is associated with reduced outflow of denser waters below.

6.3 Comparison with observations

In Sect. 5.3.3, we proposed isopycnal deepening in the outflow region as the mechanism increasing volumes of lighter waters exiting the Labrador Sea. Here (Fig. 14), we compare time series of isopycnal depths in VIKING20X with those calculated from observations of isopycnal depths in the Labrador Sea outflow region. Between the early 1990s and about 2011–2012, the low-pass-filtered depth of isopycnals in the range σ₀=27.65–27.70 kg m⁻³ more than doubles (Fig. 14) in both the model and observations. After 2013, with the resumption of deeper convection in the Labrador Sea, the annual mean isopycnal depths begin to decrease. While there is good agreement between isopycnal depth changes in the model and observations on decadal timescales, at shorter interannual to pentadal timescales the agreement decreases, in particular for the lighter layers in the Labrador Sea boundary region. There are fewer observations in this region, possibly leading to the higher interannual variability observed.

Figure 14 Time series of low-pass-filtered (2-year cutoff) mean isopycnal depth from the model (red lines) and observations (black lines) in the Labrador Sea outflow region (see Fig. 1). (a) Labrador Sea boundary region, bed depths of 750–3275 m over the Labrador Slope. (b) Sea bed deeper than 3275 m in the central Labrador Sea. Note the slightly different isopycnals plotted. Isopycnals in the range σ₀=27.65–27.70 kg m⁻³ are seen to deepen across the region between the early 2000s and the early 2010s before shallowing after 2012. The deepening of the isopycnals is associated with increased outflow of lighter upper waters from the Labrador Sea.

Using the same observational dataset, Yashayaev and Loder (2017) studied the LSW (Labrador Sea Water). The LSW typically forms the intermediate water mass (approximately σ₀=27.74–27.80 kg m⁻³) immediately below the upper layer waters we are interested in here. The detailed repeat observations described show cumulative surface heat fluxes in the Labrador Sea after the early 1990s that lead to increased temperatures in the depth range 200–2000 m (Yashayaev and Loder, 2017, Fig. 3), resulting in the 2012–2016 LSW class being one of the deepest and most persistent ever observed. The recent work (Volkov et al., 2022) also shows rising sea levels in the Labrador Sea between 1994 and 2010 associated with increased ocean heat content, due to reduced surface heat flux out of the Labrador Sea.

The preceding observational analysis, combined with the previously published observations from the region, provides convincing evidence that reduced heat loss over a period of over 10 years in the Labrador Sea led to deepening upper layers. Comparison of VIKING20X-modelled Labrador Sea isopycnal depths with these observations therefore supports that VIKING20X realistically reproduces isopycnal deepening on pentadal to decadal timescales.

7 Conclusions

We have used the high-resolution, eddy-rich, nested ocean–sea-ice model VIKING20X to explore the causes of the exceptional freshening (and cooling) event observed in the surface 500–1000 m of the eastern subpolar North Atlantic in the years from 2012 onwards. We find a major cause to be reduced surface heat loss from the Labrador Sea during the decade preceding the exceptional freshening (Fig. 15).

Figure 15 Diagram illustrating the mechanism for reduced surface heat loss leading to increased outflow from the Labrador Sea: (a) late 1990s/early 2000s for a period of higher heat loss (zig-zag red lines); (b) early 2010s for the period of reduced heat loss (small zig-zag red lines). The upper panels show water mass densities (the less-dense upper layer is shown as clear and denser intermediate and deep waters are dark grey) and mean velocities across OSNAP_W (red is northward and blue is southward). The lower panels present the number of backward tracks crossing OSNAP_W and show the contrasting transport pathways to OSNAP_E. During the two periods, the inflow of light waters along West Greenland (1) remains similar in thickness and velocity structure. Reduced surface heat loss (2) in the freshening period (b) leads to less transformation from light to deep water, deepening the less dense layer (3) and so increasing the outflow of the light waters on the Labrador Shelf (4). The lower panel shows how the proportion of water from the Labrador Sea (blue) increases and pushes eastward at OSNAP_E during the freshening event (b). The Gulf Stream origin water decreases in volume and is confined more to the eastern end of the OSNAP section. While total transport from the Labrador Sea is higher in (b) than (a), transport via the loop path southward along the eastern US has decreased in (b). Most of the increased flow from the Labrador Sea is water circulating in the SPG.

Surface heat loss in the Labrador Sea transforms lighter inflow from the West Greenland Current to denser water which flows out in the Labrador Current system. However, not all the light surface inflow water is converted to denser intermediate and deep water; some flows out of the Labrador Sea in the upper few hundred metres at densities σ₀<27.65 kg m⁻³. Reducing heat loss in the Labrador Sea can therefore lead to increased volumes of lighter waters remaining in the upper layers, deepening the isopycnals and increasing the outward transport of lighter upper-layer waters in the Labrador Current system.

The relatively fresh and cold upper layer waters flowing out of the Labrador Sea combine with warmer, saltier upper-layer waters from the subtropics. They then form the upper waters of the eastern subpolar North Atlantic together. The TS characteristics of the resulting waters are largely governed by the ratio of volumes of water from the two major source regions in the mix. In this way, reduced Labrador Sea surface heat loss, outward transport of increased volumes of lighter upper layer waters from the Labrador Sea, are found to lead to freshening and cooling downstream in the eastern subpolar North Atlantic.

The quantitative particle tracking allows us to estimate contributions to the exceptional freshening and cooling. For the salinity, we find that increasing volume transport of light waters out of the Labrador Sea accounts for 60 % of the freshening, with 27 % being due to decreasing northward transport from the Gulf Stream (due to weakening modelled AMOC). The small remainder is a combination of along-track external fluxes and changes in TS properties in the origin regions.

For temperature, the picture is more complicated. Half of the overall mean cooling of the water leaving the combined origin regions is mitigated along-track by reduced heat loss to the atmosphere from the originally cooler water. Hence, while surface salinities in the eastern subpolar North Atlantic have shown extreme freshening (the freshest in 120 years Holliday et al., 2020), temperatures, though cooler than average, have shown a relatively smaller decrease. The net cooling of water leaving the combined origin regions is, as for salinity, driven mostly by increasing volumes of water originating from the Labrador Sea (over 55 %) and decreasing volumes of water from the Gulf Stream (20 %). In contrast to the freshening, a notable contribution to the cooling is also made by the reduced mean temperature of water from the Gulf Stream at their origin (15 %).

Bryden et al. (2020), using heat and freshwater budgets, find that reduced northwards fluxes from the subtropics to the subpolar North Atlantic can explain the bulk of the observed eastern subpolar North Atlantic cooling and freshening. In contrast, while in the current analysis these processes make a notable contribution to the exceptional freshening and cooling, we find the largest contribution to be made by variability in fluxes from the Labrador Sea.

We have considered previously published suggestions of subpolar gyre-scale mechanisms driving the exceptional freshening and cooling, specifically subpolar gyre strengthening and/or expansion (see, e.g. Bersch, 2002; Sarafanov, 2009; Bersch et al., 2007; Hátún et al., 2005; Koul et al., 2020; Kenigson and Timmermans, 2020) and the rearrangement of heat and freshwater in the upper subpolar gyre caused by changing pathways of LC-Arctic water (Holliday et al., 2020) driven by wind stress curl. Our methodology, calculating the mean salinity and temperature of northward transports across the whole upper eastern SPNA (rather than area-average salinities and temperatures or focusing on a subregion), precludes pure subpolar gyre eastward expansion (with no new transport) as a cause of the overall cooling and freshening (though it could contribute to the redistribution of heat and salt transport along the OSNAP_E). Examining the subpolar gyre strength (SPGI, Fig. 5), in VIKING20X-JRA-short, subpolar gyre strength increases from 2014 after a period (2000–2013) with interannual variability but no significant trend. While the strengthening subpolar gyre coincides with the strongest freshening at OSNAP_E, time series of salinity and temperature (Fig. 7) show that this freshening and cooling starts before 2012. Time series of the variability of transports leaving the Labrador Sea and Gulf Stream origins (Fig. 5) show that the changed relative proportions of waters in the mix are set by 2013, before the subpolar gyre strengthens. So, while subpolar gyre strengthening very likely contributes to the extreme of the freshening and cooling in 2014–2016, the processes leading to most of the freshening and cooling precede the subpolar gyre strengthening.

Our work partially supports the idea that changing pathways of water of Labrador Sea origin contribute to the freshening and cooling of the ESPNA (Holliday et al., 2020), in that we see a reduction of loop path transport in favour of an increase in direct path transport. This switch to a shorter path can produce a temporary increase in the proportion of fresher, colder Labrador Sea origin water in the mix, and the shorter, more northerly, path reduces along-track heat loss.

While the above processes probably all contribute, making this exceptional freshening event the strongest for 120 years, our analysis suggests that the largest contribution is from the increased outflow of lighter upper-layer waters from the Labrador Sea due to deepening isopycnals. This outflow variability appears largely independent of both subpolar gyre strength and current pathways, relying on changing isopycnal depths rather than the surface currents measured by the sea-surface-height-based SPGI. Yashayaev and Loder (2017) show close correspondence between the annual density (and hence isopycnal depth) time series for upper intermediate waters in the central Labrador Sea and multiyear cumulative atmospheric forcing indices (of surface heat flux and winter NAO). This agreement provides support for a preconditioning role of cumulative atmospheric forcing (a tendency for more negative NAO states, typically accompanied by reduced heat loss over the SPNA) over the decade preceding the exceptional freshening, leading to the peak in outflow of lighter upper-layer waters (about 2014 in VIKING20X, Fig. 5b) and the subsequent (2016–2017) peak freshening in the eastern SPNA.

We further find the increased volumes transported from the Labrador Sea to the upper 500 m of OSNAP_E to be primarily water recirculating in the subpolar gyre rather than Arctic waters of Hudson Bay, Davis Strait, or Greenland Sea origin. These Arctic waters form the Labrador Coastal Current and the LC-Arctic outflow explored by Holliday et al. (2020) and provide a possible indirect route for fluxes out of the Arctic Ocean to influence the hydrography of the SPG (Kenigson and Timmermans, 2020). The increased recirculation in the subpolar gyre, in both absolute volume terms and as a proportion of the surface subpolar gyre circulation, implies longer residence times in the upper layers of the subpolar gyre. To first order, water crossing northwards across OSNAP-E in the upper layers does one of three things: it crosses the GSR (where it is made denser and returns in deeper layers); it is made denser within the SPG and sinks to intermediate and deep layers; or it recirculates in the upper SPG. To increase the proportion recirculating, one of the others needs to reduce. The most likely candidate is reduced conversion to denser water in the SPG, as exchanges across GSR appear fairly constant from observations. The implications of this for subpolar and overturning dynamics remain the subject of further work.

It is interesting to contrast these results with those of Robson et al. (2016), who also propose a link between Labrador Sea densities and cooling of the upper eastern subpolar North Atlantic. But, while Robson et al. (2016) link the eastern subpolar North Atlantic cooling to record low densities (lighter water) in the deep Labrador Sea (in turn reducing the transport of warm water from the south associated with decreased AMOC), here we find that the link is to increased volumes of lighter water in the upper water column.

The model results have enabled us to form a hypothesis describing a possible mechanism linking Labrador Sea surface heat fluxes with salinity and temperature changes in the eastern subpolar North Atlantic. Initial comparison of the model output with observational data adds some support for this hypothesis. VIKING20X is shown to reasonably reproduce both the strength and the spatial structure of the exceptional freshening event shown in EN4 analyses (Fig. 1). Model–data comparison of vertical structure in the Labrador Sea outflow region (Fig. 14) – the region critical to our hypothesis – shows VIKING20X reproducing observed scales and timings of isopycnal deepening and shoaling in the upper 500 m.

Deepening isopycnals in the upper 500–1000 m of the Labrador Sea outflow, and the resulting increased outflow of lighter waters, has been shown here to lead to freshening and cooling of the upper eastern subpolar North Atlantic. The increased layer thickness of surface waters in the outflow must reduce the thickness of intermediate and/or deep outflows below. Both of these deeper layers also export water which is relatively fresh and cool compared to other North Atlantic water masses of similar density. While a full freshwater and heat budget of the subpolar gyre is beyond the scope of the current work, this suggests that the exceptional freshening and cooling of the upper 500–1000 m of 2012 to 2016 may, at least in part, represent a change in the vertical distribution of heat and freshwater in the eastern subpolar North Atlantic. This would complement the horizontal redistribution within the subpolar gyre proposed by Holliday et al. (2020), and the south–north redistribution of Bryden et al. (2020). Some support for the vertical redistribution of heat and salt in the subpolar gyre can be seen in the vertical sections of model temperature and salinity anomalies presented here (Fig. 1e, f) and the OSNAP observational sections (Holliday et al., 2018, Fig. 6), which both show warming and salinification below the thermocline.

While the VIKING20X model used here contains a freshening and cooling signal in the upper eastern subpolar North Atlantic with many similarities to the observed signal, the correspondence is, of course, not exact. VIKING20X slightly overestimates the magnitude of the exceptional freshening event (by up to 25 %, Figs. 1, 2) while also having a long-term underestimate of the AMOC strength of 10 %–20 % (Fig. 9). Consider an extreme scenario where all the overestimate of freshening and cooling is due to the overestimate of the variability of transports of Labrador Sea origin, and the contribution from the variability of transports of Gulf Stream origin is 10 %–20 % too low. This gives rough estimates of a minimum contribution from Labrador Sea origin transport variability of 50 % (freshening) or 45 % (cooling) and a maximum contribution from Gulf Stream origin transport of 40 % (freshening) or 25 % (cooling). We must further caveat the results presented here based on some limitations of the Lagrangian tracking methods. We also do not assess possible biases due to eddy effects and diffusion, which are not represented in the particle tracking, or due to the streamtube assumption of volume transport remaining constant along particle tracks.