Antarctic Bottom Water and North Atlantic Deep Water in CMIP6 models

Deep and bottom water formation are crucial components of the global ocean circulation, yet they were poorly represented in the previous generation of climate models. We here quantify biases in Antarctic Bottom Water (AABW) and North Atlantic Deep Water (NADW) formation, properties, transport, and global extent in 35 climate models that participated in the latest Climate Model Intercomparison Project (CMIP6). Several CMIP6 models are correctly forming AABW via shelf processes, but 28 models in the Southern Ocean and all 35 models in the North Atlantic form deep and bottom water via open-ocean deep convection too deeply, too often, and/or over too large an area. Models that convect the least form the most accurate AABW but the least accurate NADW. The four CESM2 models with their overflow parameterisation are among the most accurate models. In the Atlantic, the colder the AABW, the stronger the abyssal overturning at 30^∘ S, and the further north the AABW layer extends. The saltier the NADW, the stronger the Atlantic Meridional Overturning Circulation (AMOC), and the further south the NADW layer extends. In the Indian and Pacific oceans in contrast, the fresher models are the ones which extend the furthest regardless of the strength of their abyssal overturning, most likely because they are also the models with the weakest fronts in the Antarctic Circumpolar Current. There are clear improvements since CMIP5: several CMIP6 models correctly represent or parameterise Antarctic shelf processes, fewer models exhibit Southern Ocean deep convection, more models convect at the right location in the Labrador Sea, bottom density biases are reduced, and abyssal overturning is more realistic. However, more improvements are required, e.g. by generalising the use of overflow parameterisations or by coupling to interactive ice sheet models, before deep and bottom water formation, and hence heat and carbon storage, are represented accurately.

1 Introduction

Bottom water formation around Antarctica and deep water formation in the North Atlantic ventilate the global abyssal and deep ocean. Ocean–ice–atmosphere interactions by the Antarctic ice shelves (Orsi, 2010; Drucker et al., 2011; Ohshima et al., 2013) or, more rarely, in open-ocean polynyas (Killworth, 1983; Campbell et al., 2019), create the coldest and densest water mass: the Antarctic Bottom Water (AABW). AABW does not stay around Antarctica but instead travels north on the sea floor as a several hundred to few thousand metre thick layer, filling all three basins (Johnson, 2008). In a substantial portion of the Atlantic, AABW spreading north is overlain by North Atlantic Deep Water (NADW) spreading south (Johnson, 2008). NADW forms in the Labrador Sea and Nordic seas because of strong winds and haline convection, respectively (Killworth, 1983). It is the saltiest of the two water masses but is also warmer and lighter than AABW and hence leaves the sea floor to continue circulating above AABW where the two meet (Johnson, 2008). NADW production has long been linked to the strength of the Atlantic Meridional Overturning Circulation (AMOC, e.g. Broecker, 1995), although observations from the recently deployed Overturning in the Subpolar North Atlantic Program (OSNAP) line (Lozier et al., 2019) suggest that this link is more complex than simply meaning more deep-water formation equals stronger AMOC. Perhaps more crucially, AABW and NADW formation provide a direct path from the atmosphere to the bottom of the ocean and as such a conduit for heat and carbon storage (Chen et al., 2019; Zanna et al., 2019). An accurate representation of deep-water formation in climate models is thus a necessary precondition for trustworthy future climate projections.

The Climate Model Intercomparison Project phase 6 (CMIP6, Eyring et al., 2016) is the latest release of CMIP, the organised effort to make global climate models comparable, notably by running them with the same forcings. In the previous instalment, CMIP5 (Taylor et al., 2012), Southern Ocean mixed layers were poorly represented (Sallée et al., 2013) and models were forming the majority of their AABW wrongly via open-ocean convection, mostly in the overly frequent Weddell Polynya (Heuzé et al., 2013). They would only stop doing so once the ocean surface had freshened enough, by 2200 (De Lavergne et al., 2014), and consequently underestimated 21st century bottom property changes (Heuzé et al., 2015). NADW formation was more accurately represented, and although this was due to overly large sea ice extents in the North Atlantic, it occurred more in the Irminger Sea than in the Labrador Sea (Menary et al., 2015; Heuzé, 2017). So far, results on CMIP6 models have shown that sea ice representation has improved in both hemispheres, but the intermodel spread remains large (Roach et al., 2020; Shu et al., 2020). The Weddell Polynya is still opening too often but only in half of the models (Mohrmann et al., 2021), i.e. those that now have the most accurate Antarctic Circumpolar Current (Meijers et al., 2012; Beadling et al., 2020). CMIP6 models also have a higher climate sensitivity than CMIP5 models (Zelinka et al., 2020), more in line with the observed sensitivity (Armour, 2017; Cox et al., 2018). Their AMOC, however, is too sensitive to the new aerosol forcings (Menary et al., 2020). CMIP6 resolution is still coarse, with most models having a horizontal resolution of 1^∘, but recent results showed that NADW formation is in fact less accurate with higher resolution (Koenigk et al., 2020). In summary, by improving these other crucial processes but keeping a low resolution, CMIP6 models should have more realistic AABW and NADW than CMIP5 models. We here investigate whether this is the case.

In this paper, we determine the characteristics of Antarctic Bottom Water (Sect. 3.1) and North Atlantic Deep Water (Sect. 3.2) in CMIP6 models, focussing first on their respective formation processes, properties, and biases. The primary objective of this paper is to quantify and discuss biases of each model, so that model users can make informed model selections. Multi-model means are also presented at the end of each subsection. We then study the global transport of these two water masses (Sect. 3.3), and specifically how their properties determine their global extent. Finally (Sect. 4), we conclude this paper by a discussion on what – if anything – has improved since CMIP5.

Ziehn et al. (2017)Wu et al. (2019)Wu et al. (2019)Rong et al. (2019)Danabasoglu et al. (2020)Danabasoglu et al. (2020)Danabasoglu et al. (2020)Danabasoglu et al. (2020)Voldoire et al. (2019)Séférian et al. (2019)Swart et al. (2019)Held et al. (2019)Kuhlbrodt et al. (2018)Volodin and Gritsun (2018)Lurton et al. (2020)Hajima et al. (2020)Tatebe et al. (2019)Mauritsen et al. (2019)Müller et al. (2018)Mauritsen et al. (2019)Yukimoto et al. (2019)Cao et al. (2018)Counillon et al. (2016)Tjiputra et al. (2020)Tjiputra et al. (2020)Park et al. (2019)Sellar et al. (2020)

Table 1 The 35 CMIP6 models used in this study, their ocean component, nominal horizontal resolution in ^∘ latitude × ^∘ longitude, vertical grid type (ρ means isopycnic, σ terrain-following, several symbols a hybrid grid) and number of vertical levels, and official reference. N/A indicates that no paper has been published yet for the CMIP6 configuration.

2 Methods

2.1 CMIP6 models and observation-based reference data

We use the 35 CMIP6 models listed in Table 1. The only criterion for choosing them was the availability of at least their seawater salinity and temperature monthly output “so” and “thetao”, respectively, for the entire historical run (January 1850 to December 2014) at the latest date of download (20 May 2020). When available, we also directly used their monthly mixed-layer depth “mlotst”; if not, we computed it from the monthly salinity and temperature as detailed in Sect. 2.2. We also made use of each model's bathymetry “deptho” and grid cell area “areacello” files to accelerate our computations. Finally, for the transport calculations, we used the monthly meridional velocity “vo”. All output data were obtained on the model's native grid, except for NorESM2-LM and -MM that submitted their temperature and salinity on an isopycnic vertical grid; for these two models, we used the regularised z level outputs.

We used only one ensemble member per model, as even by the latest date of download the majority of models had provided only one member. Furthermore, as some models are not fully independent due to sharing similar codes (Table 1), using different ensemble sizes would have accentuated the bias towards one model family. To account for this lack of independence, the correlations quoted throughout the text have been verified with different model numbers (not shown). For most models, the ensemble member we used is referred to as r1i1p1f1. It was not available for CNRM-CM6-1, CNRM-ESM2-1, MIROC-ES2L, and UKESM1-0-LL, for which we used r1i1p1f2. Neither were available for HadGEM3-GC31-LL, for which we used r1i1p1f3.

