Local versus far-field control on South Pacific Subantarctic mode water variability

In the South Pacific Subantarctic mode water (SAMW) formation region, central and eastern pools of SAMW have been found to be linked to winter mixed-layer thicknesses that vary strongly interannually and out of phase across the basin. This mixed-layer variability is associated with peaks in sea level pressure variability at a quasi-stationary anomaly situated between the two pools. To investigate how surface forcing, as well as the propagation of upstream anomalies, affects the formation of these SAMW pools, a set of adjoint sensitivity experiments with a density-following feature are conducted. Adjoint sensitivities reveal that local cooling can lead to an increase in the SAMW pool volume through mixed-layer-depth changes and the lateral movement of the northern boundary of the pool. In addition, upstream warming along the Antarctic Circumpolar Current can lead to an increase in the SAMW pool volume through lateral density surface movement shifting the southern boundary polewards. The density properties are advected from upstream to the downstream pool over 1 year. Optimal conditions for SAMW formation involve a combination of local cooling and upstream warming of SAMW formation sites. Hence, South Pacific SAMW variability is particularly sensitive to atmospheric modes which lead to a dipole in heating across the formation sites.

1 Introduction

Pools of Subantarctic mode water in the South Pacific Ocean are climatically important due to their role in the uptake and transport of heat, carbon, and nutrients around the global oceans (Williams et al., 2023). The overturning circulation in the Southern Ocean has a two-cell structure (Marshall and Speer, 2012). The upper cell involves the ventilation of circumpolar deep water (CDW) to form northward-flowing mode and intermediate waters. The deeper cell is formed via upwelling CDW being transformed through brine rejection during ice formation. The overturning circulation is a combination of a time-mean circulation and a time-varying eddy circulation.

Subantarctic mode water (SAMW) is a subsurface water mass that can be identified by its low absolute value of potential vorticity (McCartney, 1977), where a water mass with low potential vorticity is a vertically well mixed water mass. SAMW is formed in areas of deep winter mixed layers through surface ventilation and winter heat loss on the northern side of the Antarctic Circumpolar Current (ACC) (McCartney, 1982). The springtime shoaling of mixed-layer depths leaves behind a signal of well-mixed water with low potential vorticity, known as SAMW, and is subsequently advected away from the formation sites. If not subsequently obducted in later winters, SAMW circulates into the subtropics where it may remain for many decades (Jones et al., 2016). SAMW is found in the Southern Ocean in all three basins to the north of the ACC, with the Pacific Subantarctic mode waters being the densest and the Atlantic being the lightest (Sallée et al., 2010). SAMW is important as it is the main pathway of heat and carbon uptake from the Southern Ocean into the global ocean. The Southern Ocean dominates global ocean heat uptake (Williams et al., 2023) and contributes 40 % of CO₂ uptake historically (Frölicher et al., 2015). SAMW has been found to be thickening (Gao et al., 2018), which is contributing to the bulk of recent ocean warming (Li et al., 2023). As the thickness of Subantarctic mode waters is seen to be very important, density coordinates are a better fit than depth coordinates for assessing changes in the formation rates and properties of mode waters.

Due to the advent of the Argo time series, there are now orders of magnitude more ocean observations available than at the turn of the century, particularly in the Southern Ocean (Johnson et al., 2022). An Argo time series allows for interannual and internal variability in SAMW regions to be considered (Li et al., 2021). In the Pacific and Indian Subantarctic mode water formation regions, the central and eastern pools have been found to have wintertime thicknesses that vary strongly interannually and out of phase across each basin (Figs. 1 and S1 in the Supplement) (Meijers et al., 2019). This thickness variability is associated with changes in atmospheric forcing, as evident from peaks in the variability in sea level pressure between the central and eastern pools in each basin. Therefore, each basin-wide SAMW pool can be further split into a central and eastern pool (Cerovečki and Meijers, 2021). Cerovečki and Meijers (2021) found that these dipoles in the thickness of the pools within each ocean (Indian and Pacific) basin were driven by surface heat flux differences on either side of the basin, over each pool, and that these differences arose due to sea level pressure anomalies situated between the pools. For positive (negative) anomalies, warm (cool) air would be drawn from the north (south) over the central SAMW pools to the west of the sea level pressure anomaly, driving anomalous heat gain (loss) and resulting in shoaling (deeper) mixed layers and smaller (increased) SAMW layer volumes. The opposite would happen over the pool to the east of the sea level pressure anomaly, resulting in a clear mixed layer and SAMW dipole across the basin. Furthermore, Cerovečki and Meijers (2021) found that these sea level pressure anomalies tend to be correlated across the Pacific and Indian basins, as part of the wider wave-number-3 mode of pan-Antarctic sea level pressure variability. However, it was also found that this process was not always the case, particularly in the presence of strong El Niño–Southern Oscillation (ENSO) events, which primarily only affected the Southeast Pacific basin. Meijers et al. (2019) found that the large-scale atmospheric modes, Southern Annular Mode (SAM) and ENSO, mapped onto the sea level pressure anomalies and had similar fingerprints on the mixed-layer-depth patterns. This similarity suggests that SAM and ENSO anomalies may reinforce or weaken the dipole patterns depending on their relative phases. For example, positive ENSO is associated with heat gain over the central box and heat loss over the eastern box, thereby providing a west–east dipole in the thermal forcing. While these studies primarily focused on the role of the local heat fluxes in driving the SAMW formation site mixed-layer-depth variability, both Cerovečki and Meijers (2021) and Meijers et al. (2019) noted the presence of a clear 1-year lag between mixed-layer-depth properties transferring from the central Pacific pool to the eastern Pacific pool. Both studies suggest, although do not quantify, the potential for advective processes to set, or precondition, the response of the downstream SAMW pools.

Figure 1 Winter mixed-layer-depth anomalies (m) for the years 2006–2017 (shown in panels a–l respectively) from ECCOv4r4 data. The winter mixed-layer depth is calculated using the depths at which the neutral-density increase from the surface is ≥0.03 kg m⁻³. Green boxes show the central and eastern Subantarctic mode water formation pools. Black contours show sea surface height proxies (0.25, 1.0, and 1.5 m) for the fronts of the ACC. (m) Time series of detrended winter mixed-layer depth (m) within the central (blue) and eastern (red) boxes from 1992 to 2017. Dots on the x axis show years with strong differences in thickness between the central and eastern pools; magenta dots have positive central and negative eastern anomalies, whilst blue dots have negative central and positive eastern anomalies.

