Diagnosing the Atlantic Meridional Overturning Circulation in density space is critical in warmer climates

The Atlantic Meridional Overturning Circulation (AMOC) plays a crucial role in shaping the global climate system by redistributing heat and influencing large-scale climate patterns. Utilizing the AWI-CM3 model, we investigate the AMOC sensitivity to an abrupt climate change scenario (abrupt-4xCO₂) with respect to pre-industrial climate (PI), comparing the AMOC diagnosed in depth (z-AMOC) and density (ρ-AMOC) space. Water mass transformations are assessed to analyze the impact of background climate on surface-forced and interior-mixing-induced transformations. We find that both the location and magnitude of AMOC maximum are directly affected by the framework choice. In PI, the ρ-AMOC maximum is substantially stronger than that of the z-AMOC, while at 26° N the two diagnostics are nearly equivalent. Consequently, the variability of the z-AMOC maximum correlates only with that at 26° N, reflecting isopycnal flattening into constant depth levels in the subpolar North Atlantic inherent to these diagnostics. Strong AMOC weakening is observed under 4xCO₂ forcing in both frameworks until simulation year 75. Subsequently, both diagnostics reveal a weaker AMOC with an approximate strength of 7.1 Sv, although with z-AMOC displaying a slight recover towards the end of the simulation and ρ-AMOC oscillating steadily around 5 Sv. At 26° N, variability patterns remain comparable to PI, albeit with an additional ρ-AMOC weakening of approximately 2 Sv, indicating divergence between the representation of AMOC dynamics in the subtropical Atlantic within both frameworks in comparison to PI. The diagnostics in density space allow for the attribution of this further ρ-AMOC weakening to increased entrainment of fresher overflows from amplified GIN seas overturning and reduced deep convection in the Labrador and Irminger Seas. In contrast, the diagnostics in depth space only reveals reduced downwelling around the southwestern Greenland coast and along the path of the Gulf Stream, features that are more challenging to evaluate against available observations. Thus, the comparison between z-AMOC and ρ-AMOC indicates that diagnosing the AMOC in density space provides more physically meaningful information regarding the state of the water mass transformations and their contribution to ocean circulation regimes across the entire Atlantic basin, not only the subpolar North Atlantic, and especially as the climate continues to warm. These findings emphasize the importance of diagnosing AMOC in density space to better understand water mass transformations, which are concealed in depth space and to capture AMOC variability in warmer climates, across all latitudes.

1 Introduction

The Atlantic Meridional Overturning Circulation (AMOC) plays a crucial role in the global climate, accounting for approximately 25 % of global heat transport and nearly half of the deep-water formation in the ocean (Srokosz et al., 2012). The AMOC is also deemed as a global tipping point that is weakening and likely to collapse if subject to enough climate forcing, although the magnitude of this weakening or probability of collapse remain uncertain and are under heavy discussion (McKay et al., 2022; Ditlevsen and Ditlevsen, 2023; van Westen et al., 2025; Dima et al., 2025; Zimmerman et al., 2025). Such potential weakening under anthropogenic climate change has gained increasing scientific attention, given its association with potentially catastrophic climatic events, particularly affecting economically vulnerable regions (Orihuela-Pinto et al., 2022; Meccia et al., 2025; Zhang et al., 2024). Consequently, understanding the AMOC sensitivity to climates warmer than present, here defined relative to the Pre-Industrial (PI) period (1850 CE (Common Era)), is critical for assessing the potential environmental and societal impacts of future climate change and to better prepare for them (Bellomo et al., 2021).

At current background climate, the AMOC firstly consists of an upper limb of warm, saline waters flowing northward towards the subpolar North Atlantic (SPNA) (Buckley and Marshall, 2016). This upper limb is fed by two primary pathways: the “cold route”, which transports relatively cold and fresh waters from the Drake Passage in the Southern Ocean, and the “warm route”, which carries relatively warm and saline waters from the Agulhas Plateau at the boundary between the Indian and Atlantic Oceans (Rühs et al., 2019), with the latter representing the dominant contribution. The maximum AMOC overturning strength occurs within its mid-depth cell, at approximately 1000 m below the ocean surface and between 30–65° N (Matos et al., 2020). Upon reaching the northern North Atlantic, these waters lose heat to the atmosphere, become denser, sink, and return southward as a cold lower limb dominated by the North Atlantic Deep Water (NADW) (Buckley and Marshall, 2016). Furthermore, meridional overturning in the Atlantic Ocean includes the abyssal cell, sustained by the Antarctic Bottom Water (AABW) sourced from the Southern Ocean (Biló et al., 2024).

The AMOC is primarily driven by winds, buoyancy fluxes and interior mixing processes that vary over synoptic to multi-centennial timescales (Buckley and Marshall, 2016). These dynamics are strongly influenced by fluctuations in heat and freshwater fluxes throughout the Atlantic Ocean, which regulate the density of water masses transported via ocean currents (Sévellec and Fedorov, 2016). However, observing variations in buoyancy fluxes at fine resolution remains challenging due to current technological limitations, despite recent advances in AMOC fingerprint analysis and observational arrays (Frajka-Williams et al., 2019). Consequently, numerical model simulations remain indispensable for investigating AMOC variability across multiple scales. In this context, the Coupled Model Intercomparison Project (CMIP) was developed aiming at providing a standardized framework for developing climate model simulations and facilitating robust inter-model comparison. The sixth phase, CMIP6 (Eyring et al., 2016), has enabled a wide range of studies focused on evaluating the AMOC sensitivity to different forcings and scenarios and over various time and spatial scales. Concurrent with the increased number of studies concerning the AMOC response to climate change, another critical factor has gained attention: the choice of vertical coordinate system employed to derive AMOC-related output and diagnostics (Foukal and Chafik, 2024). These coordinate systems refer to the AMOC diagnostics in depth (z-AMOC) or density (ρ-AMOC) space. In depth space, the AMOC is derived from the zonally averaged streamfunction varying with latitude and depth, whereas in density space, the AMOC varies with density instead of depth. It is important to highlight, however, that ρ- and z-AMOC are not two different AMOCs but different representations of the same circulation, thus emphasizing different aspects of it.

The choice of coordinate system becomes particularly relevant when comparing the overturning estimates in the SPNA with those in the subtropical North Atlantic (STNA). In the STNA, the strong stratification and relatively flat isopycnals allow the southward limb to flow directly beneath the northward upper limb, resulting in similar estimates of magnitude and variability of ρ- and z-AMOC (Moat et al., 2025). Conversely, in the SPNA, sloped isopycnals induce horizontal separation between the upper and lower limbs, resulting in divergence between ρ- and z-AMOC in terms of both strength and variability (Foukal and Chafik, 2024). This divergence arises because zonal averaging in depth space conceals this horizontal separation, which in turn compromises the accuracy of AMOC estimates under this representation, particularly in the SPNA. Moreover, deriving the overturning circulation in density space instead of depth space is advantageous beyond studies that focus mostly on the North Atlantic and on current climate change, as constant-depth averaging can lead to spurious features such as the Deacon cell in the Southern Ocean (Döös and Webb, 1994; Stevens and Ivchenko, 1997; Speer et al., 2000) and has been linked to discrepancies between modeled and observed AMOC variability across timescales (Liu et al., 2017). Thus, diagnosing AMOC in density space yields a more continuous and physically consistent representation of the AMOC and its underlying mechanisms (Megann, 2018; Xu et al., 2018; Sidorenko et al., 2020a, 2021; Megann et al., 2021; Foukal and Chafik, 2024).