Although we used the full historical run for robustness verifications, we present only the results for the period January 1985 to December 2014, for consistency with the observational products. Note that we neither detrended the CMIP6 historical run nor subtracted the pre-industrial control run, again for consistency with observations (which feature the climate change trend). These observations are the full-depth ocean temperature and salinity climatologies from the World Ocean Atlas 2018 (Locarnini et al., 2018; Zweng et al., 2018, respectively), the annual mixed-layer depth climatology first described by de Boyer Montégut et al. (2004), and the global bathymetry GEBCO (GEBCO Compilation Group, 2019).

2.2 Computations: deep and bottom water properties, transports, and extents

To start with, when necessary, we computed the monthly mixed-layer depth (MLD) of the CMIP6 models as per the CMIP6 procedures by first computing the monthly mean potential density σ_θ from their monthly practical salinity and potential temperature. As is requested for CMIP6, the MLD is then diagnosed as the depth where σ_θ differs from that at 10 m depth by more than 0.125 kg m⁻³. Because of averaging effects, potentially combined with the non-linearity of the equation of state, this MLD computed from monthly temperature and salinity differs slightly from mlotst, which is the monthly average of the daily or higher resolution MLD outputted by the model. As shown for CanESM5 in Fig. A1 in the Supplement, both mlotst and our recomputed MLD have the same spatial patterns, but their largest values can differ by up to 300 m. As (1) the same regions are detected as having MLD exceeding the thresholds listed below and (2) the alternative is to not use the models that do not provide mlotst, we consider the difference acceptable. Furthermore, a different threshold of 0.03 kg m⁻³ is used by the observational reference (de Boyer Montégut et al., 2004), which could lead to an apparent shallow bias for the models' mixed-layer depths (as we show in Sect. 3, it does not). We could then quantify bottom water formation in the three sectors of the Southern Ocean (south of 50^∘ S with borders at 65^∘ W and 50 and 130^∘ E, orange contours in Fig. 1), in the North Atlantic subpolar gyre (SPG, 50 to 66^∘ N, 70 to 20^∘ W) and in the Nordic seas (GIN, 65 to 80^∘ N, 30^∘ W to 20^∘ E, orange contours in Fig. 3) by computing the deep mixed volume (DMV) of each region as in Brodeau and Koenigk (2016). That is, for each month and each region, we keep only those grid cells where the MLD exceeds a critical value and sum the product MLD x cell area. We work with the maximum value of each year. As in Brodeau and Koenigk (2016) and Koenigk et al. (2020), we use a critical value of 700 m in the Nordic seas as it is the depth of the sill that connects them to the rest of the world ocean, and 1000 m in the Labrador Sea. As in, e.g. Heuzé et al. (2013); De Lavergne et al. (2014), we use 2000 m in all three Southern Ocean sectors.

Figure 1 Southern Ocean reference bottom density σ_θ (top left panel, top colour bar), and for each CMIP6 model, bottom density bias (model minus reference) averaged over 1985–2014. The white number for each model is its RMSE over the entire Southern Ocean deeper than 1000 m. The thick black line indicates maximum mixed layer deeper than 2000 m. The thin grey line indicates the 2000 m isobath. Orange lines on the reference panel delineate the Weddell (W), Amery (A), and Ross (R) sectors for the DMV calculation (see Sect. 2.2.).

Table 2 Median and maximum deep mixing volume (DMV; see Sect. 2) for the Southern Ocean sectors (orange contours in Fig. 1) for each CMIP6 model over 1985–2014. Values given in 10¹³ m³, which is approximately the DMV of a 1^∘$\times \mathrm{1}{}^{\circ }$ grid cell with a 1000 m mixed layer. Number in brackets indicates how many years out of 30 that the DMV is different from zero, i.e. the number of years with deep convection.

Figure 2 Multi-model summary of the biases for the Southern Ocean. Columns (from left to right) show observational reference and multi-model mean property (same colour bar for both), mean bias (model minus reference), and standard deviation of the difference. Rows (from top to bottom) show bottom density, bottom potential temperature, bottom salinity, and mixed-layer depth. White numbers are the median value over the deep Southern Ocean, defined as in Fig. 1. Grey lines indicate the 1000 and 2000 m isobaths.

Figure 3 North Atlantic reference bottom density σ_θ (top left panel, top colour bar). For each CMIP6 model, bottom density bias (model minus reference) averaged over 1985–2014 is given. Orange lines on the reference panel delineate the subpolar gyre (SPG) and Nordic seas (GIN) areas for RMSE and DMV calculation. The yellow dotted line shows the SPG sector of Johnson (2008). White numbers for each model show their RMSE over the GIN (top) and SPG (bottom) areas for depths over 1000 m. The thick black line indicates the maximum mixed layer deeper than 1000 m. The dotted cyan line indicates the same in GIN deeper than 700 m. The thin grey line shows the 1000 m isobath.

Table 3 Median and maximum deep mixing volume (DMV; see Sect. 2) for the subpolar gyre (SPG) and Nordic seas (GIN; see Fig. 3) for each CMIP6 model over 1985–2014. Values given in 10¹³ m³, which is approximately the DMV of a 1^∘$\times \mathrm{1}{}^{\circ }$ grid cell with a 1000 m mixed layer. The number in brackets indicates how many years out of 30 that the DMV is different from zero, i.e. the number of years with deep convection.

Figure 4 Multi-model summary of the biases for the North Atlantic. Columns (from left to right) show the observational reference and multi-model mean property (same colour bar for both), mean bias (model minus reference), and standard deviation of the difference. Rows (from top to bottom) show the bottom density, bottom potential temperature, bottom salinity, and mixed-layer depth. The white numbers are the median value over the deep GIN (top) and SPG (bottom), defined as in Fig. 3. Grey lines indicate the 1000 and 2000 m isobaths.

We quantify biases in the models by computing the root-mean-square error (model minus reference) in temperature, salinity, and density σ_θ at the sea floor grid cell. To do so, all models had to be interpolated onto the reference's grid. After interpolation we also computed the multi-model mean properties, mean bias, and standard deviation of the bias for each grid cell. Note that we purposely keep σ_θ instead of σ₄, as σ_θ is the density used in the models' code to notably compute the MLD. For later calculations, we also compute the temperature and salinity of the water masses AABW and NADW by taking their average properties over a specific region. As we will show in Sect. 3.1, the AABW formation region really differs from model to model; as such, instead of using a limited region as in Johnson (2008), we detect AABW as having the temperature minimum anywhere deeper than 2000 m and south of 50^∘ S. For NADW, we produce two flavours: NADW_SPG as having the salinity maximum anywhere deeper than 1000 m in the small area of SPG defined by Johnson (2008, 53 to 63^∘ N, 55 to 54^∘ W, yellow box in Fig. 3), and NADW_GIN, the salinity maximum anywhere deeper than 1000 m in the GIN sector defined above. Note that the definition of NADW_GIN is different from that of Johnson (2008) because of extra tests not shown here and known past model biases in the vertical structure of the North Atlantic (e.g. Menary and Wood, 2018) and in the North Atlantic–Nordic Seas exchanges (e.g. Heuzé and Årthun, 2019).

We not only investigate the properties of AABW and NADW by their formation region but also their transport into the rest of the global ocean. For AABW, we hence compute each model's Southern Meridional Overturning Circulation or SMOC, using the same method as Heuzé et al. (2015) for comparison purposes. That is, in each ocean basin, we first integrate the meridional velocity vo from the west coast to the east coast at 30^∘ S, then we integrate this value from the sea floor to the surface. The SMOC then is the northward deep maximum. We use a similar method for the AMOC, computing it at 35^∘ N instead for comparison with the CMIP6 results of Menary et al. (2020). After integration from coast to coast and from sea floor to surface of the velocity vo, the AMOC is defined as the southward subsurface maximum. We could not directly use the meridional overturning circulation output “msftmz” for the following reasons:

• it is provided by only 18 of the models (from 10 families);

• it is in kg s⁻¹ instead of m³ s⁻¹, requiring division by the density, for which we only have the monthly mean;

• for most models, the Indian and Pacific oceans are provided as one joint region, so we could not have obtained the SMOC in each basin.

Having to interpolate the irregular model grids onto the sections instead of directly using the model output may have introduced some errors, but as the AMOC results of this paper and those of Menary et al. (2020) for the models and experiment we have in common are similar, we are confident in our MOC values. Note that two models, GFDL-ESM4 and NorCPM1, did not provide vo, limiting our transport analysis to 33 CMIP6 models.

Finally, to investigate in CMIP6 the link found in Heuzé et al. (2015) between the SMOCs and the northward extent of AABW layer, we chose to re-create the Johnson (2008) maps of AABW and NADW volumes in the global ocean for CMIP6 models. However, using the same approach as Johnson (2008), whereby we would have to determine the characteristics of every water mass in each basin for each model, is beyond the scope of this paper. Instead, as below the core of NADW the global ocean (excluding the Arctic) is a mixture of NADW and AABW only, we determine the NADW and AABW contents at each depth from a conservative property χ using the mixture equation of Jenkins (1999):

$$\begin{array}{}\text{(1)}& {\displaystyle }{\displaystyle }{\text{NADW}}_{\text{content}}={\displaystyle \frac{{\mathit{\chi }}_{\text{AABW}}-\mathit{\chi }}{{\mathit{\chi }}_{\text{AABW}}-{\mathit{\chi }}_{\text{NADW}}}},\end{array}$$

and

$$\begin{array}{}\text{(2)}& {\displaystyle }{\displaystyle }{\text{AABW}}_{\text{content}}={\displaystyle \frac{{\mathit{\chi }}_{\text{NADW}}-\mathit{\chi }}{{\mathit{\chi }}_{\text{NADW}}-{\mathit{\chi }}_{\text{AABW}}}}.\end{array}$$

Here, as in Johnson (2008), we consider the practical salinity and potential temperature as conservative enough to be used for these calculations. We then take the 50 % content depth as the border between the NADW and AABW layers, i.e. anything with more than 50 % AABW or less than 50 % NADW is in the AABW layer, and the AABW thickness is the difference between the depth of that border and the sea floor. We finally take the median of all the combinations: temperature or salinity, NADW properties from SPG or GIN, and AABW or NADW contents. For the NADW layer, we detect the NADW core as the maximum NADW content with the extra criterion that the maximum must be larger than 80 % NADW. Tests with values ranging from 60 % to 100 % yield similar values (not shown). Then the so-called NADW thickness is the thickness from the depth of the core to the NADW–AABW boundary (or to the sea floor if there is no AABW). By working with a mixture of two water masses only, we could not try and detect the top of the NADW layer. Note that traditional methods of using a fixed temperature and/or salinity for water mass determination cannot be applied to potentially biased climate models. The northward extent of AABW in each basin is defined as the northernmost latitude of the uninterrupted contour of thickness =2000 m that starts in the Southern Ocean. We do the same for the southward extent of NADW in the Atlantic Ocean. For this part of the analysis, we show only the results for NADW that originated in SPG; NADW that originated in GIN seems to leave the Nordic seas in no model (not shown).

3 Results

In this section, we first look at bottom water formation and properties in the Southern Ocean and then deep water formation and properties in the North Atlantic. It is only in the last section that we analyse both water masses together by determining their global transports and volumes. In this section, we only talk about CMIP6. The comparison with CMIP5 will come in Sect. 4.

3.1 Southern Ocean bottom water characteristics in CMIP6 models

Shelf overflow and open-ocean deep convection in the Southern Ocean

Presently in the real ocean, Antarctic Bottom Water is primarily formed in several locations (including the Weddell Sea, the Ross Sea, and by Adélie Land) as water is cooled, made saltier, and becomes denser on the continental shelves and then cascades down the continental slopes, entraining deep waters on its way to the sea floor (visible as the densest areas in Fig. 1). The CMIP6 models' bottom density bias on the shelves suggests that 19 out of 35 models may form dense water on the shelf: ACCESS-CM2 (Weddell and Ross), ACCESS-ESM1-5 (Ross), CAMS-CSM1-0 (Ross), the four CESM2 (Ross), CanESM5 (Weddell and Ross), GFDL-ESM4 (Weddell and Ross), the three GISS (Ross mainly), HadGEM3-GC31-LL (Weddell and Ross), INM-CM5 (Weddell and Ross), IPSL-CM6A-LR (Weddell and Ross), the two NorESM2 (Ross mainly), SAM0-UNICON (Ross mainly), and UKESM1-0-LL (Weddell and Ross). The other 16 models are too light (strong negative bias in Fig. 1). Mean biases over 30 years are not enough to determine whether the dense shelf water flows into the deep ocean; we instead created movies of the monthly bottom density over the entire historical run for these 19 models, of which 2 are available as a video supplement.

The movies let us split these 19 models into the following three groups.

• INM-CM5 and NorESM2-LM show overflowing from the Ross shelf to the deep basin (as does NorESM2-MM) and in the Amery sector (see video).

• GFDL-ESM4, HadGEM3-GC31-LL, IPSL-CM6A-LR, SAM0-UNICON, and UKESM1-0-LL may overflow in the Ross sector, but we would need higher temporal resolution data to be certain.

• The other 11 models occasionally have a plume of dense water leaving the shelf, but it is nowhere near as dense as the shelf water it originates from (see video example of ACCESS-CM2).

In summary, in no model is there any (obvious) shelf export in the Weddell Sea. INM-CM5 (terrain following, high horizontal resolution model) and the two NorESM2 (hybrid isopycnic models) are the only ones forming AABW accurately via shelf processes, in the Ross sector only. The other high resolution models are not dense enough on the shelf or not clearly exporting their shelf water. There is no clear link between shelf processes and resolution (horizontal or vertical), vertical grid type, or ocean model component.

If up to 8 models may export dense water from the shelf, how do the other 27 models of our study form their AABW? Deep convection, so far observed only once in the real ocean, in response to the 1974–1976 Weddell Polynya (Killworth, 1983), is the next obvious process to investigate. Mixed-layer depths exceeding 2000 m are most prevalent in the Weddell sector (black contours in Fig. 1, DMVs in Table 2, MLD in Table A1). Of our 35 models, 24 exhibit deep convection in the Weddell Sea, of which 19 do so for most years of our study period. Most of these models also have a too large and too frequent Weddell Polynya (Mohrmann et al., 2021), except for GFDL-CM4 and IPSL-CM6A-LR, which may be convecting under sea ice cover, and the two MIROC, which are ice-free (Mohrmann et al., 2021; Roach et al., 2020). In the Amery and Ross sectors, we need to distinguish between the models that have non-zero DMV because of open-ocean deep convection and those with coastal polynyas. In the Amery sector, aside from MCM-UA-1-0, whose deep convection area is just a continuation of the Weddell one, we have 10 out of 35 models with open-ocean deep convection: the two CNRM, the two EC-Earth3, GISS-E2-1-H, MIROC6, MPI-ESM-1-2-HAM, MPI-ESM1-2-LR, and the two NorESM2. In the Ross sector, we have 15 models with open-ocean deep convection: 14 in the Ross Sea and 1 (NESM3) in the Amundsen Sea. In the Amery and Ross sectors, there is no link anymore between DMVs and the polynya activity in Mohrmann et al. (2021), suggesting that bottom water formation occurs under ice cover. There is, however, a strong correlation of +0.47 between DMV in the Ross sector and DMV in the Weddell sector, i.e., models that convect a lot do it in both sectors. Behrens et al. (2016) suggests that a strong DMV is associated with a strong Antarctic Circumpolar Current (ACC), while Cabré et al. (2017) find that strong DMV would weaken the westerly winds, i.e. may weaken the ACC. Here the only relationship between DMV and ACC, using the values from Beadling et al. (2020), is in agreement with Cabré et al. (2017): the more deep convection in the Amery sector, the weaker the ACC (correlation of −0.54, significant at 95 % level). We find no relationship with DMV in the Weddell or Ross sectors.