The main aims of this study are to understand the effects of different surface forcing on Subantarctic mode water pool volumes. It is known that heat fluxes and wind stresses are the most important factors driving Subantarctic mode water volume change; however, the local versus far-field contributions have not been quantified, despite several studies indicating that advection of upstream water masses, as well as their modification, is likely to be important in setting SAMW mixed-layer and subducted properties. To provide this quantification, adjoint sensitivity experiments are used to quantify when and where each forcing affects each Subantarctic mode water pool over the South Pacific.

Adjoint sensitivity analysis is a technique used in meteorology and oceanography to efficiently estimate the sensitivity of numerical model outputs with respect to various inputs (Errico, 1997). This process involves starting from a quantity of interest, known as the objective function, and propagating its sensitivities backward through a linearized version of the model. This allows for efficient determination of how various model variables affect the objective function, facilitating applications in data assimilation, parameter estimation, stability analysis, and sensitivity studies.

Adjoint models have been applied to estimate ocean circulation using data assimilation (Marotzke and Wunsch, 1993), as well as for adjoint sensitivity experiments, such as in a study of the sensitivity of Atlantic heat transport (Marotzke et al., 1999; Bugnion et al., 2006). Many aspects of Atlantic overturning (Köhl, 2005; Czeschel et al., 2010, 2012; Sévellec and Fedorov, 2016), such as heat and volume transport of various limbs of this circulation system (Heimbach et al., 2011; Kostov et al., 2019), have also been studied. In addition to transport, adjoint sensitivity experiments have also been applied to quantify the sensitivity of integral quantities such as volume- and time-mean heat content to local and remote atmospheric forcing (Jones et al., 2018b, 2020b; Boland et al., 2021b; Meijers et al., 2023b). More recently, the method has been applied to the sensitivity of volume (Kostov et al., 2024) and stratification (Pimm et al., 2024). If the objective function is defined as the integral of a tracer concentration, the sensitivity fields act as an adjoint tracer, revealing the source waters for the control volume (Fukumori et al., 2004; Chhak and Di Lorenzo, 2007; Song et al., 2016). It has also been used to diagnose the causes of bottom pressure variations, leading to insights into how atmospheric variability affects large-scale circulation (Fukumori et al., 2007b, 2015). More recently, adjoint sensitivity experiments have been used to inform observing system design by quantifying the proxy potential of selected observing sites (Loose et al., 2020; Loose and Heimbach, 2021). The wide range of adjoint sensitivity applications highlights its versatility.

Firstly, the ECCOv4r4 model and adjoint approach are outlined (Sect. 2.1). Pre-existing adjoint machinery uses depth coordinates (Boland et al., 2021a; Jones et al., 2020a) and does not account for mode water thickness changes. To address the aims a new density-following adjoint model feature is created (Sect. 2.2), which allows the variability in the volume of Subantarctic mode water in the central and eastern South Pacific to be considered more explicitly. In Sect. 2.3, the experimental design, which sets out the objective functions used and the horizontal and vertical areas chosen for each South Pacific Subantarctic mode water pool experiment, is discussed. Section 3 considers the resulting impacts, which are defined as the adjoint sensitivities convolved with forcing variability, to understand when and where different surface forcings impact each Subantarctic mode water pool separately. Then, the linear causal mechanisms implied by the adjoint sensitivities are tested using forward perturbation experiments, Sect. 4. Finally, the wider context of the work is set out in the discussion and conclusions.

2 Experimental design

2.1 The ECCOv4r4 state estimate and the adjoint approach

Numerical modelling experiments are used to understand the local and far-field influences of surface forcing, including heat and freshwater fluxes, and wind stresses, on the central and eastern southeastern Pacific Subantarctic mode water pools. Forward and adjoint modelling is completed using the ECCOv4r4 (Estimating the Circulation and Climate of the Ocean version 4 release 4) state estimate (ECCO Consortium et al., 2020). ECCOv4r4, which covers a 26-year time period from 1992 to 2017, is used for all numerical experiments completed in this study. ECCOv4r4 is a physically consistent ocean state estimate, which combines ocean observations with the MIT general circulation model (MITgcm) equations to give a more accurate description of the oceanic system with insight into the underlying processes (Forget et al., 2015a). Furthermore, the ECCOv4r4 state estimate has been shown to have good agreement with observations in studies of Southern Ocean mode water properties (Jones et al., 2018a; Boland et al., 2021a; Forget et al., 2015b).

In a more traditional modelling approach, forward perturbation experiments define a quantity of interest and chosen area, and the quantity of interest is then perturbed by a small amount within the chosen area. The effects of the perturbation on the model state are seen in the difference from the unperturbed control state. For adjoint sensitivity experiments, one quantity of interest, or objective function, is chosen, and the adjoint sensitivities of this function to every independent variable at each time, latitude, longitude, and depth are then calculated backwards in time (Jones et al., 2018a; Pimm et al., 2024). The resulting adjoint sensitivities are linear gradients of the objective function to all model independent variables at all positions and times. These adjoint sensitivities are used to identify potential causal mechanisms affecting ocean properties. However, a limitation of adjoint sensitivity experiments is that the adjoint sensitivities produced are only linear approximations. Previous work has suggested that the linear approximation is useful for several years to decades when using large-scale quantities of interest (Jones et al., 2018a; Boland et al., 2021a; Pimm et al., 2024).

2.2 Density-following feature

Previously, within MITgcm, a fixed vertical mask, defined by depth levels constant in time, has been used when employing adjoint sensitivity experiments. The fixed vertical levels lead to the objective function being integrated over a fixed volume. However, for examining mode waters, this approach can be limited, as mode waters span across multiple depths over time and are physically defined by density or potential vorticity characteristics, which shift with time. This depth choice leads to adjoint sensitivity experiments not being able to capture the full mode water features properly, as they do not capture the isopycnal thickness changes that are a defining characteristic of SAMW variability. To address this limitation, a density-following feature (Appendix A) is added into MITgcm that allows the objective function to be defined as an integral between two neutral-density surfaces that move in time (Campin et al., 2022). A problem that this new feature resolves is that mode water sensitivities previously had to be interpreted as information for the specific defined mode water pool combined with other water masses (Meijers et al., 2023a). However, the new density-following feature can be used to isolate specific water masses only, as defined by their density structure. Therefore, this feature will provide a more accurate interpretation of water masses from the adjoint sensitivities.