Although diagnosing AMOC in density space provides clear advantages, the establishment of ρ-AMOC as the standard diagnostic is still hindered by the research community's long-standing familiarity with z-AMOC, built over decades of studies employing the latter definition, with some articles providing a supplementary figure of ρ-AMOC remapped onto depth coordinates to facilitate comparisons between depth and density space AMOC representations (e.g. Xu et al., 2018; Tesdal et al., 2023; Foukal and Chafik, 2024). An additional caveat includes the higher computational cost associated with diagnosing ρ-AMOC (Sidorenko et al., 2021), which can discourage its implementation in studies that require long integration periods or high-resolution output. Furthermore, while the streamfunction in density space has been requested as output in CMIP6 (Griffies et al., 2016), it was not provided consistently for all experiments by all participating modelling centers (Baker et al., 2025; Jackson and Petit, 2023), which limits model intercomparison. While these barriers remain, the scientific gain through diagnosing ρ-AMOC outweighs these challenges, driving its increasing recognition in recent decades, with a strong momentum towards diagnosing AMOC in density space either in modelling studies or observational arrays at various latitudes (e.g. Frajka-Williams et al., 2023; Jackson and Petit, 2023; Fu et al., 2023; van Westen et al., 2025). In particular, observational arrays such as the OSNAP (Overturning in the Subpolar North Atlantic Program; Lozier et al., 2017) and RAPID-MOCHA (RAPID Climate Change – Meridional Overturning Circulation and Heatflux Array; McCarthy et al., 2015), already provide ρ-AMOC output and have, since their launch, changed our view on overturning in the subpolar and subtropical North Atlantic (Lozier et al., 2019; Moat et al., 2025; Frajka-Williams et al., 2023). Consequently, this dichotomy between the studies employing ρ- and/or z-AMOC frameworks introduces increased uncertainty regarding the occurrence and timing of a substantial AMOC weakening under a warming climate. At the inception of CMIP7, where both ρ- and z-AMOC are requested (Fox-Kemper et al., 2025), we see a timely opportunity to advertise the more widespread adoption of ρ-AMOC, at least as an additional, if not even the main, overturning diagnostic.