Up to now, there are still seven models that have no open-ocean deep convection and no shelf overflow: the four CESM2, GISS-E2-1-G and GISS-E2-G-CC, and NorCPM1. GISS-E2-1-G and GISS-E2-G-CC have non-zero DMV when considering the entire historical run (1850–2014, not shown), with GISS-E2-1-G convecting once in the Weddell sector and once in the Ross sector, and GISS-E2-1-G-CC thrice in the Weddell sector. The four CESM2 models do not, but they have an overflow parameterisation that artificially moves dense water from the shelf to the deep basins (Briegleb et al., 2010). If the water on the shelf exceeds a critical density, a pipe artificially transports this dense water from the shelves to the deep basin. Without having to cascade, the dense shelf water keeps its properties. Consequently, the absence of cascade means that we cannot detect it on the overflow movies. For NorCPM1 though, we could not determine how its AABW is formed; maybe it formed before 1850. In conclusion, most models form their AABW via open-ocean deep convection. Even the models that seem to represent shelf processes accurately exhibit open-ocean deep convection. Somewhat surprisingly, the only relationship between the DMV and the climate sensitivities of Zelinka et al. (2020) is in the Weddell sector (correlation of −0.36): models that convect a lot there have a low sensitivity, which is to be expected as heat and CO₂ are sent to the deep ocean. What is surprising is that the relationship holds only in the Weddell sector. Hence, the sensitivity might be more linked to polynya activity, which is linked to deep convection only in the Weddell sector.

AABW properties

Does the way CMIP6 models form their AABW impact its characteristics, as it did in CMIP5 (Heuzé et al., 2013)? Only density biases are shown in Fig. 1, but salinity and temperature biases are provided in Figs. A2 and A3, respectively. A total of 10 models have a negligible bottom density bias (root-mean-square error, RMSE, lower than 0.05 kg m⁻³, white numbers in Fig. 1): UKESM1-0-LL, CanESM5, IPSL-CM6A-LR, CESM2-WACCM-FV2, CESM2-FV2, CESM2, CESM2-WACCM, GISS-E2-1-H, CNRM-CM6, and HadGEM3-GC31-LL. A total of 12 more models have an acceptable bottom density bias (RMSE lower than 0.1 kg m⁻³), including all the other models that are based on NEMO. That means that 22 of 35 models have acceptable biases. Let us investigate rather what may be common to the 13 models that are performing poorly.

INM-CM5 is the only model that is biased dense. Its bottom salinity is extremely high (RMSE of 0.42, Fig. A2) while its bottom temperature is rather accurate (RMSE of 0.8 ^∘C, Fig. A3). Its predecessor INM-CM4 had a similar issue, although whether this was caused by a too short spin-up or non-conservation of salt could not be determined (Alexander Gusev, personal communication, July 2014). All the other 12 models are biased light (Fig. 1). For BCC-CSM2-MR, BCC-ESM1, MCM-UA-1-0, MIROC-ES2L, MIROC6 and MRI-ESM2-0, it is because of a fresh bias (Fig. A2). The other 6 models have relatively accurate bottom salinity, but are biased warm (Fig. A3). CMIP6 models are overall biased light or biased dense in the entire deep Southern Ocean (excluding the shelves).

Finally, the multi-model mean bottom water is too light (mean bias of 0.06 kg m⁻³), too warm (0.85 ^∘C), and too fresh (0.02 psu) throughout the entire Southern Ocean (Fig. 2). Technically, as the multi-model mean bottom temperature is 0.57 ^∘C, we can even say that on average CMIP6 models do not have AABW (defined as having a temperature below 0 ^∘C) at the bottom of the Southern Ocean. The multi-model mean biases in mixed-layer depth exceed 2000 m in all sectors, especially in the Weddell Sea, but they do so over a small region so that the average bias over the entire deep Southern Ocean is low (151 m). As mentioned before, the biases in properties are moderate in the open ocean but exceed 0.1 kg m⁻³ or 0.1 psu on the shelves: the multi-model mean shelves are too light and too fresh as many models do not form high-salinity shelf water. This is also the reason why the standard deviation (Fig. 2 rightmost column) is largest on the shelves. In the open ocean, the standard deviation is larger in the Weddell Sea than in other sectors for the temperature, salinity and MLD, as it is where DMVs differ most (Table 2).

The models with low biases in bottom density also tend to have zero to low DMVs in the Weddell Sea, but the relationship does not hold for maximum DMVs larger than 200 10¹³ m³ (Table 2). NorESM2-LM and NorESM2-MM notably have low biases but very high DMVs, but they also do shelf overflows. Open-ocean deep convection leads to a warming and salinification of bottom waters (Zanowski et al., 2015); one hypothesis is then that models that hardly convect stay closer to the bottom density value they were initialised with. In the case of the CESM2 suite, the overflow parameterisation may help form accurate bottom water. Biases are not the whole story though. As expected, we do find a significant relationship (95 % level) between the actual temperature and salinity of AABW and the DMV: in the Amery and Ross sectors, more deep convection leads to warmer AABW (correlation of +0.33 and +0.29, respectively) as in Wang et al. (2017). In the Weddell sector, however, more deep convection leads to fresher AABW (correlation of −0.35), which in fact is consistent with the short-term response of the Southern Ocean to deep convection in Zanowski et al. (2015). The multi-model mean AABW salinity is 34.606±0.154; the reference value from Johnson (2008) is 34.641. The multi-model mean AABW temperature is $-\mathrm{0.45}\pm \mathrm{0.73}$ ^∘C; the reference value from Johnson (2008) is −0.88 ^∘C. Therefore, the multi-model mean AABW is warmer and fresher than the reference, and more DMV worsens these biases. Note that the values of the individual models are given in Table A2.

To summarise, in the Southern Ocean most models form their AABW by open-ocean deep convection. In the Weddell Sea, this convection seems tied to the Weddell Polynya activity and impacts the AABW salinity most: more deep convection means fresher bottom salinity. In the Amery and Ross sectors, it is linked more to the bottom temperature: more deep convection means warmer bottom salinity. Models that seem to form dense water via shelf processes also exhibit deep convection, so we cannot determine whether overflows alone would make the Southern Ocean more accurate. The models that convect the least or not at all tend to be the most accurate. For the CESM2 family, accurate bottom properties and lack of deep convection may both be the result of their overflow parameterisation (Briegleb et al., 2010; Snow et al., 2015). For another model, NorCPM1, the accuracy in all properties may come from its observation assimilation rather than accurate model physics (Counillon et al., 2016).

We will study the impact of these biases on the global transport of AABW in Sect. 3.3, but as we cannot do so without investigating the AABW–NADW balance in the Atlantic basin, let us first evaluate the representation of NADW in CMIP6 models.

3.2 North Atlantic deep water in CMIP6 models

Deep water formation in the North Atlantic subpolar gyre and Nordic seas

In the North Atlantic, all 35 CMIP6 models of our study exhibit deep convection in the subpolar gyre (black contours in Fig. 3 and Table 3). As in CMIP5 (Heuzé, 2017), a large proportion of them convect not only in the Labrador Sea as the reference but also intensely south of Iceland (Irminger Sea):

• 6 of the 35 models convect only in the Labrador Sea (Fig. 3): CNRM-CM6-1, CNRM-ESM2-1, EC-Earth3-Veg, IPSL-CM6A-LR, NorCPM1, and NorESM2-LM;

• 9 of the 35 models convect both in the Labrador and Irminger seas, but the two regions are not connected: ACCESS-CM2, BCC-ESM1, CESM2-WACCM-FV2, EC-Earth3, HadGEM3-GC31-LL, INM-CM5, MPI-ESM1-2-LR, NorESM2-MM, and UKESM1-0-LL;