This addition of the objective function integrated over density surfaces that move vertically in time is achieved using a sigmoid function. The density-following feature is used to update the vertical mask at each time step, and the linear partial differentials to this objective function are then calculated in the adjoint step of the model run.

In addition, for the density-following feature, another feature is added to MITgcm to permit an objective function of volume between the upper- and lower-bounding-density surfaces to any model independent variables (Campin et al., 2023). The objective function of volume between the two surfaces is given as follows:

$$\begin{array}{}\text{(1)}& J={\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }t}}\underset{t_{\mathrm{1}}}{\overset{t_{\mathrm{2}}}{\int }}\underset{{\mathit{\sigma }}_{\mathrm{l}}}{\overset{{\mathit{\sigma }}_{\mathrm{u}}}{\int }}\underset{A}{\int }\mathrm{d}A\mathrm{d}\mathit{\sigma }\mathrm{d}t,\end{array}$$

where J is the objective function, Δt is the time interval between t₁ and t₂ over which the objective function is applied, σ_l and σ_u are the lower- and upper-bounding-density surfaces, and A is the area over the horizontal area chosen.

After both features were created, testing was conducted using both simplified models and the full ECCOv4r4 state estimate. After these tests were completed successfully, the density-following feature (Campin et al., 2022) and the volume objective function feature (Campin et al., 2023) were added to MITgcm using GitHub.

2.3 Defining the central and eastern pools

In the South Pacific Subantarctic mode water formation region, the central and eastern Subantarctic mode water pools have been found to have wintertime thicknesses that vary strongly interannually and out of phase across the basin (Meijers et al., 2019). The South Pacific Subantarctic mode water formation region is chosen here, as it has been the focus of several previous studies (Cerovečki and Meijers, 2021; Meijers et al., 2019; Li et al., 2021) and is also uniquely subject to ENSO forcing. The South Indian SAMW formation region is also found to have central and eastern Subantarctic mode water pools which vary in thickness, so the results found here may be more widely applicable.

The central and eastern regions are identified using the local maximum winter mixed-layer depth and low absolute potential vorticity from ECCOv4r4 (Fig. 1), which are separated by the Pacific–Antarctic Ridge longitudinally. Meijers et al. (2019) found that the wintertime thickness of the central and eastern pools of Subantarctic mode water formation vary out of phase with each other. The winter mixed-layer-depth anomalies have opposing signs in each pool which alternate each year (Fig. 1a–l). The relative phases of the dominant atmospheric modes in the Southern Ocean, SAM and ENSO, drive the variability in the SAMW volume and heat content. Expanding on this research, Cerovečki and Meijers (2021) show that the local sea level pressure anomalies are the main drivers of the variability in the SAMW volume, although SAM and ENSO map onto these anomalies. When SAM and ENSO are in phase with each other, positive enhanced westerly winds and associated cool Ekman transport lead to a deepening of the more southern eastern pool which is reinforced by enhanced southerlies. When the atmospheric modes are out of phase, there is a less coherent response in the Subantarctic mode water formation pools, meaning that the interannual dipole signal is stronger in some years than others; these years are identified by magenta and blue dots along the winter mixed-layer-depth time series (Fig. 1m). The central and eastern pool variability in the South Pacific SAMW is generally only seen in the later half of the ECCOv4r4 state estimate. Furthermore, towards the start of the ECCOv4r4 state estimate, it is not as reliable due to a lack of observations in the Southern Ocean. Characteristic density surface outcrop positioning is further to the north if the Subantarctic mode water pool has a positive winter mixed-layer-depth anomaly and further to the south if the Subantarctic mode water pool has a negative winter mixed-layer-depth anomaly (Fig. 2a, b). Over the seasonal cycle of the mixed layer, the outcrop positions of the density surfaces move north in winter and south in summer.

Figure 2 (a) Time-mean neutral density at the surface over longitude and latitude shown for central region density surfaces (26.7 and 27.1 kg m⁻³). Blue contours are for when the central pool mixed-layer-depth anomaly is in the negative phase and the eastern pool mixed-layer-depth anomaly is in the positive phase. Magenta contours are when the central pool mixed-layer-depth anomaly is in the positive phase and the eastern pool mixed-layer-depth anomaly is in the negative phase. The black box shows the central Pacific pool. (b) Time-mean neutral density at the surface over longitude and latitude shown for eastern region density surfaces (26.9 and 27.15 kg m⁻³). The black box shows the eastern Pacific pool. (c) Time-mean neutral density over latitude and depth averaged over the longitudes of the central region. The dotted blue line is the mean mixed-layer depth over central longitudes when the central pool mixed-layer-depth anomaly is in the negative phase and the eastern pool mixed-layer-depth anomaly is in the positive phase. The dotted magenta line is the mean mixed layer over the central longitude depth when the central pool mixed-layer-depth anomaly is in positive phase and the eastern pool mixed-layer-depth anomaly is in negative phase. The black lines show the central Pacific pool latitudes. (d) Time-mean neutral density over latitude and depth averaged over the longitudes of eastern region. The black lines show the eastern Pacific pool latitudes.

This information, along with Figs. S2 and S3, is used to define the horizontal and density masks that are used to create the objective functions for the following adjoint sensitivity experiments. The objective function of volume in the central and eastern Subantarctic mode water pools is given by the following respective expressions:

$$\begin{array}{}\text{(2)}& {\displaystyle }{J}_{\mathrm{C}}\left(V\right)={\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }t}}\underset{t_{\mathrm{1}}}{\overset{t_{\mathrm{2}}}{\int }}\underset{{\mathit{\sigma }}_{\mathrm{lC}}}{\overset{{\mathit{\sigma }}_{\mathrm{uC}}}{\int }}\underset{A_{\mathrm{C}}}{\int }\mathrm{d}A\mathrm{d}\mathit{\sigma }\mathrm{d}t,\text{(3)}& {\displaystyle }{J}_{\mathrm{E}}\left(V\right)={\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }t}}\underset{t_{\mathrm{1}}}{\overset{t_{\mathrm{2}}}{\int }}\underset{{\mathit{\sigma }}_{\mathrm{lE}}}{\overset{{\mathit{\sigma }}_{\mathrm{uE}}}{\int }}\underset{A_{\mathrm{E}}}{\int }\mathrm{d}A\mathrm{d}\mathit{\sigma }\mathrm{d}t,\end{array}$$