Motivated by this uncertainty surrounding the importance of diagnosing AMOC in density space under amplified climate change and across different latitudes, we conducted two simulations using the AWI-CM3 model (Alfred Wegener Institute Climate Model version 3; Streffing et al., 2022): (1) the CMIP6 abrupt-4xCO₂ scenario (4xCO₂ hereafter), characterized by quadrupled atmospheric CO₂ concentration relative to pre-industrial levels (∼1120 ppmv); and (2) the CMIP6 piControl experiment (PI hereafter) as our control climate. We employ the 4xCO₂ experiment as our warm climate instead of the CMIP6 shared socioeconomic pathway scenarios (SSPs; O'Neill et al., 2016) for future climate change because the idealized and abrupt forcing in 4xCO₂ allows a clearer attribution of AMOC changes to greenhouse-gas forcing alone. In contrast to SSPs, which combine multiple forcings and are constrained by scenario-specific socioeconomic assumptions, the 4xCO₂ framework facilitates the diagnosis of isolated physical mechanisms underlying the AMOC response. Moreover, while the 4xCO₂ simulation may appear extreme relative to the current climate trajectory, prior research has demonstrated its utility in contextualizing ongoing climate change (Yang et al., 2023) and in assessing AMOC resilience to high CO₂ forcing considering natural variability (Baker et al., 2025), as this experiment provides a clearer assessment of physical mechanisms under strong radiative forcing. Our primary objective is to determine whether the framework choice (z-AMOC vs. ρ-AMOC) affects our understanding of AMOC changes that occur in the severely perturbed warm climate state of 4xCO₂. Additionally, we aim at assessing whether ρ-AMOC becomes a more relevant diagnostics also in the STNA in warmer climates, through evaluating the degree of similarity between ρ- and z-AMOC mean state and variability at high latitudes and at 26° N. Finally, we demonstrate the potential and advantages of adopting ρ-AMOC as standard output in future model intercomparison and model-observation comparison efforts, independently from the scope of the study.

2 Methodology

2.1 Description of Model and Simulations

Utilizing the AWI-CM3 model (Streffing et al., 2022), we conduct a quasi-equilibrium pre-industrial simulation (PI) that we employ as the control climate, and a perturbed warm climate simulation with an abrupt quadrupling of PI CO₂ concentration (4xCO₂). Both simulations span 200 years and were branched off from an equilibrated 1000 year PI spin-up simulation, adhering to CMIP6 guidelines (Eyring et al., 2016), as introduced in Sect. 1. Analyses of mean large-scale processes were performed utilizing the output of the final 50 years of each simulation.

AWI-CM3 comprises both an atmospheric and an ocean and sea ice component. The atmospheric component is the open-source version of the Integrated Forecast System, version 43r3v2 (OpenIFS; Buizza et al., 2017), developed by the European Center for Medium Range Weather Forecast (ECMWF). The ocean component, FESOM2.5 (Finite volumE Sea-ice Ocean Model, version 2.5; Danilov et al., 2017; Scholz et al., 2019, 2022,  FESOM2 hereafter), developed at AWI, features a triangular unstructured mesh and the embedded Finite volumE Sea Ice Model (FESIM; Danilov et al., 2015). Relevant fluxes are coupled between ocean and atmosphere via the OASIS3-MCT4 coupler (Ocean, Atmosphere, Sea Ice, Soil Model Coupling Toolkit, version 4; Craig et al., 2017).

OpenIFS was employed at approximately 100 km horizontal resolution on a cubic octahedral grid (TCo95L91) for higher computational efficiency (Malardel et al., 2016). FESOM2 was employed with an ocean mesh configuration of one-degree nominal resolution, refined near the Equator ($\mathrm{1}/\mathrm{3}\mathit{°}$), north of 50° N (∼24 km), and near coastlines. In this configuration, FESOM2 consists of approximately 127 000 surface mesh nodes and 47 depth levels implemented under an Arbitrary–Lagrangian–Eulerian (ALE) vertical coordinate scheme, which allows for flexible vertical layer configurations. For our simulations, we employed the z^∗ vertical coordinate system (Adcroft and Campin, 2004), which scales vertical layers proportionally to sea surface height, thereby reducing spurious mixing in regions with strong stratification (Petersen et al., 2015; Scholz et al., 2019). Vertical mixing in FESOM2 is parameterized using the K-Profile Parameterization (KPP; Large et al., 1994) scheme, through which locally enhanced vertical diffusivity, set to 0.01 m² s⁻¹, is implemented to represent convection arising from local static instability. Given that the mesh employed in this study is not eddy-resolving, mesoscale eddy stirring is included via the Gent–McWilliams (GM) parameterization (Gent and McWilliams, 1990), implemented according to the explicit eddy-induced stream-function algorithm of Ferrari et al. (2010), as detailed in Danilov et al. (2017) and evaluated in Scholz et al. (2019).

To enable diagnostics of AMOC in density space, we remapped the model vertical levels onto 89 uneven density bins referenced to 2000 dbar (σ₂; kg m⁻³), ranging from σ₂=30–40 kg m⁻³, following methodologies described by Megann (2018), Xu et al. (2018), and Sidorenko et al. (2020a). The specific bins used in this study are provided in Matos (2025d), and the analysis of the ρ-AMOC sensitivity to the number of bins implemented is discussed in detail in Sidorenko et al. (2020b). The model outputs standard variables in depth space while computing transports in density space during runtime, thereby optimizing storage and computational costs (Sidorenko et al., 2021). In addition to the calculation of z-AMOC and ρ-AMOC, surface-forced diapycnal (ψ_S) and interior-mixing-induced (ψ_I) water mass transformations that are the two main components of ρ-AMOC (hereafter referred to as surface and interior transformations) were computed.

2.2 Mathematical Framework and Definitions

Our study adopts the mathematical framework defined in Sidorenko et al. (2020a) and applies the algorithms developed by Sidorenko et al. (2020b) for AMOC diagnostics on unstructured meshes. In these algorithms, the overturning streamfunctions are computed directly on the model's native unstructured mesh using a binning approach that avoids interpolation and preserves mass conservation. Specifically, mesh elements (triangles) or edges intersecting a selected latitude or density bins are identified, and fluxes normal to those sections are integrated using the native finite-volume discretization, ensuring full geometric consistency with the unstructured grid. In the following, we provide the mathematical definitions used in this study:

1. AMOC in depth space (z-AMOC):

   $$\begin{array}{}\text{(1)}& {\mathit{\psi }}_z(y,z)=\underset{\text{North}}{\overset{y}{\int }}\underset{x_{\mathrm{E}}}{\overset{x_{\mathrm{W}}}{\int }}w(x^{\prime },y^{\prime },z)\mathrm{d}{x}^{\prime }\mathrm{d}{y}^{\prime },\end{array}$$

   where w denotes the vertical velocity that includes the GM bolus component w^∗. x and y represent the zonal and meridional distances, respectively, and z represents the depth. The limits from North to y indicate the cumulative integration from a reference latitude (usually the northern boundary of the domain) down to the latitude y. The limits from x_E to x_W indicate the full integration from the eastern to the western boundary of the domain.

2. AMOC in density space (ρ-AMOC):

   $$\begin{array}{}\text{(2)}& {\mathit{\psi }}_{\mathit{\sigma }}(y,\mathit{\rho })=\underset{\text{North}}{\overset{y}{\int }}\underset{x_{\mathrm{E}}}{\overset{x_{\mathrm{W}}}{\int }}{w}_{\mathit{\rho }}(x^{\prime },y^{\prime },\mathit{\rho })\mathrm{d}{x}^{\prime }\mathrm{d}{y}^{\prime },\end{array}$$

   where w_ρ represents the diapycnal velocity across the isopycnal surface ρ, reconstructed from the divergence of horizontal flow within each density bin. In our simulations, the vertical structure of the isopycnals that are binned and converted to diapycnal velocity is dynamically modified during runtime via convective adjustment induced by the KPP parameterization. Consequently, w_ρ directly incorporates diffusive water mass transformations, which facilitates the decomposition of ρ-AMOC into contributions from surface and interior transformations, and model drift.

3. Surface-forced diapycnal water mass transformations (ψ_S; surface transformations):

   $$\begin{array}{}\text{(3)}& {\mathit{\psi }}_{\mathrm{S}}(y,\mathit{\rho })={\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }\mathit{\rho }}}\underset{x_{\mathrm{E}}}{\overset{x_{\mathrm{W}}}{\int }}\underset{\text{North}}{\overset{y}{\int }}\underset{\mathit{\rho }-\frac{\mathrm{\Delta }\mathit{\rho }}{\mathrm{2}}}{\overset{\mathit{\rho }+\frac{\mathrm{\Delta }\mathit{\rho }}{\mathrm{2}}}{\int }}{F}_{\mathrm{p}}(x^{\prime },y^{\prime },{\mathit{\rho }}^{\prime })\mathrm{d}{\mathit{\rho }}^{\prime }\mathrm{d}{y}^{\prime }\mathrm{d}{x}^{\prime },\end{array}$$

   where F_p represents buoyancy flux and Δρ is the size of the density bin. This component reflects surface buoyancy fluxes that are confined to the surface mixed layer and relates the rate of density transformation in a given latitude and density class to the surface buoyancy fluxes into that density class over its outcrop area (Walin, 1982; Speer and Tziperman, 1992).

4. The model drift, $\frac{\mathrm{\Delta }V}{\mathrm{\Delta }t}$, which corresponds to the rate of volume change over time. In our simulations, the model drift was calculated and found to be negligible around the density and latitude of the upper cell maximum (Fig. A2).

5. Interior-mixing-induced water mass transformations (ψ_I; interior transformations):

   $$\begin{array}{}\text{(4)}& {\mathit{\psi }}_{\mathrm{I}}(y,\mathit{\rho })={\mathit{\psi }}_{\mathit{\sigma }}(y,\mathit{\rho })-{\displaystyle \frac{\mathrm{\Delta }V}{\mathrm{\Delta }t}}-{\mathit{\psi }}_{\mathrm{S}}(y,\mathit{\rho }).\end{array}$$

   This component accounts for transformations forced by diapycnal mixing in the ocean interior, cabbeling and thermobaricity (Megann et al., 2021), excluding mixing caused by surface buoyancy fluxes. Therefore, it accounts for water masses initially driven by surface transformations and subsequently modified in the ocean interior.

Please note that, in our study, the AMOC is computed using vertical velocity, w, rather than the conventional approach based on meridional velocity (v; see Table A1 for the equations). We make this choice because the hydrostatic Boussinesq formulation used by FESOM2 implies incompressibility (Banerjee et al., 2024), which makes both methodologies mathematically equivalent (see Sidorenko et al., 2020b, for a detailed comparison) and because not using meridional velocity is advantageous with the spatial discretization used by FESOM2. Using the meridional velocity would require integration along the boundaries of the control volumes, which is less convenient for arbitrary unstructured meshes, and it would also neglect important information about diapycnal velocities. In contrast, the vertical velocity approach naturally yields the AMOC in density space, reduces noise introduced by the beforementioned integration along the boundaries of the control volumes, and enables more efficient online diagnostics in FESOM2 (Sidorenko et al., 2020a, 2021). Additionally, since model drift is negligible during model runtime and at the density range of the upper cell (Fig. A2), surface and interior transformations constitute the only ρ-AMOC components in our simulations.