• 17 of the 35 models convect both in the Labrador and Irminger seas as one large SPG deep convection area: ACCESS-ESM1-5, BCC-CSM2-MR, CESM2, CESM2-FV2, CESM2-WACCM, GFDL-CM4, GFDL-ESM4, GISS-E2-1-G, GISS-E2-1-G-CC, GISS-E2-1-H, MCM-UA-1-0, MIROC6, MPI-ESM-1-2-HAM, MPI-ESM1-2-HR, MRI-ESM2-0, NESM3, and SAM0-UNICON;

• the final 3 of the 35 models convect only in the Irminger Sea: CAMS-CSM1-0, CanESM5, and MIROC-ES2L.

As in Koenigk et al. (2020), the higher-resolution NorESM2-MM and MPI-ESM1-2-HR have larger deep convection areas than their corresponding lower-resolution NorESM2-LM and MPI-ESM1-2-LR. Note that the difference between NorESM2-MM and NorESM2-LM is in the atmospheric component resolution. There is, however, no robust relationship across models between the horizontal resolution and the DMV. There is a relationship with the climate sensitivities of Zelinka et al. (2020) though: the larger the DMV in SPG, the less sensitive the model (correlation of −0.36), which as already discussed in the Southern Ocean part is not surprising. Consequently, there is also a strong correlation between the DMV in the Weddell Sea and in SPG (+0.57): models that convect a lot in the Weddell Sea convect a lot in SPG as well. As already mentioned, no causation can be inferred: deciphering whether global biases in DMV are responsible for the models' sensitivities or sensitivities are set by other processes and impact the DMV is beyond the scope of this analysis.

All models except INM-CM5-0 and NorCPM have deep convection in the GIN seas as well. Moreover, in GIN models convect most years, with a minimum as high as 24 of 30 years (Table 3). There is more variability in the SPG, but likewise the majority of models convect all years. Besides, they convect too deep. While in the Southern Ocean, deep convection to the sea floor can happen (Killworth, 1983), in the North Atlantic it should not go much beyond 1000 m (e.g. Våge et al., 2009). In the SPG, only the three models CanESM5, INM-CM5, and NorCPM1 have maximum mixed-layer depths just exceeding 1000 m (Table A1). An extra four models, ACCESS-CM2, CESM2, CESM2-WACCM, and MIROC-ES2L, have tolerable depths up to 2500 m. However, the vast majority convects too deep, often to the sea floor. It is the same in GIN, albeit with different models: this time, the four models CESM2, MPI-ESM1-2-LR, and NESM3 have MLDs up to 2000 m, and all the other models go to 3000 m or even the sea floor. There is no significant correlation between the DMV in SPG and that in GIN. In summary, CMIP6 models exhibit deep convection in the North Atlantic too often, too deep, and over too large an area. It is not possible to determine the single most accurate model in the North Atlantic or even in each subregion; model users must choose a compromise between correct representation of the variability, location, depth, or extent.

North Atlantic bottom properties

CMIP6 water property biases at the bottom of the North Atlantic are of the same order of magnitude as those at the bottom of the Southern Ocean (shading in Fig. 3). Three models have bottom density biases resulting in an RMSE lower than 0.05 kg m⁻³ in both SPG and GIN: CESM2-WACCM, HadGEM3-GC31-LL, and UKESM1-0-LL. An additional nine models have an RMSE lower than 0.1 kg m⁻³ in both SPG and GIN: CESM2, CESM2-FV2, CNRM-CM6-1, EC-Earth3-Veg, GFDL-CM4, GFDL-ESM4, IPSL-CM6A-LR, MRI-ESM2-0, and SAM0-UNICON. As for the other 23 models, it depends on the region.

• INM-CM5, NorCPM1, NorESM2-M, and NorESM2-MM are biased dense in both regions because they are biased salty (Fig. A4). The magnitude of the bottom cell biases in INM-CM5 is very grid cell dependent, maybe because of faulty regularisation of the sigma grid (even though it did not have this problem in the Southern Ocean).

• ACCESS-ESM1-5 is accurate in the SPG sector but biased dense (salty) in GIN.

• CESM2-WACCM-FV2, GISS-E2-1-G, and GISS-E2-1-G-CC are accurate in the SPG sector but biased light in GIN. For CESM2-WACCM-FV2, it is because of a warm bias (Fig. A5); for GISS-E2-1-G and GISS-E2-1-G-CC, it is because of a fresh bias.

• ACCESS-CM2, BCC-CSM2-MR, CNRM-ESM2-1, EC-Earth3, MPI-ESM-1-2-HAM, and MPI-ESM1-2-LR are biased light in SPG but accurate in GIN. EC-Earth3 is the only one that is biased fresh. The other models are biased salty but warm.

• The last nine models are biased light in both regions. For CanESM5, GISS-E2-1-H, MCM-UA-1-0, MIROC-ES2L and NESM3, this is caused by a salty bias; for BCC-ESM1, CAMS-CSM1-0, MIROC6, and MPI-ESM1-2-HR, this is a warm bias.

All models except CanESM5 and IPSL-CM6A-LR are in fact biased warm compared to the World Ocean Atlas 2018 bottom temperature. The evolution of the bottom properties throughout the entire historical run is complex, with significantly different variabilities depending on the model (not shown), and their analysis is beyond the scope of this paper. All that we can say that the warm bias is not a result of only the modelled climate change or any drift.

The multi-model mean biases reflect the individual model biases that we just described (Fig. 4): in the deep ocean, in both GIN and SPG, the multi-model mean is biased salty (0.07 in GIN, 0.15 in SPG) but warm (1.68 and 1.91 ^∘C respectively), and hence slightly light (0.04 and 0.07 kg m⁻³). The across-model spread, given by the standard deviation (right column in Fig. 4), is on the order of 0.1 kg m⁻³ or psu and 1 ^∘C in both regions, similar to what we found for the Southern Ocean. The multi-model mean MLD is too deep over too large an area compared to the reference (bottom row in Fig. 4), corresponding to a bias of nearly 800 m on average in both regions. The standard deviation is relatively high though, also reaching 800 m on average, indicative of model disagreement over the location of the deep MLD. The largest bias in standard deviation, i.e. the location where the across-model agreement is lowest, is in the Greenland Sea, probably a consequence of the across-model differences in sea ice representation (Shu et al., 2020).

In the SPG, there is a strong relationship between the bottom temperature RMSE and the climate sensitivity of Zelinka et al. (2020): the more sensitive the model, the less biased the temperature in SPG (correlation of −0.68). There is a somewhat significant (at the 90 % level) relationship between the DMV and the bottom density bias in SPG only: the more the model convects, the less biased it is (correlation of −0.29). There is, however, no relationship with the location of deep convection itself, e.g. MIROC-ES2L that convects only in the Irminger Sea has a similar bias (magnitude and sign) as MPI-ESM1-2-LR that convects only in the Labrador Sea; MPI-ESM1-LR in turns has a large bias and convects at the same location as UKESM1-0-LL, which has a low bias. The four CESM2 models and their overflow parameterisation by the Denmark Strait are again among the most accurate, which was in fact the original motivation for that parameterisation (Briegleb et al., 2010). NorCPM1 is somewhat disappointing; it is built on NorESM2-LM and is supposed to have improved performance thanks to data assimilation (Counillon et al., 2016). Its bottom density is indeed better than NorESM2-LM's in SPG but is biased even denser (saltier) in GIN. There is no across-model relationship between the sensitivity or the DMV and the salinity or GIN though, so the cause of NorCPM's and other models' bottom density bias in the GIN seas remains unknown. It may be linked to their respective biases in the representation of the Atlantic Multidecadal Oscillation, as suggested by Lin et al. (2019). There is, however, a strong interhemispheric correlation: the more biased the bottom density in the Southern Ocean, the more biased it is in SPG (+0.81).