where Δt is the time interval between t₁ and t₂, which spans from January to December 2017; σ_uC, σ_lC, σ_uE, and σ_lE are the upper-bounding-density (u) and lower-bounding-density (l) surfaces for each of the central (C) and eastern (E) pools; and A_C and A_E are the latitude–longitude areas of each Subantarctic mode water pool. The bounding neutral-density surfaces are defined by the following expressions:

• σ_uC=26.7 kg m⁻³, σ_lC=27.1 kg m⁻³,

• σ_uE=26.9 kg m⁻³, σ_lE=27.15 kg m⁻³.

The horizontal area of each pool is defined by the following expressions:

• A_C = 180–116° W, 40–55.3° S,

• A_E = 116–72° W, 47–62° S.

An experiment for each southeastern Pacific Subantarctic mode water pool is then run individually to give the outputs of the adjoint sensitivity of volume to each independent variable.

3 Impacts

Impacts, defined as the product of adjoint sensitivities and the variation in the input variable, are the focus of this study. For example, the surface heat flux impacts on volume, I_ℋ, are calculated from the product of the adjoint sensitivities of volume to heat flux and the heat flux anomalies:

$$\begin{array}{}\text{(4)}& I_{\mathcal{H}}={\displaystyle \frac{\partial V}{\partial \mathcal{H}}}\times (\mathcal{H}-\stackrel{\mathrm{‾}}{\mathcal{H}}),\end{array}$$

where ℋ represents surface heat flux, $\partial V/\partial \mathcal{H}$ represents the sensitivity of volume to heat flux, and $(\mathcal{H}-\stackrel{\mathrm{‾}}{\mathcal{H}})$ represents surface heat flux anomaly relative to the multiyear mean annual average, which is calculated including the seasonal cycle. In ECCOv4r4, a positive surface heat flux represents a heat loss from the ocean, whereas a negative surface heat flux represents a heat gain by the ocean. The impacts of freshwater flux, eastward wind stress, and northward wind stress are calculated in the same way using their respective variables. An advantage of using impacts over adjoint sensitivities is that the influence of different surface forcings can be more easily compared to each other, as they are (1) of the same unit and (2) more “realistic” in terms of relative magnitude. Both adjoint sensitivities (Figs. S4 and S5) and impacts (Figs. S6 and S7) are useful to understand processes and provide complementary information as to the intrinsic sensitivity of volume and the extent to which surface forcing affects volume.

3.1 Quantifying the relative importance of surface forcing

To determine the relative importance of the upstream versus local impacts on each control volume, the impacts are separated by region based on whether the impacts are within the adjoint control volume or anywhere else in the global ocean, referred to as local and nonlocal respectively, such that

$$\begin{array}{}\text{(5)}& I_{\mathrm{total}}=I_{\mathrm{local}}+I_{\mathrm{nonlocal}},\end{array}$$

for each different forcing variable considered. Separating the impacts into local and nonlocal components helps to quantify which regions are most important for each surface forcing (Fig. 3). For Fig. 3, both I_local and I_nonlocal are integrated over their respective domains, to give a time series showing the average evolution of the impact over time for each region. Total impacts, I_total, may be seen in Fig. S8 in order to better compare the relative magnitude of each separate impact.

Positive heat flux impacts represent areas where either cooling leads to an increase in Subantarctic mode water volume or warming leads to a decrease in Subantarctic mode water volume. Negative heat flux impacts represent areas where either cooling leads to a decrease in Subantarctic mode water volume or warming leads to an increase in Subantarctic mode water volume. These signs carry over for the other independent variables, freshwater flux and wind stresses. The size of the volume change that can be inferred from the heat flux impacts can be considered over one grid cell. Impacts are given per grid cell for a snapshot of a 14 d period with units of the objective function per units of the independent variable. It is assumed that a snapshot of a 14 d period is similar to the average of the 14 d period. In order to determine the volume change that is ultimately caused by a change in surface forcing at a specific lead time, impacts are considered at each grid point. The net effect is found by integrating over all grid cells and lead times.

Surface heat flux has the largest impact on Subantarctic mode water volume in the central Pacific region initially (Fig. 3a). Preceding winter (July–August–September) local heat flux is most important, with the far-field, or upstream, contribution only providing a small impact (Fig. 3a). Surface heat flux is also the most important impact on Subantarctic mode water volume in the eastern Pacific region (Fig. 3e). In contrast to the central pool, the nonlocal contribution becomes very important after a 1-year lead time for the eastern pool. At a 1.5-year lead time, the most important impact on the eastern Subantarctic mode water pool is upstream heat flux, which is positioned within the central Subantarctic mode water pool (Fig. 4e). At a 1.5-year lead time, the central pool has local impacts of size that are approximately $+\mathrm{5}\times {\mathrm{10}}^{\mathrm{12}}$, whilst the nonlocal impacts are less than 𝒪(10¹²). Compared to the interannual and seasonal variability in the volume of the central mode water pool, the impacts are of the same order of magnitude, showing that heat flux is very important in the formation of the mode water pool. In contrast to the central pool at a 1.5-year lead time, the eastern pool has local impacts of less than 𝒪(10¹²) but has nonlocal impacts of size that are approximately $-\mathrm{2}\times {\mathrm{10}}^{\mathrm{12}}$. The negative sign of the nonlocal impacts indicates that the central region is where warming can lead to an increase in volume in the eastern pool; the underlying mechanisms of this process are explored in Sect. 4.2. The interannual and seasonal variability in the volume of the eastern mode water pool is of 𝒪(10¹²), showing that other impacts are also important in the formation of the eastern pool. This analysis helps to quantify the previously suggested role of advective processes in preconditioning the layer thickness of the downstream eastern Subantarctic mode water pool (Cerovečki and Meijers, 2021; Meijers et al., 2019).