3 AMOC in Density and Depth Space

Figure 1a and b display the meridional overturning in the Atlantic basin in density space. In PI (Fig. 1a), four cells that are associated with different water masses can be identified:

1. the mid-depth clockwise cell, located within the density classes of σ₂=36–36.9 kg m⁻³: this is the cell that we define as ρ-AMOC, that is associated with NADW formation in the subpolar North Atlantic. Its maximum overturning recirculation, located between the latitudes of 40–60° N, will be referred throughout this study as the ρ-AMOC maximum for simplicity;

2. the shallower clockwise cell above the ρ-AMOC, located within the density classes of σ₂=32–36 kg m⁻³: this cell will be referred as the subtropical gyre cell and is centered at 25° N, representing the diapycnal component of the subtropical gyre associated with Subtropical Mode Water formation (Megann, 2018; Groeskamp et al., 2019);

3. the abyssal anticlockwise cell bellow the ρ-AMOC, located within the density classes of σ₂=36.9–37.05 kg m⁻³: this cell represent the Antarctic Bottom Water (AABW) entrainment in the Atlantic Ocean, from the Southern Ocean (Orsi et al., 1999).

4. the weaker clockwise cell north of the ρ-AMOC, located within the density classes of σ₂=36.4–37.05 kg m⁻³: this cell, confined to the latitudes of 60–80° N, represent the Nordic Seas Overflow Waters (NSOW), which are formed in the Greenland, Iceland, and Norwegian (GIN) Seas and exported southward across the Greenland–Scotland Ridge through topographically constrained passages, and contribute to NADW formation via entrainment into the Atlantic Ocean (Saunders, 2001).

Figure 1 Mean AMOC strength in units of Sverdrup (1 Sv ≡10⁶ m³ s⁻¹) averaged over the last 50 years of the PI and 4xCO₂ simulations. Panels show (a, b) ρ-AMOC, (c, d) ρ-AMOC remapped into depth coordinates, and (e, f) z-AMOC for (a, c, e) PI and (b, d, f) 4xCO₂. In (c, d), ρ-AMOC is remapped into depth coordinates by loading the mean layer thickness of each density class, cumulatively summing these thicknesses to obtain the bottom depth of each class, then shifting the cumulative sum down one index (with the surface layer reset to zero) so that each transport bin appears at the depth of its upper boundary.

Conversely, Fig. 1e and f display the meridional overturning in the Atlantic basin in depth space. In comparison with Fig. 1, the mid-depth clockwise cell (z-AMOC) in Fig. 1e extends from 750–1500 m up to 55° N, with its maximum located between 30–45° N. The abyssal cell is centered near 4 km depth and reaches up to 45° N, with strengths ranging from 2–6 Sv. The subtropycal gyre cell is confined to the upper 500 m and centered around 20° N, in contrast to its more poleward location in density space. Another difference between the frameworks is that the AMOC and abyssal cell maxima appear stronger in density space, while the subtropical gyre maximum is weaker. Overall, the apparent southward displacement of all three cells in depth space relative to density space arises from the averaging of sloping isopycnals into constant-depth intervals in the subpolar North Atlantic. Specifically for the AMOC maximum, this subtropical location is also evident in the z-AMOC representation of many ocean, climate, and Earth System models (Baker et al., 2025), highlighting a systematic bias that arises from the z-AMOC diagnostics.

To facilitate direct comparison between ρ-AMOC and z-AMOC and visualizing not only the better latitudinal coherence of the ρ-AMOC maximum in PI but also of the other before-mentioned cells, we also remapped ρ-AMOC onto depth levels (Fig. 1c and d). Moreover, we provide a detailed comparison of the location and strength of the ρ- and z-AMOC maxima under PI and 4xCO₂ conditions in Table 1.

Table 1 Location and strength of the ρ- and z-AMOC maxima, computed as the average of the last 50 years of each simulation.

In 4xCO₂, the ρ-AMOC maximum weakens significantly to 7.66 Sv in comparison to its PI strength of 21.11 Sv, accompanied by a shift in density from σ₂=36.68–35.87 kg m⁻³. Such an upward shift of the AMOC maximum can also be referred to as AMOC shoaling and is a common feature in modelling studies that employ elevated radiative forcing with respect to the pre-industrial climate (Matos et al., 2020; Baker et al., 2025). The latitude of the maximum undergoes a northward shift, from 55–58° N (Fig. 1a and b), directly followed by the NSOW signature that shifts to the same density class. Both of these shifts are associated with sea ice loss and subsequent increase of deep convection in the GIN Seas that modify the density and rate of entrainment of NSOW into the North Atlantic (Chafik and Rossby, 2019). Similarly, the z-AMOC maximum strength declines from 16.32 Sv in PI to 7.61 Sv under 4xCO₂ conditions, with a corresponding shoaling from 910 to 790 m but with a southward shift in latitude from 39–35° N (Fig. 1e and f, Table 1). Overall, Fig. 1 highlights two major consequences of quadrupling atmospheric CO₂ concentrations to the AMOC in both frameworks: (1) the AMOC weakening and shoaling, and (2) the AABW weakening and upward expansion. These findings regarding the AMOC mean state in PI and under a warmer climate, both in depth and density space, corroborate those from other models with similar configurations (Xu et al., 2018; Sidorenko et al., 2020a, 2021; Megann, 2018; Megann et al., 2021; Baker et al., 2025), which increases confidence that the simulated mechanisms are robust features of the AMOC response to strong radiative forcing.

To evaluate the consistency of these phenomena across our 200 year integration period, we define two AMOC indices in both density and depth spaces, derived from the streamfunction of each model year: (1) AMOC_max, denoting the annual maximum overturning between north of 50° N, representing subpolar AMOC; and (2) AMOC₂₆, denoting the annual maximum overturning at 26° N, representing the subtropical AMOC. The indices were generated within the entire density range for ρ-AMOC, and below 500 m for z-AMOC, generating continuous 200-point time series for each index.

Figure 2 displays the time series of both indices, with annual means represented by thin lines and 15 year running means superimposed as thick lines to attenuate interannual fluctuations and accentuate multi-decadal variability. A first-degree polynomial trend was fitted and subtracted from the multi-decadal time series to facilitate correlation analysis using Pearson's correlation test (Table A2). Furthermore, the magnitude of AMOC variability was quantified by its standard deviation (σ). Please note that 50 year averages previously referenced relate exclusively to the climatological fields (mean state) depicted in Fig. 1.