Finally, looking at the NADW properties instead of the biases (individual values in Table A2), we find that the multi-model mean NADW_SPG, which Johnson (2008) refers to as upper NADW (UNADW) or Labrador Sea Water (LSW) is too warm and too salty: 4.86±0.81 ^∘C and 35.163±0.143 in CMIP6, instead of 3.32 ^∘C and 34.894 in Johnson (2008). The multi-model mean NADW_GIN in contrast is accurate: 0.77±0.99 ^∘C and 35.001±0.169 compared to 1.30 ^∘C and 34.878 for the water mass called lower NADW (LNADW) or Iceland–Scotland Overflow Water (ISOW) in Johnson (2008), despite five models having a temperature below 0 ^∘C. Again in SPG, models with a high sensitivity have a lower salinity and temperature (correlations of −0.31 and −0.45 respectively), which is consistent with the links previously found with the DMV. But again no relationship can be found in GIN. In GIN, we found no relationship between the biases or properties and the horizontal or vertical resolution or between the grid type and ocean model component. In CMIP5, Heuzé and Årthun (2019) had found strong across-model biases in the inflow to the Nordic seas caused by the large-scale oceanic and atmospheric circulations, as well as the bathymetry, while Lin et al. (2019) showed that GIN property biases can be linked to the representation of multidecadal variability. Investigating the exact cause of the biases in GIN is beyond the scope of this paper, not least because in the next section we will show that NADW_GIN does not contribute to the global NADW in CMIP6 models. For now, we can conclude that the bottom property biases in GIN are not related to deep water formation in the region.

3.3 Global transport of NADW and AABW in CMIP6 models

In this last section, we shall determine the global fate of NADW and AABW once they leave their source regions. For NADW, this fate is tied to the strength of the Atlantic Meridional Overturning Circulation. The mean AMOC value lies at 18 Sv ($\mathrm{1}\mathrm{Sv}={\mathrm{10}}^{\mathrm{6}}{\mathrm{m}}^{\mathrm{3}}{\mathrm{s}}^{-\mathrm{1}}$), although observations both at the RAPID/MOCHA-array at 26.5^∘ N (e.g. Duchez et al., 2016) and in the more recently deployed OSNAP-lines in the subpolar North Atlantic (Lozier et al., 2019) reveal a strong interannual variability of up to 5 Sv. Aside from INM-CM5 and its AMOC of 63±19 Sv, all models fall in that range (individual values in Table B1), resulting in a multi-model median of 19.5±9.5 Sv. Observations of the southern MOC at 30^∘ S are rarer. From box inverse modelling, Lumpkin and Speer (2007) estimated the Atlantic SMOC at 5.6±3 Sv; apart from the three GISS models and MIROC-ES2L that are too weak, all models are in that range (Table B2), giving a multi-model median of $\mathrm{2.8}\pm \mathrm{1.4}.$ Observational values in the Indian Ocean range between 3 and 27 Sv (Huussen et al., 2012); thus unsurprisingly, from the weak MCM-UA-1-0 (1.5±1.6 Sv) to the strong GFDL-CM4 (11±18 Sv), all models are in that range and the multi-model median is 3.0±2.5 Sv. This is a remarkable improvement from CMIP5, where a majority of models had an Indian SMOC close to 0 (Heuzé et al., 2015). In the Pacific finally, Lumpkin and Speer (2007) estimated the MOC to be 11±5 Sv. MCM-UA-1-0 is again the weakest (3.9±1.9 Sv), and the only model that falls out of the observational range, resulting in a multi-model median of 5.9±3.0 Sv. In summary, the AMOC and southern MOCs are rather accurately represented in CMIP6 models.

The across-model correlations among the transports are strong and significant (95 % level): the stronger the SMOC in the Indian Ocean, the stronger it is in the Pacific Ocean as well (correlation of +0.37). In contrast, a strong SMOC in either of these basins corresponds to a weak SMOC in the Atlantic (Atlantic–Indian, correlation of −0.45; Atlantic–Pacific, correlation of −0.34). A weak SMOC in the Atlantic corresponds to a strong AMOC (correlation of −0.30), as previously found by Patara and Böning (2014) in the NEMO model. We are not implying causation from the correlations, but it is interesting to find relationships between the biases quantified in Sects. 3.1 and 3.2 and the transports. In agreement with Patara and Böning (2014), a stronger Atlantic SMOC is associated with lower temperature biases (correlation of 0.29), i.e. colder AABW (−0.35), whereas a stronger Pacific SMOC is associated with stronger density biases (+0.36). A stronger AMOC is associated with larger biases in temperature and salinity in SPG (correlations of +0.33 and +0.37 respectively) and in particular a saltier NADW_SPG (+0.34, as in the paleoclimate simulations of Menviel et al., 2020). The Atlantic SMOC is the only transport that is linked to the strength of the Antarctic Circumpolar Current (ACC, values from Beadling et al., 2020): the stronger the ACC, the stronger the Atlantic SMOC (correlation of +0.37). There is no significant direct relationship between the transports and the DMV, which at least for the AMOC is in agreement with the recent observations of Lozier et al. (2019) and modelling work of Årthun et al. (2019). Unlike, e.g. Menary et al. (2015) or Koenigk et al. (2020), we find no link between the MOCs and the models' horizontal resolution.

In line with Heuzé et al. (2015), we expect the transports to impact the interbasin spread of NADW and AABW, i.e. that the stronger the transport, the further from its source the water mass will travel. To investigate this, we recreated the AABW and NADW thickness maps of Johnson (2008) as Figs. 5 and 6, respectively. For Fig. 6, we show only NADW_SPG; surprisingly, in no model could we find NADW_GIN beyond the Nordic seas (not shown). In agreement with observations and Johnson (2008), AABW occupies the majority of the water column in most of the Indian and Pacific oceans, but its northward extent is limited in the Atlantic Ocean. Said extent is highly model dependent in the Atlantic, whereas it extends as far north as the basin limits permit in most models in the Indian and Pacific oceans. Finally, in most models, AABW in the Indian Ocean seems to come from the Pacific. The NADW southward expansion in the Atlantic is also model dependent, with some reaching to the Antarctic Circumpolar Current (e.g. BCC-ESM1) and others not even leaving the North Atlantic subpolar gyre (e.g. UKESM1-0-LL). As explained in the methods, the NADW layer in the Indian and Pacific oceans is most likely biased by our calculation method that takes into account only two water masses and thus shall not be discussed further.

After extracting the southernmost extent of NADW and northernmost extents of AABW for each model (see Tables B1 and B2), we do find, as expected, that the stronger the AMOC, the further south NADW extends in the Atlantic (correlation of 0.32). The stronger the Atlantic SMOC, the further north AABW extends in the Atlantic (correlation of 0.40). As we previously found an anticorrelation between the AMOC and the Atlantic SMOC across CMIP6 models, the Atlantic balance is complete: models with strong AMOC and weak SMOC have their Atlantic dominated by NADW (e.g. CESM2-WACCM), whereas those with a weak AMOC and strong SMOC are filled with AABW (e.g. IPSL-CM6A-LR). Although there was no significant correlation between the DMV and the transports, we do find that the larger the DMV, the further the extent of NADW (DMV SPG, correlation of +0.34) or AABW (DMV Weddell, +0.51). We found no significant correlation between the northward extent in the Indian or Pacific oceans and either the SMOCs or DMVs, or with the strength of the ACC. There are, however, relationships with their bottom properties: the northward extent of salty models is less than that of fresh models (correlations of −0.31 in the Indian and −0.44 in the Pacific). As we also find a strong positive relationship (correlation of +0.72) between the salinity of AABW and the salinity gradient across the ACC computed by Beadling et al. (2020), i.e., we find that the fresh models have a weak gradient to overcome, this result is not surprising. We can even speculate that in the absence of NADW, AABW would expand further north in the fresher models regardless of their SMOC.

In conclusion, in CMIP6 models as in the real ocean, deep convection impacts bottom water characteristics and biases: in the Southern Ocean, deep convection seems associated with more biased bottom waters; in the North Atlantic, the more the models convect, the less biased they are. Either way, these biases then impact the deep and bottom water transport: a saltier NADW is associated with a stronger AMOC, and a colder AABW is associated with a stronger Atlantic SMOC. These transports then impact the location of the “NADW–AABW border” in the Atlantic: stronger AMOC and weaker Atlantic SMOC (the two transports are anticorrelated), further southward extent of NADW, and less northward extent of AABW. In the Indian and Pacific oceans, the northward extent is larger in the fresher models, which are the ones with weak fronts in the ACC. To summarise, deep and bottom water formation are crucial for an accurate representation of the global deep ocean. We conclude this paper with a discussion of changes in deep and bottom water modelling since CMIP5 and what we can expect from the next generation(s) of simulations.