Eastward and northward wind stresses also play a role at lead times of 3 years or more (Fig. 3b, c). Local and nonlocal eastward wind stress are both equally important in impacting central Pacific Subantarctic mode water volume at up to a 3-year lead time when nonlocal eastward wind stress dominates (Fig. 3b). At an approximately 3-year lead time, nonlocal eastward stress also dominates over the surface heat fluxes. Local and nonlocal northward wind stress plays a small role in impacting the volume of the central Subantarctic mode water pool (Fig. 3c). Similarly, the eastward wind stress impacts are larger than the northward wind stress impacts for the eastern Subantarctic mode water pool (Fig. 3f, g). After a 3-year lead time, the total eastward wind stress is a more dominant contribution over nonlocal heat fluxes for the eastern Subantarctic mode water pool (Fig. 3f). Freshwater flux has a very small impact on the volume of both the central and eastern Subantarctic mode water pools (Fig. 3d, h).

Initially, local heat flux is the dominant impact on volume, as seen by the amplitude of its impacts. At lead times longer than 1.5 years, the central pool continues to be dominated by local heat fluxes, but the eastern pool is dominated by upstream heat fluxes. This shows the importance of upstream heat fluxes in the eastern Subantarctic mode water pool due to advection from the upstream central pool to the downstream eastern pool via the mean flow. The fact that no such nonlocal dominance occurs within the central pool indicates the relative dynamic isolation of the Pacific from the Indian formation sites (Cerovečki and Meijers, 2021). In contrast, the volume of the eastern pool is strongly impacted by the propagation of previously formed SAMW and mixed layers from the central Pacific in the previous winter. This 1.5-year lead time aligns with the timescales hypothesized by Meijers et al. (2019) and represents the first quantification of this process. After a 3-year lead time, local and nonlocal eastward wind stresses play an important role in impacting the volume of both Subantarctic mode water pools.

During the first year lead time in the central region, the maximum overall impact is made up of the following contributions: 88 % from surface heat flux, 6 % from eastward wind stress, and 5 % from northward wind stress. Of this total surface heat flux forcing for the central region, 94 % is from local sources and 5 % is from nonlocal sources. During the first year of lead time in the eastern region, the maximum impact is made up of the following contributions: 43 % from surface heat flux, 43 % from eastward wind stress, and 13 % from northward wind stress. Of this total surface heat flux forcing for the eastern pool, 85 % is from local sources and 14 % is from nonlocal sources. Of this total eastward wind stress forcing for the eastern pool, 61 % is from local sources and 38 % is from nonlocal sources. Hence, the impact is most strongly affected by local forcing, and the contribution of nonlocal forcing varies with the type of forcing and the region.

Figure 3 (a) Impacts of surface heat flux for the central pool at the surface (anomalies multiplied by sensitivities, m³) split into local and nonlocal contributions. Panel (b) is the same but for eastward wind stress (m³). Panel (c) is the same but for northward wind stress (m³). Panel (d) is the same but for surface freshwater flux (m³). Panels (e–h) are the same but for the eastern pool. Note that the y axis for each separate impact is different.

3.2 Effect of surface heat flux on the volume of waters within density outcrops

As surface heat flux is revealed to be the dominant impact on the volumes of both the central and eastern Subantarctic mode water pools, the present subsection focuses on the role of the surface heat flux to reveal the underlying mechanism.

Latitude–longitude wintertime snapshots of heat flux impacts are considered at the surface for each Subantarctic mode water pool (Fig. 4). Figure 4 compares the heat flux impacts with the position of the bounding density surface outcrops. Initially, the impacts are mainly positive along the outcrop of the upper-bounding-density surfaces and negative along the outcrop of the lower-bounding-density surfaces (Fig. 4a, b). On the northern side of the boxes, ocean heat loss drives the upper isopycnal northward, increasing Subantarctic mode water volumes. On the southern side of the boxes, ocean heat loss also drives the lower isopycnal northward, but this makes the volume smaller within the fixed horizontal box. This suggests a dominance of the lateral movement of upper and lower bounding isopycnals rather than vertical movement involving mixed-layer depths, which is investigated more fully in Sect. 4.

After a 1.5-year lead time for the central Subantarctic mode water pool, the impacts are mainly positive and are either to the north of the upper bounding surface or in between the upper and lower bounding surfaces (Fig. 4c), with some small negative impacts mainly to the south of the outcrop of the lower-bounding-density surface. These impacts follow along the outcropping contours, along the ACC. At lead times of 1.5 years or more for the eastern Subantarctic mode water pool, the impacts are positive locally and to the north of the upper-bounding-density surface and negative upstream on the southern side of the central region (Fig. 4d). Due to the more southerly location of the eastern Subantarctic mode water pool, the outcropping density surfaces are closer together.

At a 2.5-year lead time in the central region, there are still mainly positive sensitivities to the north of the upper bounding surface, a mix of impacts in between the upper- and lower-bounding-density surfaces, and negative sensitivities to the south of the lower-bounding-density surface (Fig. 4e). For the central region, impacts within the northern flow pathway are part of the wider subtropical gyre and are, therefore, generally slower than the ACC flow. The forcing in the southern section of the horizontal box leads to the formation of Subantarctic mode waters that are seasonally capped and insulated from the atmosphere to a greater extent than the northern pathway. Therefore, the northern pathway anomalies are more likely to be atmospherically damped. Combining these factors means that longer pathways of sensitivities to the south of each horizontal box are expected. After a 2.5-year lead time in the eastern region, the impacts are mainly small when compared to the central region's impacts (Fig. 4f), and other forcings are similarly important with respect to impacting the volume of the eastern region at this time.

Local and upstream heat fluxes are the dominant contribution to the volume of both Subantarctic mode water pools for up to 3 years in the central pool and up to 2 years in the eastern pool. A longer time series of heat flux impacts may be seen in Figs. S6 and S7 for the central and eastern Subantarctic mode water pools respectively.

Figure 4 (a, c, e) Surface heat flux impacts (m³) at the surface for lead times of 0.5, 1.5, and 2.5 years. The central adjoint sensitivity target region (green) and the outcrops of the bounding density surfaces (26.7 and 27.1 kg m⁻³, black contours) are included. (b, d, f) Surface heat flux impacts (m³) at the surface for lead times of 0.5, 1.5, and 2.5 years. The eastern adjoint sensitivity target region (green) and the outcrops of the bounding density surfaces (26.9 and 27.15 kg m⁻³, black contours) are included.