Figure 2 Annual mean time-series (thin lines) of ρ-AMOC (dark blue, royal blue) and z-AMOC (green, yellowgreen) in units of Sverdrup (Sv) under PI (dark blue, green) and 4xCO₂ (royal blue, yellowgreen) forcing with superimposed 15 year rolling mean (thick lines). (a) AMOC strength at 26° N (AMOC₂₆ index), (b) maximum AMOC strength north of 50° N (AMOC_max index).

The AMOC₂₆ and AMOC_max indices (Fig. 2) emphasize that the choice of framework significantly influences the interpretation of AMOC changes to different forcings across latitudes. In PI (Fig. 2a), ρ-AMOC_max is about 13 Sv stronger than z-AMOC_max, with both displaying similar intra- to multi-decadal variability, but not correlated (Table A2). This difference in magnitude arises from the distinct latitudinal location of the AMOC maximum in density and depth space and has implications for the comparison of z-AMOC with AMOC observations in the SPNA, potentially creating model-observation discord and adding further uncertainty to assessments of the AMOC state as the climate continues to change. Under 4xCO₂ conditions, the magnitude of AMOC weakening is similar in both diagnostics until simulation year 50. Thereafter, ρ-AMOC_max and z-AMOC_max diverge in both magnitude and variability, with z-AMOC_max displaying a slight recovery toward the end of the simulation and ρ-AMOC_max oscillating around 7 Sv, suggesting a mechanism that inhibits AMOC recovery that is concealed in depth space diagnostics. Moreover, z-AMOC_max in PI and in 4xCO₂ are weaker than at 26° N, which highlights the fact that z-AMOC does not display overturning in the subpolar North Atlantic, as its maximum is located at 39° N and not at 50° N.

The AMOC₂₆ time series (Fig. 2b) further underscore the influence of framework choice to the interpretation of the AMOC strength and variability. In PI, both ρ-AMOC₂₆ and z-AMOC₂₆ exhibit interannual to multi-decadal variability, with an average strength of approximately 13 Sv. The equivalent ρ-AMOC₂₆ and z-AMOC₂₆ strengths in PI are attributable to the less steep isopycnals in subtropical regions, which are more closely aligned with z-coordinates, thus enhancing the similarity of z- and ρ-AMOC (Wett et al., 2025; Moat et al., 2025). In contrast, in 4xCO₂, ρ-AMOC₂₆ and z-AMOC₂₆ diverge substantially in magnitude and variability, suggesting that continued diagnosis of the AMOC in depth space, even in subtropical-focused studies predicated on their similarity under pre-industrial conditions, may lead to biased conclusions. In summary, Figs. 1 and 2 together demonstrate that the two frameworks, despite targeting the same circulation, can yield fundamentally different representations and interpretations of AMOC sensitivity to climate change. Ultimately, the broader question of which framework is more physically meaningful and advantageous, whether in model-based studies or in observational analyses, remains open, particularly when assessing the AMOC in climates warmer than present.

Regarding the influence of background climate on the AMOC, previous modelling studies have indicated that the extent of AMOC decline under increased radiative forcing depends on initial AMOC strength and subsequent changes in large-scale climate patterns, including surface air temperature, precipitation and sea-ice concentration (Bellomo et al., 2021; Zhao et al., 2024; Sigmond et al., 2020). These investigations indicate a tendency for models exhibiting stronger (weaker) AMOC in the control state to manifest a greater (smaller) degree of AMOC weakening. However, our modelling effort reveals that these conclusions might require revision once more research institutes adopt the ρ-AMOC diagnostics in their simulations, as our PI ρ-AMOC_max (stronger) and z-AMOC_max (weaker) reach a comparable strength in at the end of the 4xCO₂ simulation. Because these two magnitudes are derived within the same simulation, the large-scale climatic consequences of quadrupling CO₂ remain identical, even though the interpretation of the AMOC's role in driving these changes depends on the chosen diagnostics. We certainly do not aim at questioning the discussions arising from AMOC analyses in these or other studies, but simply refer to the fact that the AMOC is usually stronger in the subpolar North Atlantic in density space in comparison to depth space, therefore, the degree of ρ-AMOC weakening tends to be greater than in z-AMOC (Baker et al., 2025). Furthermore, the implementation of ρ-AMOC as a standard output may lead to convergence between models previously exhibiting diverging AMOC strengths, as the magnitude of strength and weakening may be framework-dependent. Therefore, our findings suggest that calculating AMOC solely in depth space may underestimate both the magnitude of AMOC weakening and its long-term impacts.

In addition to the relation between AMOC initial strength and subsequent weakening, Bellomo et al. (2021) demonstrated AMOC recovery occurring after surface air temperature stabilizes. However, even upon reaching quasi-equilibrium, our simulations show that this recovery only holds true for z-AMOC, with ρ-AMOC showing damped oscillations around a nearly collapsed state, which contradicts the recent hypotheses of AMOC collapse in warmer climates (Ditlevsen and Ditlevsen, 2023). These damped oscillations can be attributed to feedbacks involving interactions between AMOC and Arctic and Antarctic winter sea-ice throughout the simulation period (Zhao et al., 2024; Nobre et al., 2023). Additionally, previous studies suggest that winter sea-ice loss in the Arctic leads to increased freshwater input into the Atlantic and Pacific basins via Arctic gateways, which directly influences interior cross-isopycnal mixing through entrainment of overflows with fresher and warmer characteristics, thereby reducing or disrupting deep-water formation (Ionita et al., 2016; Mecking et al., 2017; Liu and Fedorov, 2022). An intensified hydrological cycle is also suggested as a potential driver of AMOC weakening due to freshwater forcing from runoff and the melting of land and sea-ice (Sévellec et al., 2017). In fact, the Arctic winter sea-ice area through our simulation (Fig. 3) reveals a similar multidecadal variability superimposed onto the intra-decadal variability of both ρ-AMOC indices, suggesting that Arctic winter sea ice loss and subsequent increased entrainment of lighter overflow waters constitute one of the mechanisms driving the AMOC weakening in our 4xCO₂ simulation.

Figure 3 Winter Sea ice area in the Arctic and Southern Oceans in PI and 4xCO₂ simulations. Annual means are represented in thin lines and a superimposed 15 year rolling mean is represented in thick lines to highlight multidecadal variability.

However, such a relationship is not apparent upon comparing z-AMOC and Arctic winter sea ice area curves, indicating that stable sea-ice loss resulting from a 4xCO₂ climate establishes a new regime of persistently weak AMOC and may preclude AMOC recovery, at least on interannual to centennial timescales. On the other hand, some studies indicate a z-AMOC recovery due to the salinization of the surface ocean layers driven by increased net-evaporation over the North Atlantic in warmer climates (e.g. Weijer et al., 2020; Ackermann et al., 2020), which contradicts our ρ-AMOC results and supports our z-AMOC results. However, the magnitude and timescale of this recovery vary greatly among the model simulations (Bonan et al., 2022; Zhao et al., 2024; Nobre et al., 2023), which increases the uncertainty in the AMOC state in warm climates when diagnosed in depth space. Additional climatic phenomena have been suggested as drivers or contributors to AMOC weakening, especially in the common era, such as atmospheric blocking, meso- to submesoscale dynamics, and increase in atmospheric CO₂ driving sea-ice and ice sheet melting (Ionita et al., 2016; Mecking et al., 2017; Liu and Fedorov, 2022). However, to derive such analyses, we must have model output at finer time and spatial resolutions, and transient simulations, respectively. Our simulations, however, were performed in a quasi-equilibrium formulation without dynamical ice sheets, deriving annual ρ-AMOC output and at a relatively coarse resolution. Further studies tackling the AMOC sensitivity to multiple forcings, resolutions, and components of the Earth system, while comparing both diagnostics would be beneficial to increase the confidence in our understand of the AMOC dynamics across multiple scales.