Figure 5 Thickness of the Antarctic Bottom Water layer in observations (top left panel) and in each CMIP6 model. See Sect. 2 for more information.

Figure 6 Thickness of the North Atlantic Deep Water layer in observations (top left) and in each CMIP6 model from the NADW core to its bottom. See Sect. 2 for more information.

4 Discussion: changes since CMIP5 and way towards CMIP7

In CMIP5 models, no model assessed by Heuzé et al. (2013) could represent dense shelf overflows correctly. Consequently, models relied on open-ocean deep convection for their bottom water formation. The right amount of deep convection in the Weddell Sea was required for accurate bottom properties; models that convected too little or too much were the most biased. This relationship does not hold for CMIP6 anymore, and it is the models that convect the least that tend to be the most accurate (Fig. 1 and Table 2). It may be because many models are now artificially prevented from opening polynyas and convecting in the Weddell Sea (Mohrmann et al., 2021). However, as the Weddell Polynya has now reopened in the real ocean (Campbell et al., 2019), future models may remove their “polynya-prevention” schemes again. Another reason for CMIP6 models seemingly not needing Southern Ocean deep convection to have accurate bottom properties may be that, as we showed in this paper, several CMIP6 models successfully represent shelf processes. This was an unexpected result considering that horizontal resolutions have not increased much since CMIP5, suggesting that models have improved their parameterisations instead (Danek et al., 2019). Regardless of the formation process, bottom density biases are smaller in CMIP6 than they were in CMIP5 (RMSEs in Fig. 2 of Heuzé et al., 2013). The new version of the models that performed well in CMIP5 also performs well in CMIP6 (e.g. the IPSL and NorESM families), and the others have improved (CanESM4 had a bias of 0.17 kg m⁻³; CanESM5 had a bias of 0.03 kg m⁻³). The worst performing model of CMIP5 was INMCM4. The worst performing model of CMIP6 with respect to Southern Ocean bottom properties is its successor, INM-CM5-0, but even this model saw its bias halve. INM-CM5-0 has both shelf processes and open-ocean deep convection, whereas INMCM4 had neither, which probably contributed to ridding the model of its cold bottom bias (Zanowski et al., 2015).

In the North Atlantic, to the best of our knowledge, most CMIP5 studies focussed on the relationship between deep water formation and the AMOC or the warming hole (e.g. Menary and Wood, 2018) but did not investigate bottom property biases. The one exception is Ba et al. (2014), who found a recurrent cold bias, with the World Ocean Atlas 2018 as reference, we find in contrast that most CMIP6 models have a warm bias at the bottom of the North Atlantic. Deep water formation in the North Atlantic in the majority of CMIP5 models occurred too often, too deep, and over too large an area (Heuzé, 2017). These findings are still valid for CMIP6 (Fig. 3 and Table 3). One noticeable improvement is that the models whose CMIP5 predecessor convected only in the Irminger Sea now convect in the entire subpolar gyre, including the Labrador Sea. Unfortunately, some of the models that performed well in CMIP5 when considering the location of deep convection in the SPG, i.e. had a relatively small area in the Labrador Sea, have also expanded to the entire SPG (e.g. the CNRM family). Therefore, the inaccurate models may be on the way to improvement. This is most likely because the Arctic sea ice is better represented in CMIP6 than in CMIP5 (Shu et al., 2020), but the ones that were relatively accurate have degraded. The same holds for the Nordic seas: CMIP6 models are convecting even more than CMIP5 models did, and they were already convecting too much. In an increasingly warmer and ice-free climate, Lique and Thomas (2018) predict that deep water formation would migrate from the North Atlantic subpolar gyre to its subtropical gyre and from the Nordic seas to the Arctic. Liu et al. (2019) adds that this will depend on whether meltwaters will most strongly impact the stratification, shutting down deep convection, or the horizontal gradients and hence the winds, pushing meltwater away from convection areas. For now, we observe that from the very icy CMIP5 to the more accurately de-iced CMIP6 models, deep water formation regions just expanded to occupy most of the space available in SPG and GIN. It is unclear whether increasing the resolution of future models would solve this issue: Danek et al. (2019) dramatically reduced mixed-layer depths in SPG by using an adaptive mesh with 5–15 km resolution, while Koenigk et al. (2020) finds that DMVs in the SPG become even larger in the high resolution versions of the models that participated in HighResMIP. Without changing the horizontal resolution, a more systematic inclusion and better representation of the stratosphere may be enough to reduce deep convection in the North Atlantic (Haase et al., 2018).

Regarding the transports, as noted by Menary et al. (2020) the AMOC is stronger in CMIP6 than in CMIP5, which they attribute to the aerosol forcing. Except for INM-CM5, which is now way too strong or which uploaded incorrect velocity fields, this increase is not that strong and most models are in the observational range. In the case of the CNRM family, a stronger AMOC is in fact a much more accurate AMOC (from 12 Sv in CMIP5 to 19 Sv in CMIP6). The NorESM models have a weaker AMOC in CMIP6, which is more accurate than their CMIP5 version (from 32 Sv in CMIP5 to 21 Sv in CMIP6). The two highest-resolution models have weakened so much that their AMOC is too low (GFDL-CM4 and MPI-ESM1-2-HR). This seems in contradiction to Koenigk et al. (2020), who found that increased resolution in HighResMIP leads to a stronger AMOC, but their result is mostly true when the models reach an eddy-resolving resolution, which they do not here in CMIP6. It is harder to determine whether the southern MOCs at 30^∘ S have improved since the values from inverse modelling (Lumpkin and Speer, 2007) and observations (Huussen et al., 2012) have very large uncertainties. All that we can say is that the Atlantic SMOC is stronger in CMIP6, and thus only the GISS family continues to have an Atlantic SMOC around 0 Sv. In the Indian Ocean, no model has a transport of 0 anymore, which resulted in a doubling of the multi-model mean from 1.6 Sv in CMIP5 to 3 Sv in CMIP6, giving it the same importance as the Atlantic SMOC. The Pacific SMOC remains the strongest of the three and sees no significant difference between CMIP5 and CMIP6 except for the two models that used to be around 0, INMCM4 (INM-CM5-0 is now at 10 Sv) and GISS-E2-H (GISS-E2-1-H now at 7 Sv). As in CMIP6, the Southern Ocean representation from the sea floor to the surface (Beadling et al., 2020) has improved, as well as the ACC (see also Beadling et al., 2020); it is no surprise that more models are now capable of exporting AABW to the rest of the world ocean. To the best of our knowledge, the global extent of AABW and NADW, presented here for CMIP6 in Figs. 5 and 6, respectively, was not assessed in CMIP5, so we cannot determine whether improved Southern Ocean characteristics lead to an improved global water mass distribution.