The heat flux impacts are separated by their position relative to the outcrops of the bounding density surfaces, in order to show the emphasis of flow patterns on the dipole in sensitivities (Fig. 5). For the central region, heat flux impacts to the north of 26.7 kg m⁻³ and in between 26.7 and 27.1 kg m⁻³ are mainly positive. The heat flux impacts are also influenced by the seasonal cycle and are much stronger in austral winter (Fig. 5a). The negative impacts in the central region are mostly to the south of the 27.1 kg m⁻³ outcrop, and they are mainly much smaller than the positive impacts. For the eastern region, the total impacts are smaller than the central region (Fig. 5b). Positive impacts are to the north of the upper-bounding-density surface, 26.9 kg m⁻³. Negative impacts are mainly in between 26.9 and 27.15 kg m⁻³ and to the south of 27.15 kg m⁻³. The northern flow pathways are part of the wider subtropical gyre and are, therefore, generally slower than the ACC.

Figure 5 (a) Impacts of surface heat flux for the central pool at the surface (anomalies multiplied by sensitivities, m³) averaged over latitude and longitude for the components north of the outcropping density surface 26.7 (thick red), in between the density surface outcrops 26.7 and 27.1 (dashed blue), and to the south of the outcropping density surface 27.1 (thin green) over lead time (year). (b) Impacts of surface heat flux for the eastern pool at the surface (anomalies multiplied by sensitivities, m³) averaged over latitude and longitude for the components north of the outcropping density surface 26.9 (thick red), in between the density surface outcrops 26.9 and 27.15 (dashed blue), and to the south of the outcropping density surface 27.15 (thin green) over lead time (year).

The positive and negative heat flux impacts are separated over the bounding density surfaces. Positive impacts are mostly to the north of the upper-bounding-density surface outcrop for both Subantarctic mode water pools. For the central region, positive impacts are also found in between the upper- and lower-bounding-density surface outcrops. For the eastern region, negative impacts are found in between the upper- and lower-bounding-density surface outcrops. Negative impacts are mostly found to the south of the lower-bounding-density surface outcrop for both Subantarctic mode water pools. Also, positive impacts are mostly local, and negative impacts are mostly upstream of the central or eastern adjoint sensitivity region. The local and upstream pattern of the impacts is used to inform the design of the forward perturbation experiments in Sect. 4.

4 Forward perturbation experiments

Forward perturbation experiments are completed in both local and upstream regions of the central and eastern mode water pools to help understand the different physical processes at play. Impacts are mainly positive to the north of the local region, implying that if there is a cooling with a positive anomaly, there will be an increase in the mode water pool volume. Negative impacts are found upstream of both of the control regions, implying that if there is a warming with a negative anomaly, there will be an increase in the mode water pool volume. Hence, the combination of local cooling and upstream warming is the ideal condition to maximize mode water volume. These forward perturbation experiments will help explain how the mode water pool volume increases in both cases.

Forward perturbation experiments are performed in regions with reference to the sensitivity of volume in both the central and eastern pools to heat flux maps. The forward perturbation experiments shown here reference the sensitivity of the volume of the central pool to heat flux only, as the results are mostly similar. Firstly, warming and cooling perturbation experiments are conducted in an area local to the central mode water pool. This area is chosen so that it is also an area of positive sensitivity at a 1.5-year lead time (Fig. 4c). The positive adjoint sensitivities indicate that a warming anomaly will lead to a decrease in volume, whereas a cooling anomaly will lead to an increase in volume. Next, warming and cooling perturbation experiments are performed in an area upstream of the central mode water pool, located within the ACC. This area is chosen so that it is an area of negative sensitivity at a 1.5-year lead time (Fig. 4c). Negative adjoint sensitivities indicate that a warming anomaly will lead to an increase in mode water volume, whereas a cooling anomaly will lead to a decrease in mode water volume. The combination of local cooling and upstream warming increases the mode water volume. All four perturbations are of a magnitude of 100 W m⁻² and last for a 3-month time period starting on the 1 June 2016. These local and upstream responses are explored next in more detail.

4.1 Local forward perturbation experiments

At a 1.5-year lead time, the sensitivities of volume to surface heat flux (Fig. S3) are both positive and negative, with the positive sensitivities mainly to the north of upper-bounding-density surface outcrop. The perturbed region is thus defined towards the north of the central box, and the response to this perturbation remains within the northern part of the central box, as the mean circulation is weak and the ACC acts as a barrier to the advection of the response. Mode waters are generally formed on the northern edge of the ACC, so this type of perturbation is important for the short-term formation of Subantarctic mode waters.

Initially, a local warming perturbation leads to a decrease in neutral density (Fig. 6a). The mixed-layer depth also begins to shoal and the upper-bounding-density surface moves further to the south (Fig. 6c), creating a smaller volume of mode water in the central pool. The results are consistent with the expectations from the adjoint sensitivity experiments and suggest that nonlinear processes are not important on this timescale.

A cooling perturbation leads to an increase in neutral density when compared to the control experiment (Fig. 6b). Cooling leads to the mixed layer deepening and the upper-bounding-density surface moving further to the north (Fig. 6d), and the volume of mode water between the upper- and lower-bounding-density surfaces increases.

The movement of the upper-bounding-density surfaces does not persist for a long period of time (Fig. 6e, f), although this local process continues to happen over multiple winters. Here, it can be seen that local cooling over austral winter will lead to a short-term increase in the mode water pool volume, as expected. This increase in volume is due to convection as well as the lateral movement of bounding density surfaces.

Figure 6 Evolution of a negative heat flux perturbation (surface warming) and a positive heat flux perturbation (surface cooling) in a local area of positive sensitivity in the central adjoint experiment for neutral density minus the control neutral density (kg m⁻³) over depth (m) and latitude (degrees) at the chosen longitude 148° W (the midpoint of the central pool) at (a, b) the start of the warming perturbation, (c, d) 3 months after the perturbation, and (e, f) 6 months after the perturbation. Dotted magenta lines are the perturbed mixed-layer depth. Dotted green lines are the control mixed-layer depth. Thin magenta lines are the perturbed neutral-density surfaces from 26.7 to 27.1 kg m⁻³. Thick green lines are the control neutral-density surfaces from 26.7 to 27.1 kg m⁻³.

4.2 Upstream forward perturbation experiments

The upstream forward experiments are designed to reveal how upstream warming leads to an increase in the downstream volume of mode waters. At a 1.5-year lead time, the sensitivities of volume to surface heat flux (Fig. S3) are both positive and negative, with the negative sensitivities mainly to the south of the upper-bounding-density surface and within the ACC, suggesting that advection from upstream source waters may be important.