4 Large-scale and regional patterns of water mass transformations

As detailed in Sect. 2, the diapycnal velocity we diagnose inherently reflects water mass transformations. This provides a natural framework for investigating the respective contributions of surface and interior transformations to the AMOC response to abrupt climate change. These analyses are not performed in depth space here for two main reasons. First, in depth coordinates, changes in density at a given depth from both surface buoyancy forcing and interior mixing are superimposed, making it challenging to disentangle their individual contributions. Additionally, buoyancy fluxes, driven by heat and freshwater exchanges, act primarily along density surfaces rather than depth levels, therefore diagnosing them in depth-space may lead to more biases arising from averaging across constant depth intervals (Walin, 1982). In contrast, in density space, interior transformations can be obtained directly from the difference between the overturning streamfunction, model drift, and surface transformations (Eq. 3). This approach is less convoluted and allows for a clearer attribution of density changes in the mixed layer to modifications in the thermohaline properties of water masses and mixing with waters advected from, and to, other regions.

Surface transformations (ψ_S; Fig. 4a and c), are not intrinsically suited for water mass tracking since they measure only the magnitude and sign of diapycnal transformations (Xu et al., 2018). Nevertheless, they indicate upwelling or downwelling of water masses and provide attribution of their latitudinal circulation regimes based on the density range of each transformation cell. In the pre-industrial simulation (Fig. 4a), three main surface transformation cells can be identified within the latitudinal limits of the AMOC upper limb:

1. a tropical cell centered at σ₂=31.5 kg m⁻³, that enters the basin from the south and upwells at around 10° N, were the peak anticlockwise overturning is reached;

2. a subtropical cell centered at σ₂=35 kg m⁻³, that becomes progressively denser due to buoyancy loss in the subtropics, thereby reaching a maximum of 40° N; and

3. a subpolar cell, that occupies the density classes between σ₂=36.75–37 kg m⁻³, flowing northward until ∼65° N, where it becomes denser, and return southwards as the NADW.

The Antarctic Intermediate Water (AAIW) can also be identified within the subpolar cell at intermediate density (σ₂≈36.55 kg m⁻³) and is originated in the Southern Ocean at intermediate density levels due to its low temperature and salinity contrast (Santoso and England, 2004). In summary, surface transformations occur toward lighter waters in the tropics, and toward denser waters in the subtropical and subpolar North Atlantic in PI, but do not fully resembles the ρ-AMOC structure as depicted in Fig. 1a, which arises from the sum of surface and interior transformations in our simulation (as model drift is negligible). Therefore, the three cells identified in Fig. 4a are modified by interior transformations that either reinforce or counteract their contribution of surface transformations to the overturning circulation. Within the tropical cell, interior transformations act towards denser waters, while in the subtropical cell they act towards lighter waters, thereby counteracting surface transformations. This interplay results, in Fig. 1a, in the clockwise (downwelling) circulation for σ₂<32 kg m⁻³ south of the Equator and anticlockwise (upwelling) circulation before the shallower secondary maxima in the subtropical North Atlantic. At higher densities, interior transformations drive the AMOC maximum recirculation cell, indicating the entrainment of water masses and processes associated with the non-linearity of the equation of state (e.g. cabbeling and thermobaricity), that contribute to NADW formation in the North Atlantic at σ₂=36.8–36.9 kg m⁻³. The abyssal cell associated with Antarctic Bottom Water (AABW) in Fig. 1 appears in Fig. 4b as an anticlockwise cell counteracting surface transformations for σ₂=36.8–37 kg m⁻³. The signature of NSOW north of 60° N and σ₂=36.4–37.05 kg m⁻³ is also present, which indicates that its contribution to ρ-AMOC is primarily driven by interior transformations.

Figure 4 Mean surface (ψ_S; a, c) and interior (ψ_I; b, d) water mass transformation streamfunctions in (a, b) PI and (c, d) 4xCO₂, respectively. Blue contours indicate upwelling (i.e., transformations toward lighter waters), while red contours indicate downwelling (i.e., transformations toward denser waters).

In 4xCO₂ (Fig. 4c and d), a similar counterbalance between surface and interior transformations occurs in the tropics and subtropics. While their individual contribution to water mass transformations in these regions is stronger, the net result is a weaker secondary maximum than in Fig. 1a. In the subpolar cell, both surface and interior transformations occur at lighter density classes, reflecting shoaling of the ρ-AMOC. This shift is accompanied by a large-scale reduction in meridional heat transport to the subpolar North Atlantic of approximately 0.2 PW between 30 to 65° N (Fig. 5), which is consistent with the decrease of about 4 Sv in surface transformations in this region (Fig. 4c). Concurrently, interior transformations within the subpolar cell act toward lighter rather than denser waters, yielding a net ρ-AMOC weakening. This suggests that the NSOW contributes with lighter waters to NADW formation in comparison to PI, occupying the same density classes of the subpolar cell in Fig. 1b in comparison to Fig. 1a. Moreover, the abyssal cell remains identifiable as an anticlockwise deep cell in Fig. 4d albeit at much reduced magnitude and occupying density classes from σ₂=36.4–37 kg m⁻³, consistent with the weakened and upward-expanded AABW in Fig. 1b. Overall, this reversal in contribution of interior transformations to the AMOC strength, concurrently with its pronounced multidecadal variability (Fig. 6f and h) driven by modifications in the entrainment of lighter overflows attributed to extensive Arctic winter sea ice loss (Fig. 3), yields the dampened oscillation of ρ-AMOC_max and ρ-AMOC₂₆ in 4xCO₂ (Fig. 2). This suggests that the variability of interior transformations in the subpolar North Atlantic is the main driver of ρ-AMOC variability. Nevertheless, surface transformations remain the dominant driver of ρ-AMOC strength throughout the simulation (Fig. 6a, c, e, and g), predominantly forced by heat flux contributions toward denser waters in both PI and 4xCO₂ (Fig. A1). This implies that the cooling of the upper ocean within the subpolar North Atlantic triggers deep convection, which is either reinforced (in the PI case) or dampened (in the 4xCO₂ case) by interior transformations, modulated by the variability of the latter.

Figure 5 Atlantic Ocean meridional heat transport for PI and 4xCO₂.