What can we expect from a hypothetical CMIP7? Higher resolution, most likely, although that was already expected from CMIP6 and did not happen. As explained above and by Koenigk et al. (2020) or Danek et al. (2019), a higher resolution would not necessarily improve deep water formation. Holt et al. (2017) goes as far as stating that shelf processes will not be correctly represented until the horizontal resolution remains lower than 1∕72^∘, which they expect might be reachable by the most advanced computers within 10 years. Unfortunately, we do not all have access to these computers, so that even now computing the global monthly mixed-layer depth of the highest-resolution model (GFDL-CM4, 1∕4^∘) required over 600 core hours for the 165 years of the historical run. Higher-resolution output will be impossible to manage, unless cloud-computing solutions such as PANGEO become the norm (Odaka et al., 2020). Instead of increasing the resolution, a seemingly easier solution would be to improve parameterisations (Holt et al., 2017), especially overflow parameterisations (Snow et al., 2015). Briegleb et al. (2010) first showed that an overflow parameterisation to transport water from the Nordic seas to the rest of the North Atlantic resulted in an improved representation of the ocean there. In CMIP6, the CESM2 models with their “pipes” in the North Atlantic and Antarctic shelves were among the most accurate models, especially for AABW. It would be interesting to see whether such a parameterisation on a different model would yield the same results, or whether the CESM2 models are just very accurate. Efforts could also concentrate on improving other components of the climate model, for example the atmosphere, as an improved representation of the stratosphere would supposedly decrease unrealistic deep and bottom water formation (Haase et al., 2018). However, where most progress can probably be made is in the cryosphere. As deep water formation is tied to the sea ice behaviour in both hemispheres, efforts such as sea ice MIP (SIMIP, Notz et al., 2016) dedicated to the modelling and coupling of sea ice may be the way forward. Likewise, the results of ice sheet MIP (ISMIP6, Nowicki et al., 2016) may shed a light on the debated impact of glacial meltwater on deep and bottom water formation (De Lavergne et al., 2014; Liu et al., 2019).

5 Conclusions

In this paper, we determined the characteristics of Antarctic Bottom Water and North Atlantic Deep Water in 35 models that participated in the latest instalment of the Climate Model Intercomparison Project, CMIP6: their formation, properties, transports, and extent in the global ocean. We focussed on the last 30 years of the historical run, January 1985 to December 2014. In the Southern Ocean (Sect. 3.1), bottom water formation is now more accurate, with several models representing shelf processes. Open-ocean deep convection in the Weddell Polynya still happens in more than half of the models, but it is not a requirement for accurate bottom water properties. In fact, the most accurate models were the ones with little to no open-ocean convection, especially the CESM2 family that has an overflow parameterisation. In the North Atlantic (Sect. 3.2), models convect too often, too deep, and over too large an area, but in the subpolar gyre that area has migrated from the Irminger Sea (in CMIP5 models) to the more accurate Labrador Sea. The models that convect the most in the North Atlantic subpolar gyre also have the least biased NADW. NADW formed in the subpolar gyre of the models clearly spreads southward, but the signature of the portions formed in the Nordic seas is less evident. The saltier the NADW, the stronger the AMOC and the further south the extent of NADW (Sect. 3.3). That extent is limited by the strength of the abyssal overturning in the southern Atlantic or SMOC, with stronger Atlantic SMOC (caused by colder AABW) resulting in a further northward extent of AABW. In the Indian and Pacific oceans, the extent is directly related to the AABW properties, not the SMOCs: models with a comparatively fresh AABW are also the ones with weak fronts across the Antarctic Circumpolar Current, and hence those that can travel the furthest north. In summary, for the deep and bottom water masses in CMIP6, their formation impacts their properties, which impact their transport and global extent, which in turn will have large impacts on global predictions of thermal expansion and sea level rise (Zickfeld et al., 2017), carbon storage (Tatebe et al., 2019), ecosystem changes (Sweetman et al., 2017), etc. Although CMIP6 models represent AABW and NADW more accurately than CMIP5 models did, a lot still need to be improved, especially deep and bottom water formation (Sect. 4).

How to improve deep water formation in climate models then? A higher horizontal resolution may not be the answer as, depending on the model, it either reduces (Danek et al., 2019) or increases deep convection even further (Koenigk et al., 2020). In the ocean component, one solution could be a more systematic inclusion of overflow parameterisation (Snow et al., 2015); in this study, it seems very effective for CESM2. The one data-assimilating model, NorCPM (Counillon et al., 2016), also proposes an interesting option. In the rest of the model, improving the representation of the stratosphere seems effective at reducing open-ocean deep convection (Haase et al., 2018). Whatever the future holds, we hope it will feature a more systematic archiving of useful parameters. The situation has improved since CMIP5, but there are still CMIP6 models that do not provide their monthly mixed-layer depth, and overturning streamfunctions (especially in density space) are a rarity. Making output directly available on cloud-computing-based systems such as PANGEO (Odaka et al., 2020) should also be a priority to let researchers work on heavy CMIP data as soon as they are released, regardless of their computing and storage capacities.

Appendix A: Bottom properties

In this appendix the following information is presented.

• Fig. A1 is a comparison between “mlotst” and the mixed-layer depth computed from “thetao” and “so” for one model.

• Figs. A2 and A3 show the Southern Ocean bottom salinity and temperature, respectively, to complement the bias discussion of Sect. 3.1.

• Figs. A4 and A5 show the North Atlantic bottom salinity and temperature, respectively, to complement the bias discussion of Sect. 3.2.

• Table A1 presents the maximum mixed-layer depth and convective area for each model and each region, corresponding to the DMVs discussed in Sects. 3.1 and 3.2.

• Table A2 presents the salinity and temperature of NADW and AABW for each model, briefly discussed in Sects. 3.1 and 3.2 and used for the thickness computations for Figs. 5 and 6.

Figure A1 Maximum monthly mixed-layer depth in the North Atlantic over 1985–2014 for the model CanESM5: (a) using the model output “mlotst” and (b) when computed from the monthly temperature and salinity. Over the entire 30 year period, the root-mean-square error in the Nordic seas is 305 m, and in the subpolar gyre it is 21 m.

Figure A2 Southern Ocean reference bottom practical salinity (top left panel, top colour bar). For each CMIP6 model, bottom practical salinity bias (model minus reference) averaged over 1985–2014 is shown. The white number for each model is its RMSE over the entire Southern Ocean deeper than 1000 m. The thick black line indicates maximum mixed layer deeper than 2000 m. The thin grey line shows the 2000 m isobath.

Figure A3 Southern Ocean reference bottom potential temperature (top left panel, top colour bar). For each CMIP6 model, bottom potential temperature bias (model minus reference) averaged over 1985–2014 is shown. The white number for each model is its RMSE over the entire Southern Ocean deeper than 1000 m. The thick black line indicates maximum mixed-layer deeper than 2000 m. The thin grey line shows the 2000 m isobath.

Figure A4 North Atlantic reference bottom practical salinity (top left panel, top colour bar). For each CMIP6 model, bottom practical salinity bias (model minus reference) averaged over 1985–2014 is shown. The white numbers for each model is its RMSE over the GIN (top) and SPG (bottom) areas for depths over 1000 m. The thick black line indicates maximum mixed-layer deeper than 1000 m. The dotted cyan line indicates the same in GIN deeper than 700 m. The thin grey line shows the 1000 m isobath.

Figure A5 North Atlantic reference bottom potential temperature (top left panel, top colour bar). For each CMIP6 model, bottom potential temperature bias (model minus reference) averaged over 1985–2014 is shown. The white numbers for each model is its RMSE over the GIN (top) and SPG (bottom) areas for depths over 1000 m. Thick black line indicates maximum mixed-layer deeper than 1000 m. The dotted cyan line indicates the same in GIN deeper than 700 m. The thin grey line shows the 1000 m isobath.

Table A1 Supplementary version of Tables 2 and 3 showing the 30 year max MLD (m) and max area (in 10 000 km², which is the approximate area of a 1^∘ cell).

Table A2 For each CMIP6 model, 30 year median practical salinity (S) and potential temperature (θ, ^∘C) of the NADW formed in the subpolar gyre (SPG) or Nordic seas (GIN) and of the AABW are given.

Appendix B: Transports

In this section, you will find two tables to complement Sect. 3.3.

• Table B1 presents the AMOC and southernmost extent of NADW in each model.

• Table B2 presents the SMOC and northernmost extent of AABW in the Atlantic, Indian, and Pacific oceans in each model.

Table B1 For each CMIP6 model, the 30 year median AMOC at 35^∘ N (in Sv) and the southernmost latitude (in degrees north) of the 2000 m thick NADW layer in the Atlantic from Fig. 6 are given.

Table B2 For each CMIP6 model, the 30 year median southern MOC at 30^∘ S (SMOC, in Sv) and the northernmost latitude (in degrees north) of the 2000 m thick AABW layer in each ocean from Fig. 5 are given.