This warming perturbation is positioned upstream of the central mode water pool within the fronts of the ACC. Warming initially leads to a decrease in neutral density at the site of the warming perturbation (Fig. 7a). Warming also leads to a small shoaling of the mixed layer and laterally moves the outcrop of the lower-density surface further to south (Fig. 7c). The volume of the central mode water pool is increased as a result. Following the volume response to the perturbation as it is advected along the ACC, the mixed-layer depth quickly deepens back to the normal state, but the lateral movement of the lower-density surface persists.

The identically positioned cooling perturbation leads to an increase in neutral density, which is of an opposing sign to the warming perturbation (Fig. 7b). A total of 3 months after the perturbation, the mixed-layer depth deepens and the outcrop of the lower-bounding-density surface moves laterally northwards (Fig. 7d). The control and perturbed mixed-layer depths shoal over the seasonal cycle to the same depths, having a short-lived and local effect on the mixed layer. However, the movement of the density surfaces lasts longer and advects along the ACC, with the neutral-density changes in the surface ocean.

Hence, through lateral density surface movement, forward perturbations with opposing-sign perturbations lead to similar results with opposing signs, in accord with the adjoint sensitivity and impact analysis.

Figure 7 Evolution of a negative heat flux perturbation (surface warming) and a positive heat flux perturbation (surface cooling) in an upstream area of negative sensitivity in the central adjoint experiment for neutral density minus the control neutral density (kg m⁻³) over depth (m) and latitude (degrees) at the chosen longitude 140^∘E (upstream of the central pool) at (a, b) the start of the warming perturbation, (c, d) 3 months after the perturbation, and (e, f) 6 months after the perturbation. Dotted magenta lines are the perturbed mixed-layer depth. Dotted green lines are the control mixed-layer depth. Thin magenta lines are the perturbed neutral-density surfaces from 26.7 to 27.1 kg m³. Thick green lines are the control neutral-density surfaces from 26.7 to 27.1 kg m⁻³.

5 Discussion and conclusions

This study addresses how surface forcing controls the Subantarctic mode water volume in the Pacific sector.

Adjoint sensitivities indicate that local cooling leads to an increase in the mode water pool volume via mixed-layer deepening. In addition, local cooling also leads to the upper density surface outcrop laterally moving to the north, further increasing the mode water pool volume. The position of the local cooling is asymmetric, being located over the northern half of the central pool. The central pool is located in a more northerly position than the eastern pool, so this local, northern boundary regime is more important with respect to forcing the central pool.

Adjoint sensitivities also indicate that upstream warming along the ACC leads to an increase in the mode water pool volume. The mechanism that leads to this nonintuitive result is the southward, lateral movement of the lower-density surface outcrop. The mixed layer shoals from this warming upstream of the SAMW formation region as expected, but this property is not transferred along the ACC and the mixed-layer memory in the upstream region is weak after local seasonal shoaling occurs. This southern boundary regime is more important in the eastern mode water pool than the central mode water pool, as upstream impacts are of a magnitude of 2×10¹² for the eastern pool, whereas upstream impacts are of a magnitude of less than 1×10¹² for the central pool. The eastern pool is located further south than the central pool, and more of the eastern pool is controlled by the ACC, which explains the eastern pool impacts being controlled from upstream at lead times longer than 1.5 years. In contrast, the central pool, with its weaker connection to the ACC, is controlled by perturbations more locally. Whilst defining a water mass class by density boundaries will always make its volume sensitive to processes that move those outcropping positions, these bounding surfaces are chosen based on the observed and modelled positions of the deep winter mixed layers associated with SAMW formation. Additionally, this work demonstrates that whilst volumes are most impacted via mixed-layer-depth changes due to instantaneous forcings consistent with present literature (Meijers et al., 2019; Cerovečki and Meijers, 2021), the associated lateral movement of density surfaces is the most impactful driver of volumes due to upstream forcing.

Optimal conditions for mode water formation involve local, instantaneous cooling and a preconditioning upstream warming 1 year prior at mode water formation sites. The asymmetric relative positions of the two mode water pools relative to the ACC mean that the major atmospheric modes of variability in the Southern Ocean may exert differing impacts on the two pools. The main atmospheric modes in the South Pacific are the El Niño–Southern Oscillation (ENSO) and the Southern Annular Mode (SAM). Meijers et al. (2019) considers the instantaneous effects of SAM and ENSO on heat fluxes over the central and eastern mode water pools. Both atmospheric modes lead to a dipole in heat fluxes over the two pools. A positive SAM leads to instantaneous warming over the central box (Meijers et al., 2019), which has been shown here to lead to the lower-density surface moving laterally to the south. Over 1 year, this lateral density surface movement is transferred along the ACC and leads to an increase in the mode water pool volume in the eastern pool. A positive SAM (or sea level pressure anomaly between the pools) is also associated with cooling over the eastern pool, which leads to a local increase in volume. This preconditions the eastern pool to already be thicker before any upstream advective processes reach the pool. The positive SAM mode affects the southern boundary more than the northern boundary (Meijers et al., 2019), which leads to SAM being the dominant mode over the eastern pool, as it is located further south and more within the ACC than the central pool. Recently, SAM has been largely in a positive phase, which reinforces previous work that states that SAMW is getting thicker (Gao et al., 2018). When SAM and ENSO are in phase with each other, ENSO will also reinforce the pattern of SAMW pools becoming thicker. Therefore, this response will be stronger in the South Pacific than in other Southern Ocean mode water formation regions when both ENSO and SAM are in phase with each other.

This study reveals that the southern boundary of the eastern box is influenced significantly by the heat flux and properties of the central box approximately 1–1.5 years earlier. This advective timescale is consistent with that suggested by Meijers et al. (2019). Advection may potentially interact with the previously described dipole in mixed-layer depth to enhance SAMW formation in the eastern pool. During a dipole year, heat fluxes may be high over the central pool, driving a mixed-layer-depth shoaling and a dipole deepening of the mixed-layer depth in the eastern pool. The warming of the central pool is then advected downstream to the eastern pool the following winter, where that water mass finds a preconditioned “deep” SAMW layer and acts to enhance the SAMW volume by shifting the isopycnals polewards. In this manner, we may expect the eastern pool to be influenced by the central pool, driving deeper mixed-layer depths and, hence, stronger SAMW export. The quantification of this impact is beyond the scope of this study but is supported by the presence of deeper dipole mixed layers in the eastern pool than in the southern pool.