Figure 6 Hovmöller diagrams of surface and interior transformations in density space for PI (a–d) and 4xCO₂ (e–h). First two columns (a, b, e, f) display the surface transformations, ψ_S, at (a, e) the ρ-AMOC maximum latitude and (b, f) 26° N, while last two columns (c, d, g, h) display the interior transformations, ψ_I, at (c, g) the ρ-AMOC maximum latitude and (d, h) 26° N. Positive (red) values indicate transformations toward lighter waters (buoyancy gain/upwelling), and negative (blue) values indicate transformations toward denser waters (buoyancy loss/downwelling). Red dashed lines denote the upper and lower density bounds of the ρ-AMOC upper limb, computed from the mean over the last 50 years of each simulation. The ρ-AMOC maximum latitude for both simulations is displayed in Table 1.

Figures 4 and 6 present a latitudinal density-based depiction of the ρ-AMOC components. However, acknowledging that a key factor contributing to the discrepancy between the ρ- and z-AMOC is the loss of information on the horizontal separation of the AMOC upper and lower limbs inherent in the diagnosis of the z-AMOC, the spatial distribution of vertical transport and diapycnal transformations within the North Atlantic is investigated (Fig. 7), specifically at the density and depth levels corresponding to the maximum of the AMOC upper limb (Table 1). These characteristic densities differ between experiments, as each represents the overturning maximum specific to its climate state. It should be noted that, within Fig. 7, negative values (depicted in blue) indicate downward fluxes (buoyancy loss or downwelling), while positive values (depicted in red) indicate upward fluxes (buoyancy gain or upwelling).

In PI (Fig. 7a), pronounced downwelling occurs in the Labrador (LS) and Irminger Seas (IS), and along the Norwegian coast, while comparatively weaker upwelling is confined to areas north of Iceland. The location of these transformations, occurring at the isopycnal surface of σ₂=36.68 kg m⁻³, generally align with the diapycnal velocity field (Fig. 7c), but the magnitude of upward and downward fluxes reveals the intensified contribution from interior transformations to ρ-AMOC. Additionally, in Fig. 7c, upward fluxes emerge in the Baffin Bay, along the Gulf Stream, along the path of the North Atlantic Current (NAC), and in the western GIN seas, whereas downward fluxes are stronger in the IS compared to the LS, which is consistent with observations (Lozier et al., 2019; Fu et al., 2023).

Conversely, vertical velocities at ∼910 m in PI (Fig. 7e) exhibit much noisier patterns of upwelling and downwelling compared to diapycnal velocities (Fig. 7c). A significant downward flux is evident south of Greenland, contributing to z-AMOC, while upward fluxes prevail in the LS and IS. Deep convection at these sites is still represented by a dipole pattern, indicating upwelling (downwelling) in the southern (northern) LS and IS. However, z-AMOC lacks a subpolar maximum (Fig. 1e), as opposing fluxes when integrated zonally result in weaker overturning strength. In contrast, due to the predominance of downwelling as the form of vertical motion from 30–50° N, the zonally integrated northward overturning is amplified in this region, explaining the subtropical location of the z-AMOC maximum.

Under 4xCO₂ forcing, the diapycnal velocity field (Fig. 7d) exhibits significantly reduced downward flux in the LS coupled with the emergence of upward fluxes south of Greenland. Strong downward fluxes remain in the IS, primarily driven by surface transformations (Fig. 7b). These fluxes emerge from the eastern GIN seas towards the northeastern North Atlantic, reflecting amplified GIN seas deep convection and subsequent entrainment of lighter NSOW in the North Atlantic at this isopycnal surface, consistent with the shift of the overturning maximum to a lighter isopycnal in 4xCO₂. Concurrently, upward fluxes are reduced along the path of the Gulf Stream, which is coherent with the reduced ocean meridional heat transport that is more pronounced north of 10° N (Fig. 5). In contrast, the vertical velocity field (Fig. 7f), indicates reduced downward and upward fluxes across all regions, except the eastern North Atlantic and Greenland Sea. Overall, the reduced downward flux near the latitudinal location of both the ρ- and z-AMOC maxima highlights that both diagnostics capture the reduced deep convection that drives AMOC weakening in 4xCO₂, albeit at different latitudinal bands (e.g. SPNA in density space and STNA in depth space).

Given that deep convection in the North Atlantic occur in the SPNA, we argue that only the diagnostics in density space provide physically meaningful information on the AMOC response to abrupt climate change in this region, by directly accounting for the contributions of surface and interior transformations to deep water formation. Comparing the diapycnal velocity field with observations in the SPNA reveals that this metric is able to capture the North Atlantic Ocean circulation regime in the current climate, like the dominant contribution of the IS to deep convection in the SPNA (Lozier et al., 2019; Sanchez-Franks et al., 2024) and the horizontal separation of the AMOC upper and lower limbs (Foukal and Chafik, 2024). Moreover, Årthun et al. (2023) suggest that not only does the Irminger Sea remain the dominant deep convection site in the North Atlantic, but also the GIN Seas experience amplified overturning under future climate projections, as suggested by the 4xCO₂ diapycnal velocity field (Fig. 7d). With respect to oceanic heat transport, AMOC weakening is associated with a decrease in Gulf Stream volume transport towards the subpolar North Atlantic (Piecuch and Beal, 2023), a phenomenon evident in Figs. 5 and 7d. However, the vertical velocity field (Fig. 7e) contradicts the aforementioned studies, displaying a dominant contribution from the Labrador Sea, increased upwelling in the Greenland-Iceland-Norwegian Seas, and a convoluted picture of the upper and lower AMOC limbs (Fig. 7e). Therefore, our findings corroborate the evidence from literature that ρ-AMOC can provide a substantially more accurate representation of circulation changes, and AMOC resilience and natural variability under warmer climates (e.g. Xu et al., 2018; Megann, 2018; Megann et al., 2021; Sidorenko et al., 2020a, 2021; Frajka-Williams et al., 2023; Foukal and Chafik, 2024).

Figure 7 Mean surface-forced diapycnal water mass transformations (ψ_S; a, b), diapycnal velocity (c, d) and vertical velocity (e, f) under (a, c, e) PI and (b, d, f) 4xCO₂ conditions. Depth and density levels differ: PI (910 m, 36.68 kg m⁻³) vs. 4xCO₂ (790 m, 35.87 kg m⁻³). Negative values (blue) indicate downward fluxes (buoyancy loss), whereas positive values (red) indicate upward fluxes (buoyancy gain).

5 Conclusions

We have examined the patterns and variability of the AMOC and associated water mass transformations in depth (z-AMOC) and density (ρ-AMOC) spaces using the AWI-CM3 model. Building upon previous studies that compare these frameworks in different model configurations and experiments, and through model intercomparison projects, we argue that more modelling efforts like this study are critical to evaluate the model skill in representing the basin-scale AMOC sensitivity different boundary conditions and forcings. Additionally, as surface and interior water mass transformations play a crucial role in setting AMOC strength and variability, accounting for these transformations as we do in our study is essential to advance our understanding of overturning regimes across various scales. Nevertheless, the majority of models with simulations included in the CMIP6 project at various resolutions for past, present and future climate only derive z-AMOC as output, thus a community-wide adoption of ρ-AMOC is necessary to reduce model uncertainty in the future (Foukal and Chafik, 2024). Conversely, observational arrays like the RAPID and OSNAP arrays, already provide AMOC estimates in density space, even when the difference between ρ-AMOC and z-AMOC appears negligible, like at 26.5° N (Frajka-Williams et al., 2023; Moat et al., 2025). Therefore, continuing to study the AMOC only in depth space risks increasing the model-data discord, particularly at higher latitudes. Such discrepancies raise important questions regarding the reliability of our current understanding the AMOC, and our ability to project the response of large-scale ocean circulation under the impact of climate change. To this extent, the warm climate experiment that we employed in this study provides valuable insights into AMOC sensitivity to abrupt climate change, highlighting the interplay between sea-ice formation (Fig. 3), meridional heat transport (Fig. 5), and water mass transformations (Figs. 4 and 6) in driving the AMOC dynamics (Nobre et al., 2023; Baker et al., 2025). Moreover, even though one might argue that a substantially warmer world might seem far-fetched, our findings in density space suggest that its effect on the ocean circulation is not unrealistic and quite similar to the changes in overturning circulation occurring in the past decade and projected to occur in the future (Lozier et al., 2019; Chafik et al., 2022; Sanchez-Franks et al., 2024; Årthun et al., 2023; Baker et al., 2025).

Our results demonstrate that z- and ρ-AMOC, despite targeting the same circulation, can yield fundamentally different representations and interpretations of AMOC dynamics. In our PI simulation, the first discrepancy between the frameworks concerns the latitude of the maximum overturning strength. While the ρ-AMOC maximum is located at the subpolar North Atlantic, where deep convection occurs (Lozier et al., 2019), the z-AMOC maximum is located in the subtropical North Atlantic. This reflects the loss of information on the horizontal separation of the AMOC upper and lower limbs that occur as an effect of averaging density levels into constant-depth intervals during zonal integration (Foukal and Chafik, 2024). Consequently, z-AMOC indices reflect subtropical rather than subpolar dynamics, underestimating AMOC strength by approximately 13 Sv compared to ρ-AMOC. Conversely, both diagnostics appear interchangeable at 26° N under PI conditions, implying that density space is not essential for present-climate studies utilizing subtropical AMOC strength as a basin-scale index. This equivalence, however, dissipates in warmer climates, with substantial divergence in both magnitude and variability between frameworks, not only in the subpolar North Atlantic but also at 26° N. Nevertheless, both AMOC diagnostics only agree to the AMOC strength of about 7 Sv in the final simulation years, and to the overall shoaling of the AMOC, and AABW northward-upward shift.

Furthermore, our findings suggest that the divergence between the results of our PI and 4xCO₂ simulations stems from the distinct role of interior and surface transformations in driving the AMOC, that can only be directly accounted for in density space. In 4xCO₂, overall, surface transformations drive the AMOC strength while interior transformations drive the variability of the AMOC. In this sense, deep convection is maintained in the Irminger Sea and is intensified in the GIN Seas, driven by amplified winter sea-ice melt in the Arctic. This signal is evident through the analysis of the diapycnal velocity field (Fig. 7d) and from the shoaling of the NSOW to the same density class of the mid-depth cell that represents ρ-AMOC (Fig. 1b). The deep convection in the IS and GIN seas is maintained mostly by surface transformations in 4xCO₂, contrary to PI, where interior transformations contribute at a similar magnitude. On the other hand, the effect of interior transformations in the ρ-AMOC in 4xCO₂ continues to be substantial, particularly to its variability, as the entrainment of lighter overflows occur concurrently with shifts in ρ-AMOC strength. Finally, by showing that ρ- and z-AMOC diverge not only in the subpolar North Atlantic but also in the subtropical North Atlantic under warmer climates, our work provides a novel perspective on the implications of framework choice. Additionally, we have provided substantial evidence that the role of surface and interior transformations in driving the AMOC must be explicitly assessed.

Despite the valuable insights gained through our findings we acknowledge that ρ-AMOC computations require higher computational costs and are not as readily available in large model ensembles as z-AMOC output (Sidorenko et al., 2021; Baker et al., 2025), which may delay the community-wide transition between diagnostics. However, advances towards resolving computational costs including online diagnostics and evaluation of appropriate time averaging periods may facilitate this transition (Sidorenko et al., 2020a; Megann, 2024). Indubitably, z-AMOC has provided meaningful advances in our understanding of the mechanism driving AMOC variability, and we do not advocate for its exclusion. Rather, we advocate for including ρ-AMOC as standard model output in the forthcoming CMIP phase (CMIP7), alongside with more studies that employ ρ-AMOC as the standard diagnostics. The present results, coupled with existing observational, modelling, and comprehensive studies suggest that computing AMOC in density space provides a more physically coherent depiction of the mechanisms governing AMOC internal variability and stability (e.g. Roberts et al., 2020; Sidorenko et al., 2020a, 2021; Frajka-Williams et al., 2023; Foukal and Chafik, 2024) and that both inter-model and model-observation discrepancies can be reduced through a further community-wide effort towards diagnosing the AMOC in density space. As the Earth continues to warm, diagnosing the AMOC in density space becomes even more critical given the potential divergence in basin-wide AMOC representations and their respective interpretations across different latitudes, and time- and spatial scales, as both models and observations inform political, economic, and societal decision-making under climate change scenarios.

Appendix A

Table A1 z- and ρ-AMOC equations with horizontal velocities.

Table A2 Correlation coefficients for z- and ρ-AMOC timeseries. NS = not significant.

Figure A1 Hovmöller diagrams of the contributions of heat fluxes (FH) and freshwater fluxes (FW) to surface transformations for PI (a–d) and 4xCO₂ (e–h). First two columns (a, b, e, f) display the FH at (a, e) the ρ-AMOC maximum latitude and (b, f) 26° N, while last two columns (c, d, g, h) display the FW at (c, g) the ρ-AMOC maximum latitude and (d, h) 26° N. Positive (red) values indicate transformations toward lighter waters (buoyancy gain/upwelling), and negative (blue) values indicate transformations toward denser waters (buoyancy loss/downwelling). Red dashed lines denote the upper and lower density bounds of the ρ-AMOC upper limb, computed from the mean over the last 50 years of each simulation. The ρ-AMOC maximum latitude for both simulations is displayed in Table 1.

Figure A2 Hovmöller diagrams of ρ-AMOC (a, c) and model drift (b, d) in density space for PI (a, b) and 4xCO₂ (c, d). Red dashed lines denote the upper and lower density bounds of the ρ-AMOC upper limb, computed from the mean over the last 50 years of each simulation. The ρ-AMOC maximum latitude for both simulations is displayed in Table 1.