A limitation of this study is the linearity of the adjoint sensitivities. The adjoint sensitivities are linear, as only the first derivatives are retained in the calculation. Thus, any nonlinear processes at play in mode water formation are not captured by the adjoint sensitivities (Forget et al., 2015a). However, the adjoint sensitivities and impacts are still useful for highlighting important regions and comparing the imprint of different forcings as presented here (Sect. 3) and in other works (Fukumori et al., 2007a; Pillar et al., 2016; Jones et al., 2020a).

In principle, adjoint calculations similar to the ones shown here could be carried out in the eddy-permitting resolution Southern Ocean State Estimate (Mazloff et al., 2010). However, the higher resolution will likely shorten the timescale on which the linear adjoint is applicable, potentially limiting the tool's ability to quantify sensitivities to remote forcing. A comparable study using a high-resolution model would be a welcome addition to the literature.

The benefits of using the density-following adjoint feature is that each mode water pool can be better identified using its bounding density surfaces. Previously, adjoint sensitivity experiments (Boland et al., 2021a; Jones et al., 2020a) did not properly capture isopycnal thickness changes, which are a defining characteristic of SAMW variability. The density-following feature is used to properly isolate specific water masses only, as defined by their density structure. This prior limitation was particularly an issue in the work by Jones et al. (2016), who attempted to assess changes in exported SAMW properties and struggled to quantify the impacts of changes in properties versus any changes in the density surface positions.

A further caveat of this study is the use of a fixed horizontal box. Over the seasonal cycle, the upper- and lower-bounding-density surface outcrops can leave the fixed horizontal box, but this does not occur very often. This movement out of the box leads to a fixed volume being used, at which point the only mechanism for changing the volume is mixed-layer-depth changes, which are very small. Some aspects may still be further improved to represent changes in the SAMW pool.

In terms of other mode water formations regions, our analysis may be relevant to the South Indian mode water formation region, as central and eastern mode water pools have also been identified (Cerovečki and Meijers, 2021) and as the ACC is important in its formation (Pimm et al., 2024). The South Indian mode water pool would be a good place to test the ideas put forward here involving preconditioning, the role of asymmetry in the pool latitude when compared to the ACC position, and the role of different atmospheric modes (such as SAM).

While mode waters in the Northern Hemisphere are formed in a similar way to mode waters in the Southern Ocean, the Northern Hemisphere mode waters are generally formed in regions that are influenced by different flow regimes, such as gyres and boundary currents. As a result, these Northern Hemisphere mode waters may not be as influenced from upstream sources, although the boundary currents in the Northern Hemisphere may act in a similar way to the ACC and permit the influence of upstream forcing.

In summary, we find that SAMW formation sites are strongly influenced by local instantaneous air–sea heat fluxes. Also, upstream heat fluxes that are subsequently advected into the region play an important role via the lateral movement of density surfaces rather than changes in the mixed-layer depth. The exact location of the mixed-layer pools relative to the fronts of the ACC and air–sea flux dipoles plays an important role in their response to variability in atmospheric modes, such as SAM or ENSO, and the upstream mode water pool has the potential to significantly influence the downstream pool in subsequent years, leading to a potential net asymmetry in SAMW formation and export.

Appendix A: Density-following adjoint feature

As of checkpoint 68i (Campin et al., 2022), a new feature has been added to MITgcm so that the objective function can be integrated over density surfaces that move in the vertical in time. This process is completed using a sigmoid function:

C ---- If density mask is enabled, use it here ----
            IF ( maskC(i,j,k,bi,bj).EQ.oneRS .AND.
     &           gencost_useDensityMask(kgen) ) THEN
C            - first, calculate the scalar density
             CALL FIND_RHO_SCALAR(
     I              theta(i,j,k,bi,bj),
     I              salt(i,j,k,bi,bj),
     I              gencost_refPressure(kgen),
     O              tmpsig,
     I              myThid )
C            - subtract 1000 to get sigma
             tmpsig = tmpsig - 1000. _d 0
C            - now, tmpmsk is sigmoid times this value
             tmpsig_lower = 0.5 + 0.5*tanh(gencost_tanhScale(kgen)
     &           *(tmpsig-gencost_sigmaLow(kgen)))
             tmpsig_upper = 0.5 - 0.5*tanh(gencost_tanhScale(kgen)
     &           *(tmpsig-gencost_sigmaHigh(kgen)))
C             - update mask value based on the sigmoid function
             tmp_sigmsk = tmpsig_lower*tmpsig_upper
             tmpmsk = tmpmsk*tmp_sigmsk
            ENDIF
C ---- end of density mask (but tmpmsk is used below)

Lines 2 and 3 of the code check if the density-following function (gencost_useDensityMask) has been turned on (or not) by the user in the ECCO namelist file (data.ecco) and also check that there is a horizontal mask (maskC) present. The mask (maskC) is the chosen horizontal area of interest used in the objective function that is given by the user. Here, kgen is a reference number to make sure that all variables are taken from the same group in the ECCO namelist file, as there can be many different variable groups to give different adjoint equations. If both of these checks are true, the code continues into the next step (FIND_RHO_SCALAR) to calculate the scalar density in lines 5–10. If either of the checks fails, this section is skipped and the density-following feature is not used. The function FIND_RHO_SCALAR is a pre-existing function in the MITgcm code, and it is used at other times within the model code to calculate scalar density at different time steps. The inputs theta and salt are found at every latitude, longitude, and depth for the specific time step and are used in calculating density. The input gencost_refPressure is prescribed in data.ecco by the user. The output (tmpsig) is the temporary sigma value that is used in the sigmoid function to find the upper- and lower-density surfaces. The final input, myThid, is the time step that the model is at. In line 12, tmpsig has 1000 subtracted so that sigma is in the same form as the user inputs of gencost_sigmaLow and gencost_sigmaHigh. Lines 14–15 and 16–17 use a sigmoid function to calculate where the lower-density (tmpsig_lower) and upper-density (tmpsig_upper) surface is located respectively. Lines 19–20 update the vertical mask values (tmpmsk) that are then used by the adjoint model equations to calculate the sensitivity of the objective function to the independent variables at this time step.