2.1 General configuration changes

In NorESM2, the ocean–ice model adopts a tripolar grid instead of the bipolar grid of NorESM1. The tripolar grid was chosen because it is more isotropic at northern high latitudes and is therefore more efficient with respect to time integration (longer time step) (Guo et al., 2019). Compared to NorESM1, NorESM2 has enabled a higher-frequency ocean coupling to resolve the diurnal cycle in the flux and state exchange between the ocean and atmosphere (i.e. hourly instead of daily atmosphere–ocean coupling). We also applied a full leapfrog time stepping in both the physical and ocean biogeochemical models to improve the conservation of biogeochemical tracers (Schwinger et al., 2016).

There are several configuration options for NorESM2. For CMIP6, two configurations will mostly be used, differing only in the horizontal resolution of the atmospheric and land models with otherwise identical model components. They are named NorESM2-LM and NorESM2-MM and have an atmosphere–land resolution of approximately 2^∘ and 1^∘, respectively. In addition, both versions have a different parameter tuning in the atmospheric component, necessary to achieve a top-of-atmosphere radiative balance under pre-industrial conditions. Both NorESM2-LM and NorESM2-MM share the same ocean physical ocean model (BLOM) and the same biogeochemical component (iHAMOCC) with a horizontal resolution of approximately 1^∘. The ocean model adopts 53 isopycnal layers with two non-isopycnic surface layers, which represent (1) the top 10 m of depth and (2) the rest of the surface mixed layer. A more detailed description of the ocean physical component is provided in Bentsen et al. (2020). The ocean biogeochemistry performs similarly in both model versions. Unless explicitly specified, our analysis refers to the results from NorESM2-MM (hereafter referred to as NorESM2).

2.2 Physical parameterizations

The simulated Atlantic Meridional Overturning Circulation (AMOC) strength in NorESM1 is on the strong side (30.8 Sv at 26.5^∘ N; 1 Sv =10⁶ m³ s⁻¹) compared to observation-based estimates and other CMIP5 models (Bentsen et al., 2013). The vigorous AMOC leads to an Atlantic Ocean that is too warm and saline at depth. Lack of upper-ocean stratification at high latitudes is thought to contribute to this AMOC bias. Reformulating the eddy-induced transport (Gent and McWilliams, 1990) in order to make it more efficient in the upper-ocean non-isopycnic regime of the model and adjustments to the parameterized re-stratification by mesoscale eddies (Fox-Kemper et al., 2008) lead to improved high-latitude stratification and a more realistic mixed layer depth (MLD), particularly at high latitudes.

NorESM1 has a cold bias in the depth range 200–1000 m between 50^∘ S and 50^∘ N (Bentsen et al., 2013) that could be related to insufficient upper-ocean vertical mixing. Further, the new hourly atmosphere–ocean coupling and associated tuning of MLD parameterizations lead to a cold bias in the Pacific equatorial thermocline that also indicates too little vertical mixing.

To improve this, two different approaches have been used. First, we allow surface turbulent kinetic energy (TKE) originating from wind stirring to penetrate below the mixed layer and be added to the TKE reservoir of the k-ϵ turbulence model handling the shear instability mixing of the isopycnic interior. Secondly, with hourly atmosphere–ocean coupling, high-frequency winds are available to parameterize wind work on near-inertial motions. Inspired by the approach of Jochum et al. (2013), a fraction of this energy source is added to the mixed layer TKE reservoir, modifying the MLD estimation, while a fraction of the remaining energy is used to drive upper-ocean mixing through assumed excitation and breaking of internal waves. Combining these approaches, we were able to reduce the above-mentioned upper-ocean biases significantly. The mid-depth warm bias seen in NorESM1 is also reduced. More details on the physical improvements are described in Bentsen et al. (2020). As a result of these improvements, NorESM2 simulates a more reasonable AMOC strength of approximately 21 Sv, which is within the range of the observational estimates of 17.9±3.3 Sv (Srokosz et al., 2012).

2.3 Dimethyl sulfide

NorESM2 prognostically simulates the DMS concentration according to the formulation of Kloster et al. (2006), as follows:

where the first term on the right-hand side denotes the physical transport, the second term denotes DMS production, and the third, fourth, and fifth terms denote sinks due to bacterial activity, photolysis, and flux to the atmosphere. The production term is computed as a function of temperature, pH, and simulated detritus export production (opal and CaCO₃, implicitly computed as a function of silicate concentration):

The last term in Eq. (2) describes a pH dependence of DMS production following Six et al. (2013), which is parameterized based on the deviation of actual pH from the climatological monthly mean surface pH (pH_pi) calculated from the pre-industrial control. Due to uncertainties in the pH control on DMS production, we have omitted the pH dependency term by setting α to zero in the presented simulations. Parameters γ_o and γ_c represent scaling factors for sulfur (S:Si and S:C) and are set to 0.025 and 0.125, respectively, such that DMS production is more sensitive to CaCO₃ production in the model. This assumption is due to the fact that haptophyte-class phytoplankton generally has a higher dimethylsulfoniopropionate-to-cell-carbon ratio than diatoms (Keller et al., 1989). Since the opal and CaCO₃ portions of export production are implicitly computed as a function of silicate concentration, e.g. opal (CaCO₃) production is high (low) when surface silicate is abundant and vice versa, the DMS production is also sensitive to the long-term trend in the silicate concentration.

The DMS sink due to bacteria consumption is estimated as a function of temperature (^∘C), DMS concentration, consumption rate (γ_bac=0.0864 d^−1 ∘C⁻¹), and a half-saturation constant (k_cb=10 nmol L⁻¹), as follows:

The DMS photolysis term is defined as a function of photolysis rate (γ_uv=0.0011 m² (W d)⁻¹), incoming shortwave radiation attenuated by water and phytoplankton, I_z, and DMS concentration:

Lastly, in NorESM2, the atmospheric model receives the DMS emissions (S_gas) that are prognostically computed by the ocean biogeochemistry model (see Sect. 2.8).

2.4 Riverine input

The influx of carbon and nutrients from over 6000 rivers to the coastal oceans has been implemented based on previous work of Bernard et al. (2011) with modifications. Bernard et al. (2011) have implemented the riverine fluxes based on levels in the year 1995 provided by an early version of the Global-NEWS model (Seitzinger et al., 2005). In addition to the riverine DIC, dissolved organic carbon (DOC), inorganic nitrogen and phosphorus, particulate organic carbon (POC), particulate nitrogen and phosphorus, and dissolved silicate, we have also included riverine alkalinity (ALK) and iron (Fe). Except for DIC, ALK, and Fe, all data are provided by the more recent Global-NEWS2 model (Mayorga et al., 2010), which is a hybrid of empirical, statistical, and mechanistic model components that simulate steady-state annual riverine fluxes of carbon and nutrients. The DIC and ALK fluxes are taken from the work by Hartmann (2009). The riverine Fe flux is calculated as a proportion of global gross dissolved Fe input of 1.45 Tg yr⁻¹ (Chester, 1990), weighted by the water runoff of each river. Only 1 % of the riverine dissolved Fe is added to the oceanic dissolved Fe, since approximately 80 % to 99 % of the fluvial gross dissolved Fe is removed during estuarine mixing (e.g. Sholkovitz et al., 1981; Shiller and Boyle, 1991; Bruland and Lohan, 2014). Instead of releasing the riverine fluxes to the nearest ocean grid cell (Bernard et al., 2011), we have interpolated all fluxes at river mouths in the same way as the freshwater runoff (distributed as a function of river mouth distance with an e-folding length scale of 1000 km and cutoff of 300 km) to the ocean grid. The fluxes are specified to be constant over time at contemporary levels (the year 2000). In our model, the organic form of dissolved carbon (DOC) is connected to nutrients through the Redfield ratio ($\mathrm{C}:\mathrm{N}:\mathrm{P}=\mathrm{122}:\mathrm{16}:\mathrm{1}$), and therefore other forms of dissolved organic matter (e.g. DON and DOP) are not explicitly simulated. Since Global-NEWS provides estimates of dissolved organic matter in carbon, nitrogen, and phosphate forms (DOC_riv, DON_riv, and DOP_riv, respectively) separately, only the minimum of the three riverine dissolved organic constituents is added to the DOC term in the model (i.e. DOC = DOC + min(DOC_riv, $r_{\mathrm{C}:\mathrm{N}}\cdot {\mathrm{DON}}_{\mathrm{riv}},r_{\mathrm{C}:\mathrm{P}}\cdot {\mathrm{DOP}}_{\mathrm{riv}})$). Any excess or remaining organic matter of the three constituents is then added to the corresponding inorganic pools (DIC, NO₃, or PO₄). The same concept also applies to riverine inputs of particulate organic carbon (POC) (see also Bernard et al., 2011).

2.5 Atmospheric nitrogen deposition

We apply anthropogenic nitrogen deposition fluxes that potentially affect the simulated ocean carbon uptake (through an increase in the available nitrate for biological production) to the ocean biogeochemical model. The monthly input fields, spanning the years 1850–2014, are simulated by chemistry transport models and provided by CMIP6 (Jones et al., 2016). Conservative remapping is used to interpolate the input data from a regular $\mathrm{2.5}{}^{\circ }\times \mathrm{1.9}{}^{\circ }$ latitude–longitude grid to the tripolar ocean model grid of NorESM2. Four species of wet or dry and oxidized or reduced nitrogen deposition rates are included in the input fields. All of them are added to the nitrate pool in the topmost ocean layer.

2.6 Particle export

In NorESM1, the export of particulate organic matter is parameterized with a constant sinking speed and a constant remineralization rate at depth when sufficient oxygen is available. This simplistic formulation has been shown to have difficulty in accurately representing the observed particle fluxes at depths (Kriest and Oschlies, 2008) and leads to biases in the interior distributions of biogeochemical tracers (e.g. as seen from remineralized phosphate concentrations). In NorESM2, we have developed two additional particle flux parameterizations that can be selected. They are referred to as the WLIN and AGG schemes. The WLIN scheme is similar to the standard scheme but the sinking speed is a linear function of depth: $w_{\mathrm{POC}}=min(w_{min}+az$, w_max). Here, we use $w_{min}=\mathrm{1}$ m d⁻¹, $w_{max}=\mathrm{60}$ m d⁻¹, and $a=\mathrm{60}/\mathrm{2400}$ m d⁻¹ m⁻¹. When w_min and w_max are set to zero and infinity, respectively, this scheme is equivalent to the widely used Martin curve formulation (Martin et al., 1987).

The AGG scheme implements a prognostic sinking speed calculated according to a size distribution of sinking particles. The total concentration of particulate matter, formulated as a function of phytoplankton and detritus, forms sinking aggregates with an explicitly computed size distribution that follows a power-law formula (Kriest and Evans, 1999; Kriest, 2002):

Here, d represents the particle diameter, and A and β represent parameters of the size distribution. The minimum particle diameter l is set to 0.002 cm. The sinking speed of aggregates is computed according to the diameter of each aggregate. More details of the implementation of the AGG formulation are described in Schwinger et al. (2016). Based on the performance of the different schemes, we have set WLIN to be the default scheme and it is used in all experiments presented here. More details on the differences in the simulated interior biogeochemistry using the different particle export schemes are available in Schwinger et al. (2016).

2.7 Nitrogen fixation

Nitrogen fixation by cyanobacteria is computed implicitly and directly converted into nitrate concentration. In NorESM1, nitrogen fixation is implemented such that it can occur anytime and anywhere in the surface ocean as long as the nitrate concentration is lower than phosphate (multiplied by the stoichiometric nitrate-to-phosphate ratio R_N:P following Maier-Reimer et al., 2005). However, Breitbarth et al. (2007) provided some evidence that the large-scale distribution of Trichodesmium, a type of cyanobacteria, is well constrained by seawater temperature between 20 and 30 ^∘C. Based on this, a temperature-dependent function f(T) according to Kriest and Oschlies (2015) has been added to the nitrogen fixation formulation in NorESM2. The rate of change of nitrate owing to nitrogen fixation is formulated as follows:

with

where μ is the maximum N₂ fixation rate (0.005 d⁻¹). In addition, for each mole of N fixed to nitrate, 1.25 mol of dissolved oxygen and 1 mol of alkalinity are consumed and lost in the surface layer. The corrected formula essentially limits the occurrence of N fixation to warm low latitudes.

2.8 Air–sea gas exchange

In NorESM2, the air–sea gas exchange of CO₂, O₂, N₂, N₂O, DMS, CFC-11, CFC-12, and SF₆ is computed prognostically according to the updated formulation of Wanninkhof (2014). Here, the air–sea flux F of gas X is computed as a function of surface wind speed U, Schmidt number Sc, gas solubility K_o, and the partial pressure difference of gas X based on a bulk formulation:

The net fluxes are computed for the top layer of the ocean model (i.e. 10 m of depth). The updated formulation includes a refitted temperature-dependent Schmidt number that can be applied for a temperature range of −2 to 40 ^∘C. The solubility of all gases is computed as a function of surface temperature and salinity following Weiss (1970) for O₂ and N₂, Weiss and Price (1980) for CO₂ and N₂O, Warner and Weiss (1985) for CFC-11 and CFC-12, and Bullister et al. (2002) for SF₆. The flux of DMS is only unidirectional (outgassing; pDMS_air=0).

In the historical simulation, we use the annual atmospheric concentrations of CFC-11, CFC-12, and SF₆ divided into northern and southern hemispheric values according to Bullister (2014). For the atmospheric CO₂ concentration, monthly global mean values from the CMIP6 dataset are used. In the pre-industrial simulation, a constant CO₂ concentration of 284.32 ppm is used.

2.9 Dissolved iron parameterization

Adjustments in the iron parameterization have been made to ensure that NorESM2 simulates the observed iron-limited primary production in the Southern Ocean, equatorial Pacific, and North Pacific. In the model, the consumption and release of dissolved iron associated with biological activities are determined using a fixed stoichiometry ratio ($R_{\mathrm{Fe}:\mathrm{P}}=\mathrm{6.1}\times {\mathrm{10}}^{-\mathrm{4}}$ mol Fe mol P⁻¹, previously $\mathrm{3.7}\times {\mathrm{10}}^{-\mathrm{4}}$ mol Fe mol P⁻¹). In the surface layer, a constant fraction ($\mathrm{3.13}\times {\mathrm{10}}^{-\mathrm{6}}$ kmol Fe kg⁻¹ d⁻¹, previously $\mathrm{6.27}\times {\mathrm{10}}^{-\mathrm{6}}$ kmol Fe kg⁻¹ d⁻¹) of aerial dust deposition (Mahowald et al., 2005) is instantaneously converted into bioavailable iron. The updated global input of bioavailable iron from the atmosphere is 2.8 Gmol Fe yr⁻¹, which is within the range of values used by other models (Tagliabue et al., 2016). Finally, throughout the water column, the complexation of iron by organic substances (ligands) is assumed:

where the strength of this complexation λ_Fe is set to 0.05∕365 (previously 0.005∕365) d⁻¹ and the threshold Fe_o has been reduced from 0.6 to 0.5 nmol Fe m⁻³. The latter is motivated by the fact that the observed iron concentration in the deep Southern Ocean is lower than 0.6 nmol Fe m⁻³ (Boyd and Ellwood, 2010). The higher complexation rate leads to a faster relaxation toward this value.

2.10 Ecosystem parameterization updates

The underlying marine ecosystem parameterization in NorESM2 remains the same as in NorESM1, but many of the parameters have been adjusted to reduce biases that are present in NorESM1. The main deficiencies of the spatial annual mean primary productivity pattern in NorESM1 are primary production (PP) that is too strong in the Southern Ocean, the eastern tropical Pacific, and to a lesser degree the North Atlantic, contrasted by very low PP in the subtropical gyres (most pronounced in the Pacific). The high productivity in the high latitudes is accompanied by an exaggerated annual cycle of PP, showing a spring bloom that is too strong and a decline of PP that is too fast afterwards. Tuning of the ecosystem parameterization during the development of NorESM2 was focused on improving these regional shortcomings. The changes in ecosystem parameters listed in Table 1 mainly serve two purposes: (i) to increase top-down limitation by zooplankton grazing during phytoplankton peak bloom (but not before and after) and (ii) to increase the fraction of nutrients that is routed through dissolved organic matter (DOM).

Table 1Parameter values of the ecosystem parameterizations that have been changed between NorESM1, NorESM-OC1.2, and NorESM2. NorESM-OC1.2 was an intermediate model version (Schwinger et al., 2016) and is included here for completeness and traceability. The naming of parameters follows Ilyina et al. (2013).

Download Print Version | Download XLSX

The first point is achieved through a reduction in zooplankton mortality and increases in the maximum grazing rate, assimilation efficiency, and half-saturation constant for zooplankton grazing. The latter point is mainly achieved by reducing the production of detritus by zooplankton (fecal pellets, 1-ϵ_zoo of grazing) and increasing exudation and excretion rates. This allows, in combination with a decrease in the DOM remineralization constant, more nutrients to be laterally transported out of regions where nutrient trapping occurs, i.e. the tropical Pacific.

Additionally, CaCO₃-to-phosphate and silicate-to-phosphate ratios were tuned to remove the high alkalinity bias of NorESM1 (see Schwinger et al., 2016, for details). Note that some adjustments of the parameters given in Table 1 were also necessary because of the new sinking scheme (increasing sinking speed with depth; see Sect. 2.6) and the different physical circulation fields between NorESM1 and NorESM2. An updated schematic diagram of the marine ecosystem module in NorESM2 is presented in Fig. 1.

Figure 1Schematic diagram of the ecosystem module in the ocean biogeochemical component of the NorESM2 model. The diagram is an updated version of a similar diagram in Six and Maier-Reimer (2006). Blue depicts components, processes, and parameters that have been modified in NorESM2 relative to the previous model version NorESM1. Dashed arrows indicate air–sea interactions, while wavy arrows depict riverine inputs.

Download

In NorESM1, biological productivity was computed only for the top 100 m of depth, which presumably represents the averaged euphotic layer. In an isopycnic model, there is no static interface separating this depth level. In cases in which the bulk MLD extends below 100 m, we virtually split this at 100 m, and the biological production is simulated only down to this interface. In NorESM2, we have omitted this virtual layer splitting, and the biological production is allowed to occur below 100 m as long as it is within the bulk mixed layer depth and has sufficient light (attenuated with depth and chlorophyll concentration) for phytoplankton growth.

3.1 Preformed tracers

Four new preformed tracers have been implemented in NorESM2, namely preformed dissolved oxygen (${\mathrm{O}}_{\mathrm{2}}^{\mathrm{pre}}$), phosphate (${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{pre}}$), alkalinity (ALK^pre), and dissolved inorganic carbon (DIC^pre). They are initialized to zero at the beginning of the spin-up. During the model integration, in the bulk mixed layer of the model (i.e. the top two levels), the preformed tracers are set to the respective total, non-preformed values at each time step. Below the mixed layer, the preformed tracers are advected as passive tracers by the physical processes and are hence a measure of transport-induced (e.g. circulation, ventilation, etc.) changes. The preformed tracers are included for diagnostic purposes, such as identifying the sources of model–data misfits (Duteil et al., 2012). The preformed oxygen can be used to quantify the total and apparent oxygen utilization (TOU and AOU) as well as O₂ disequilibrium (ΔO₂) in the interior ocean (Eqs. 10–11; Ito et al., 2004):

and

Here, saturated oxygen (${\mathrm{O}}_{\mathrm{2}}^{\mathrm{sat}}$) is determined as a function of temperature and salinity (Garcia and Gordon, 1992). Due to its closed coupling with ocean circulation, the preformed oxygen can be used to separate biologically from physically induced biases in the simulated interior biogeochemical tracers.

Following Bernardello et al. (2014), ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{pre}}$ is used to quantify the organic biologically mediated carbon pump (soft tissue, DIC^soft), whereas both ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{pre}}$ and ALK^pre are used to quantify the inorganic (carbonate, DIC^carb) biologically mediated carbon pump (Eqs. 12–14). We use a stoichiometry ratio of $\mathrm{P}:\mathrm{N}:\mathrm{C}=\mathrm{1}:\mathrm{16}:\mathrm{122}$ in the model such that DIC^soft represents remineralized particulate and dissolved organic materials, and DIC^carb represents the dissolution of particulate calcium carbonate in the water column:

and

Equation (14) is slightly different from that of Bernardello et al. (2014); i.e. it uses the term $r_{\mathrm{N}:\mathrm{P}}+\mathrm{1}$ instead of r_N:P because, in our updated code, we do not neglect the contribution of phosphate to alkalinity changes during biological production and remineralization. For instance, both nitrate and phosphate produced during organic remineralization increase the concentration of proton and therefore reduce the alkalinity (see also Sect. 3.1.3 of Paulmier et al., 2009).

The difference between total DIC and biologically altered DIC (DIC^soft + DIC^carb, or the biological carbon pump) therefore represents the physical carbon pump (DIC^pre), which comprises saturated and disequilibrium components (Eq. 15). Here, a tracer of saturated DIC (DIC^sat) has been implemented. It is computed at the surface layer as a function of atmospheric pCO₂, total alkalinity, sea surface salinity (SSS), and sea surface temperature (SST). Below the first layer, DIC^sat is treated as a passive tracer advected by the circulation. Since the equilibration timescale of the air–sea CO₂ fluxes is typically longer than the surface water mass residence time (e.g. in the deepwater formation regions), a disequilibrium term (DIC^diss) is computed by subtracting the saturated from the preformed DIC components. The DIC disequilibrium component is used to diagnose the importance of ventilation variability in the physical solubility pump and in the overall DIC storage (Eggleston and Galbraith, 2018):

3.2 Natural inorganic carbon tracers

In order to comply with the Ocean Model Intercomparison Project (OMIP) of CMIP6 (Orr et al., 2017), we have implemented a natural tracer of DIC (DIC^nat), which is formulated in the same manner as DIC^tot, except that the air–sea gas exchange is computed with a fixed pre-industrial atmospheric CO₂ concentration of 284.32 ppm. In a transient simulation with increasing atmospheric CO₂ concentrations, the difference between DIC^tot and DIC^nat represents anthropogenic DIC (DIC^ant):

The inclusion of a DIC^nat tracer also requires the implementation of respective “natural” components for both the alkalinity and the particulate inorganic carbon (CaCO₃) tracers. This is because the anthropogenic CO₂ uptake alters the carbonate system and therefore the dissolution of CaCO₃. We note that anthropogenic carbon does not influence the biological production (e.g. nutrient concentrations and phytoplankton growth rate) in our model. Similarly, a natural air–sea CO₂ flux has been added to the model output. These natural tracers are only activated in simulations with atmospheric CO₂ that departs from the pre-industrial value. Here, DIC^nat, TALK^nat, and ${\mathrm{CaCO}}_{\mathrm{3}}^{\mathrm{nat}}$ are initialized in the same manner as DIC, TALK, and CaCO₃, respectively. In transient historical and future scenario simulations, these tracers provide insight into the natural carbon evolution under anthropogenic climate change.

In addition, DIC^nat also provides a more precise estimate of DIC^ant entering the ocean since the pre-industrial period than simply taking the difference between DIC^tot at a specific time and its pre-industrial mean value. Alternatively, these natural tracers could be computed with a parallel transient simulation but with fixed pre-industrial atmospheric CO₂ as the ocean carbon cycle boundary condition. The inclusion of natural tracers saves substantial computational time, especially for high-resolution simulations. We note that the DIC^nat tracer has not been implemented in the sediment module. Hence, in very long-timescale transient simulations in which the sediment changes become substantial, the DIC^nat tracer may include some uncertainties.

3.3 Carbon isotopes

The ¹³C and ¹⁴C carbon isotope tracers have now been implemented in NorESM2. Naturally occurring stable carbon isotopes (¹²C, ¹³C, and ¹⁴C) and their relative abundances provide valuable information about both past and present climate. Changes in the ¹³C:¹²C ratio are used e.g. to (i) study glacial–interglacial atmospheric CO₂ changes from ice core air samples (Broecker and McGee, 2013; Bauska et al., 2016), (ii) reconstruct bottom water oxygen concentrations (Hoogakker et al., 2015), (iii) reconstruct paleo-deepwater circulation (Toggweiler, 1999; Curry and Oppo, 2005; Crucifix, 2005), (iv) investigate oceanic anthropogenic CO₂ uptake involving the Suess effect (Gruber and Keeling, 2001; Quay et al., 2003; Holden et al., 2013), and (v) evaluate model sensitivity and performance (Braconnot et al., 2012; Schmittner et al., 2013). The ¹⁴C:¹²C ratio is used as a circulation and age tracer (e.g. Skinner et al., 2017). Globally, the ${}^{\mathrm{12}}\mathrm{C}{:}^{\mathrm{13}}\mathrm{C}{:}^{\mathrm{14}}\mathrm{C}$ isotope ratio is about $\mathrm{99}:\mathrm{1}:\mathrm{10}\times {\mathrm{10}}^{-\mathrm{12}}$. In order to allow for comparison between different carbon isotope studies, the ¹³C:¹²C ratio is standardized and expressed in per mille as δ¹³C (Eq. 24) (Zeebe and Wolf-Gladrow, 2001). Our implementation of ¹³C follows the OMIP guidelines (Orr et al., 2017). The ¹⁴C is standardized as Δ¹⁴C (see Sect. 3.3.3). Variations in δ¹³C are due to isotopic fractionation processes during air–sea gas exchange, photosynthesis, and CaCO₃ formation, but the latter is neglected in our implementation (see also Sect. 3.3.2). For ¹⁴C, each fractionation factor is the quadratic of the respective ¹³C value (i.e. κ¹⁴C=κ¹³C²). Lastly, ¹⁴C is radioactive and decays with a half-life of 5730 years to ¹⁴N following

where

X represents any ¹⁴C state variable. In the model, all ¹⁴C tracers are normalized to prevent numerical errors from carrying values close to the precision of the model.

The newly implemented marine carbon isotope code parallels the respective total DIC code, and in addition includes fractionation during photosynthesis and air–sea gas exchange (as well as decay for the ¹⁴C tracers). In addition to the dissolved inorganic tracers, the following ¹³C and ¹⁴C state variables have also been implemented: dissolved organic carbon, particulate organic carbon, calcium carbonate, phytoplankton, and zooplankton. Therefore, 12 new tracers are added (isotopic DOC, POC, CaCO₃, zooplankton, and phytoplankton). Due to the long equilibration time, the isotopic tracers are only activated in the computationally efficient configurations, such as the ocean carbon cycle stand-alone configuration (NorESM2-OC) or the low-resolution version of the coupled model. Equilibrium times of the carbon isotopes are long due to the slow air–sea gas equilibration (Broecker and Peng, 1974) and long-term transient effects from the balance between sediment burial and weathering (Roth et al., 2014). Realistic equilibration times are therefore currently only obtained when bypassing the sediment module of the model when running the carbon isotopes, as well as omitting carbon isotope input from rivers. When bypassing the sediment, mass balance is maintained by redistributing POC fluxes to the sediment over the entire overlying water column and by dissolving inorganic carbon fluxes and opal fluxes at the bottom of the model immediately. Ongoing work will add the possibility for a fast sediment spin-up for use in future versions of the biogeochemical ocean model. Since Δ¹⁴C is governed mainly by radioactive decay and circulation, we focus on δ¹³C in our description of isotopic fractionation effects (Sect. 3.3.1 and 3.3.2).

5.1 Statistical performance summary

We assess, relative to the observational estimates and the earlier version NorESM1, the simulated long-term annual mean of global hydrography and biogeochemical tracer distributions. Where relevant, the mean seasonal cycles are also evaluated for the surface values. The list of parameters, the averaging periods, and the respective observational data used to assess the model performance are given in Table 2. We use spatial correlation and normalized RMSE (root mean squared error) metrics to measure the model–data difference and determine whether or not the current model version has improved and performed better than the earlier version. The normalization was done by dividing the RMSE with the standard deviation of the respective observations.

Table 2List of the simulated variables and the respective historical simulation periods over which their climatological values have been averaged. The last column indicates the observational reference used to validate the model output.

Download Print Version | Download XLSX

Figure 2b shows that, except for the surface layer, the NRMSE of most of the NorESM2 variables' mean climatology is either noticeably improved or relatively similar to that of NorESM1. At 500 m of depth, all variables but salinity have lower RMSE values. For many of the biogeochemical tracers (phosphate, nitrate, oxygen, and DIC), NorESM2 shows better agreement with data throughout the water column. Similarly, the respective spatial correlations with the observations are considerably improved for most variables, especially for nitrate, DIC, and alkalinity (Fig. 2a). For surface seasonality, NorESM2 performs fairly comparable to NorESM1, with improvements in all seasonal metrics (NRMSE and spatial pattern) for surface salinity and net primary production. NorESM2 simulates slightly larger NRMSE in its surface nutrients (phosphate, nitrate, and silicate) for both their annual and seasonal average values relative to NorESM1. This is attributed to the anomalously high concentration in the Arctic Ocean and parts of the Southern Ocean (just north of the Antarctic Circumpolar Current). More details on the performance of each specific variable are discussed in the following subsections.

Figure 2Summary of (a, c) spatial correlation coefficients and (b, d) normalized RMSE between the observations and models (NorESM1-ME, NorESM2-LM, and NorESM2-MM). In all (a, b) 3D fields, the metric values are computed over a 2D horizontal global domain ranging from the surface to 3 km of depth, while for (c, d) seasonal metrics, only surface values are evaluated. Each square contains three metric scores, i.e. for NorESM1-ME (top), NorESM2-LM (bottom-left), and NorESM2-MM (bottom right). Black dots indicate which model performs the best.

Download

5.2 Temperature

NorESM2 simulates a warm bias for both pre-industrial conditions and the contemporary surface ocean (Fig. 3d), where a warm bias as high as 5 ^∘C is simulated in some regions of the Southern Ocean. Cold biases are simulated in parts of the north and equatorial Pacific, North Atlantic subpolar gyre, and the Arctic Ocean. The climatological mean temperatures are broadly similar between the two NorESM versions, with NorESM2 having fewer spatial correlations at 1 and 1.5 km depths (Fig. 2a, b). At depths below 2 km, with the exception of the North Pacific, NorESM2 simulates fewer biases (Fig. 3e, f, h, i). This improvement is related to the more accurately simulated Antarctic Bottom Water (AABW; see also Supplement Fig. S1) in both the Pacific and Atlantic basins, as also depicted in Fig. 3. While NorESM1 simulates an AMOC strength that is too strong (roughly 31 Sv at 26^∘ N), NorESM2 has a more reasonable strength of approximately 21 Sv. Nevertheless, the warm bias of 1–2 ^∘C in the North Atlantic Deep Water (NADW) remains when compared to observations. The cold bias seen previously in the intermediate water mass in the Southern Ocean and North Pacific has been improved in the new version. In the Atlantic sector just south of 30^∘ N at 1 km of depth, there is a warm bias of up to 5 ^∘C (Fig. 3e), which could be related to the anomalously strong overflow water mass from the Mediterranean Sea. Nevertheless, the cold bias in the Atlantic and Pacific tropical–subtropical thermocline simulated in NorESM1 is now improved in NorESM2. Temperature values from NorESM2-LM are shown in Supplement Fig. S2.

Figure 3Climatological temperature values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

5.3 Salinity

The RMSE in surface salinity is reduced in NorESM2, mostly owing to improvements in the Arctic, where the previous model version simulates water that is too saline (Fig. 4d, g). Also at the surface, NorESM2 simulates biases that are anomalously too fresh and too saline in the Pacific and the Atlantic basins, respectively. In the subtropical South Pacific, the negative bias is as high as −2 psu and is consistent with a precipitation rate that is too strong in this region simulated by the atmospheric component (not shown). Positive biases are simulated in the northern Indian Ocean and along the Benguela current extension. In the interior, the salinity bias in the Pacific and the Atlantic is of the order of ±0.11 and +0.4 psu, respectively (Fig. 4e, f). In the NADW, NorESM2 displays a positive bias, which is, however, smaller than that of NorESM1 (Fig. 4h; bias greater than 0.5 psu). Similar to temperature bias, NorESM2 produces too much saline water around 30^∘ N at 1 km of depth in the Atlantic (Fig. 4e) due to overflow of high-salinity Mediterranean water masses. In the interior Pacific and Southern Ocean, both models' performances are generally comparable (Fig. 4e, f, h, i). Except for the surface salinity, the simulation with lower atmospheric resolution (NorESM2-LM) depicts a similar salinity pattern throughout the water column (Supplement Fig. S3).

Figure 4Climatological salinity values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

5.4 Mixed layer depth

Seasonally averaged mixed layer depths (MLDs) from NorESM2 and observations are shown in Fig. 5. To be consistent with the observational estimates (de Boyer Montégut et al., 2004), we have computed our MLD using the σ_T (density) criteria, which is the first depth at which the change from surface σ_T of 0.03 has occurred. In NorESM1, the MLD is generally too deep throughout most ocean regions. The new improvements in the MLD parameterization in NorESM2 reduce this bias considerably, especially in the low and middle latitudes as well as the North Pacific regions. In the subpolar North Atlantic and parts of the Southern Ocean, the simulated MLD is still deeper than the observation-based estimates. In the Weddell Sea, NorESM2 also persistently simulates a deep MLD throughout the year, a feature not seen in the observations. There is no significant difference between MLD simulated in NorESM2-MM and NorESM2-LM (Supplement Fig. S4).

Figure 5Climatological seasonal mixed layer depth as (left column) simulated in NorESM2-MM and (right column) estimated from observations.

5.5 Ocean ventilation

To assess the simulated ventilation, we compare the passive chlorofluorocarbon (CFC-11) tracer distribution in the interior Atlantic and Pacific Ocean as simulated by NorESM2 with that from observations, as shown in Fig. 6. We compare the simulated CFC-11 from calendar year 2000 with climatological estimates of Key et al. (2004). Similar to the observations, high values of CFC-11 are generally simulated in the upper 1 km in both the Atlantic and the Pacific, with the exception of the North Atlantic, where up to 0.5 nmol CFC-11 m⁻³ penetrates down to 5 km of depth. In the mid-latitude North Atlantic (i.e. 30^∘ N), the model is unable to simulate the observed high values at around 1 km of depth. This is likely due to the discrepancy in the water mass origin between the model and observations. Here, the water mass in the model is too old, as also seen in the overly high apparent oxygen utilization and dissolved inorganic carbon (see the subsections below on oxygen and DIC). The model simulates deepwater ventilation that is too strong in both the Atlantic and Pacific sectors of the Southern Ocean as well as in the high latitudes of the North Atlantic. Similar ventilation patterns are also produced by NorESM2-LM (Supplement Fig. S5), with slightly less ventilated North Atlantic.

Figure 6Concentration of CFC-11 across the vertical sections of (a) the Atlantic and (b) Pacific from NorESM2-MM (colour shadings) and observations (contour lines). Differences between NorESM2-MM and observations are shown in panels (c)–(d).

5.6 Nutrients

In general, the climatological nutrient distributions (phosphate, nitrate, and silicate) in NorESM2 are either comparable or improved compared to those of NorESM1 in both spatial correlation and normalized RMSE (Fig. 2a, b). At depth, these improvements are due to the new particulate organic carbon sinking scheme, whereby the sinking velocity increases linearly with depth. In NorESM1, the sinking speed is constant with depth, leading to more organic material being remineralized in the upper 1.5 km. The sinking scheme in NorESM1 produces anomalously high regenerated phosphate at intermediate depth (approximately between 100 and 1500 m of depth; Fig. 7h, i). A consequence of this is an oxygen minimum zone that is too strong (OMZ; see also the next subsection) and denitrification in the low latitudes, where the nitrate concentration becomes overly depleted. On the other hand, at depths greater than 1500 m, NorESM1 exhibits concentrations of all nutrients that are too low relative to the observations (panels h and i in Figs. 7–9).

Figure 7Climatological phosphate values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

Figure 8Climatological nitrate values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

Figure 9Climatological silicate values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

With the new sinking scheme, more organic material reaches the deep ocean and gets remineralized there. As a result, model–data biases in phosphate and nitrate in the interior are noticeably improved (Figs. 2b, 7e, f, 8e, f). In the Pacific OMZ, reduced export production in NorESM2 also contributes to the reduced nutrient bias. In the interior North Atlantic, despite strong positive spatial correlation with the observations, NorESM2 shows a positive bias at depth below 3000 m and north of 30^∘ S (Figs. 7e, 8e, 9e). Through analysing the quasi-conservative “PO” tracer (Broecker, 1974), we are able to attribute this to the bias in the simulated water mass (Supplement Fig. S1). In NorESM2, this region is dominated by Antarctic Bottom Water, whilst the observations indicate NADW as the main water mass. As with NorESM1, the current model still simulates nutrient concentrations that are too low in the North Pacific Intermediate Water mass (Figs. 7f, 8f, 9f). This is likely associated with the excessively low surface primary production simulated in this region (see also Sect. 5.8), leading to limited organic matter available for remineralization at depth. We note that the circulation bias in this region could also contribute to the nutrient bias. In the Southern Ocean, overly strong ventilation potentially leads to negative nutrient biases in both the Atlantic and Pacific sectors.

In NorESM2, less export production and deeper distribution of POC through the changed sinking scheme lead to a reduced rate of denitrification in the Pacific oxygen minimum zone and a greatly improved nitrate distribution (Fig. 8f). In NorESM1, far too much nitrate was consumed for denitrification (Fig. 8i). The silicate spatial distribution in NorESM2 closely resembles the phosphate distribution (Figs. 9a–c and 7a–c), with similar biases across the vertical sections of the Atlantic and Pacific (Figs. 9e, f and 7e, f). At high latitudes, NorESM2 overestimates the surface silicate concentration (Fig. 9d), which could be attributed to the reduced opal export sinking speed relative to our earlier model version (30 instead of 60 m d⁻¹ previously).

In NorESM, the bulk phytoplankton growth rate is limited by multiple nutrients (i.e. phosphate, nitrate, and dissolved iron), in addition to temperature and light (see also Sect. 2.3.2 of Schwinger et al., 2016). The model does not explicitly simulate diatom and calcifier classes. The inclusion of iron is critical to simulate the observed HNLC (high-nutrient low-chlorophyll) regions in the world ocean where year-long elevated concentrations of both nitrate and phosphate are not exhausted by phytoplankton (de Baar et al., 1995; Martin and Fitzwater, 1988). The three major HNLC regions are the subarctic North Pacific, the eastern and central equatorial Pacific, and the Southern Ocean. Here, low levels of bioavailable iron concentrations limit phytoplankton growth. In NorESM1, only the subarctic North Pacific is simulated to be iron-limited, with the other two HNLC regions being mostly limited by nitrate, as shown in Fig. 10. Iron limitation is also incorrectly simulated in parts of the subtropical North Pacific, South Pacific, and most parts of the western Pacific. Following improvements in the iron parameterization, the three HNLC regions are now shown to be iron-limited in NorESM2. Nevertheless, the iron limitation bias in the western equatorial Pacific remains.

Figure 10Maps of limiting nutrients for biological productivity as simulated in (a) NorESM1 and (b) NorESM2. The limiting nutrient is determined as the minimum between PO₄, $R_{\mathrm{P}:\mathrm{N}}\cdot {\mathrm{NO}}_{\mathrm{3}}$, and $R_{\mathrm{P}:\mathrm{Fe}}\cdot \mathrm{Fe}$.

5.7 Dissolved oxygen

At the surface level, dissolved oxygen is mostly constrained by solubility and hence sea surface temperature. Both NorESM1 and NorESM2 perform similarly for dissolved oxygen at the surface, with subtle improvements in the seasonal spatial correlation coefficient of NorESM2 (Fig. 2a, c). As with nutrients, the climatological distributions of interior oxygen improved considerably when compared to NorESM1 (Fig. 11), which is one of the CMIP5 ESMs with an overly strong OMZ (Cabré et al., 2015). In NorESM2, the OMZ volume is reduced considerably in both the equatorial Atlantic and Pacific, and the absolute dissolved oxygen values are also much improved (Fig. 11e, f). Here the remaining bias is related to the nutrient trapping issue (Six and Maier-Reimer, 1996), as can also be seen from excessively high apparent oxygen utilization (AOU; Fig. 12d). In the interior Pacific and Atlantic roughly at 2 km and deeper, NorESM2 simulates a considerably higher oxygen concentration than the observations by as much as 50 % (North Pacific Deep Water; Fig. 11f). In the North Pacific Deep Water, this is primarily due to too little organic matter flux reaching below the mesopelagic zone, leading to apparent oxygen utilization that is too low (Fig. 12d). The uncertainty in the dissolved organic carbon (not shown), which is currently not well constrained due to lack of observational data, could also play a role in this. In addition to biological constraints, the remaining interior oxygen bias can be attributed to the overly strong Antarctic Bottom Water (AABW) ventilation, which carries a colder (Fig. 3e, f) and higher oxygen-saturated water mass into both the Atlantic and Pacific bottom waters (Fig. 11e, f).

Figure 11Climatological dissolved oxygen values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

5.8 Biological production

In NorESM1, the primary production is generally too strong in the equatorial Pacific and during the summer months at high latitudes. In the subtropical oceans, the simulated production is too low relative to the observational estimates. With the newly tuned ecosystem parameters of NorESM2, the spatial productivity patterns are improved considerably, with an average absolute bias in the open ocean of less than 5 mol C m⁻² yr⁻¹ (Fig. 13a). In NorESM1 these biases reach more than 10 mol C m⁻² yr⁻¹, especially in the equatorial Pacific, the South Atlantic, and the Southern Ocean. Seasonally, the spatial patterns are also improved when compared to the observations, with a correlation of approximately r=0.5 (Figs. 2c, 13c). The NRMSE has also been reduced considerably, particularly in the October–November–December months when NorESM1 simulates a spring bloom that is much too strong in parts of the Southern Ocean, which is overestimated by as much as a factor of 8 (e.g. Figs. 13e, 14e and Nevison et al., 2015). Similarly, the spring blooms in the North Atlantic and North Pacific are overestimated in NorESM1 (Figs. 13e, 14a). The improved spring bloom characteristics at high latitudes, however, come at the cost of annual mean productivity that is too low in the North Atlantic and North Pacific.

Figure 12Climatology of apparent oxygen utilization (AOU) values across the vertical sections of (a) the Atlantic and (b) Pacific in NorESM2 (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (c)–(f).

Figure 13Annual mean (a) vertically integrated contemporary primary production and (b) transfer efficiency (T_eff-1 km) as simulated in NorESM2 (colour shadings) and from observations (contour lines). Respective difference between NorESM2-MM (NorESM1-ME) and observations are shown as colour shadings (contour lines) in panels (c)–(d).

Figure 14Mean seasonal cycle of vertically integrated contemporary primary production in the (a) North Atlantic (north of 20^∘ N), (b) North Pacific (north of 20^∘ N), (c) tropics (20^∘ S–20^∘ N), (d) mid-latitude Southern Hemisphere (20^∘ S–40^∘ S), (e) Southern Ocean (40^∘ S–60^∘ S), and (f) polar Southern Ocean (south of 60^∘ S). Shown are values from NorESM1-ME (black), NorESM2-MM (red), and observational estimates (blue). Numbers within each panel depict the regional mean values averaged over all months.

Due to unresolved physical dynamics, the model is still unable to reproduce the observed high productivity in coastal regions. Globally, the contemporary total primary production in NorESM2 (NorESM1) is 35.3 (39.9) Pg C yr⁻¹. These estimates are lower than the observational (MODIS, excluding coastal areas) estimates of 46.1 Pg C yr⁻¹ but within the broad range of CMIP5 models (Bopp et al., 2013). In our model, export production is estimated as the flux of particulate organic carbon exiting the base of the euphotic zone (here assumed to be 100 m of depth) as a function of zooplankton and phytoplankton mortalities and zooplankton faecal materials. Table 3 shows that the export production in NorESM2 is 5.39 Pg C yr⁻¹ compared to 7.90 Pg C yr⁻¹ simulated by NorESM1. Current estimates of the global export production remain highly uncertain, from 5 to > 12 Pg C yr⁻¹, with a more recent satellite-based approach that accounts for food-web processes revealing values of 6±1.2 Pg C yr⁻¹ (Siegel et al., 2014).

Table 3Annual mean biology-related metrics simulated by NorESM1-ME, NorESM2-LM, and NorESM2-MM.

Download Print Version | Download XLSX

Table 3 also summarizes other export production metrics from both model versions. As a consequence of the new particle sinking scheme, there is more POC exported into the deep ocean. This is especially reflected by the values of transfer efficiency (T_eff-1 km), which is calculated as a ratio between POC fluxes at a depth of 1000 and 100 m. The T_eff-1 km in NorESM2 (0.24) is increased by a factor of 4 relative to that in NorESM1 (0.06). Compared to a recent estimate of transfer efficiency reconstructed from observed interior biogeochemistry (Fig. 13b; Weber et al., 2016), NorESM2 still simulates transfer efficiencies of >50 % that are too high in the eastern equatorial Pacific (Fig. 13d). On the other hand, the transfer efficiency in the Southern Ocean and northern high latitudes is comparable with observations. This represents a significant improvement when compared to NorESM1 (Fig. 13d, f), which simulates nearly uniform T_eff-1 km values of less than 0.1 in all ocean regions with the exception of the eastern equatorial Pacific and equatorial Atlantic.

Figure 14 shows the seasonal cycle of biological production rates as simulated by NorESM1 and NorESM2 averaged over different ocean regions together with observation-based estimates. The regional mean values in NorESM2 are generally lower than those in NorESM1, except in the North Pacific and the tropics (Fig. 14b, c). Nevertheless, Fig. 2d shows that NorESM2 simulates a lower normalized RMSE than NorESM1 in all seasons. As stated above, this bias could be reduced even more if we excluded the continental shelf regions. Nevertheless, Fig. 14 shows that the seasonal phase and amplitude of biological production rates in NorESM2 in the North Atlantic, North Pacific, and the Southern Hemisphere regions are closer to the observations when compared to NorESM1. In the Southern Hemisphere (south of 20^∘ S), NorESM2 no longer produces the large summer bias seen in NorESM1.

5.9 DMS production and fluxes

In NorESM2, DMS production is tuned towards the climatological observations from Lana et al. (2011). Figure 15 shows the simulated surface DMS concentration for each season together with the observation-based estimates. To first order, the DMS distribution follows the spatial pattern of the simulated primary production, with higher concentrations during spring–summer periods in high-latitude regions. In the equatorial Pacific and Atlantic upwelling regions, the model overestimates DMS concentrations. Similarly, in the North Pacific, where the model underestimates productivity during the boreal spring bloom (see also Fig. 14b), the DMS concentration is lower than observed (Fig. 15b, f).

Figure 15Seasonally averaged surface concentration of DMS from observations (a–d) and as simulated in NorESM2-MM (e–h).

The simulated DMS concentration in the North Atlantic is comparable with observations. In the Southern Ocean, the model produces the observed high concentrations during summer months (JFM) along the 40–45^∘ S latitude band well (Fig. 15e). In winter periods, the DMS concentration approaches zero in the model owing to the excessively low productivity, whilst the observations still indicate values of approximately 1 µmol S m⁻³. For the present-day period (1971–2000), the simulated global DMS flux is 19.76±0.19 Tg S yr⁻¹, which is at the lower end of observation-based estimates (17.6–34.4 Tg S yr⁻¹).

5.10 Dissolved inorganic carbon and alkalinity

Except for the polar Southern Ocean, the surface alkalinity concentrations in NorESM1 are approximately 100 µmol L⁻¹ higher than the observational estimate. And since the spin-up was forced with constant pre-industrial atmospheric CO₂ concentrations, the surface DIC concentration is adjusted (also anomalously high bias; Figs. 16g and 17g) to yield the “correct” surface pCO₂. Consequently, biases in surface DIC and alkalinity compensate for one another to yield the approximately correct carbonate ion concentration (i.e. alkalinity minus DIC), seawater CO₂ buffer capacity, and air–sea CO₂ fluxes.

Figure 16Climatological dissolved inorganic carbon (DIC) values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

Figure 17Climatological alkalinity values at the (a) surface and across the vertical sections of (b) the Atlantic and (c) Pacific in NorESM2-MM (colour shadings) and from observations (contour lines). Differences between NorESM2-MM (or NorESM1-ME) and observations are shown in panels (d)–(i).

In NorESM2, we alleviate the surface alkalinity bias by increasing the surface sink associated with calcium carbonate production. This is done by increasing the CaCO₃-to-phosphate uptake ratio ($R_{{\mathrm{CaCO}}_{\mathrm{3}}:\mathrm{P}}$ in Table 1). As a result, the export of particulate CaCO₃ increases from 0.49 Pg C yr⁻¹ in NorESM1 to 0.66 Pg C yr⁻¹, with both estimates still within the range of the observation-based synthesis of 0.52±0.15 Pg C yr⁻¹ (Dunne et al., 2007). This modification increases the PIC:POC export ratio to 0.13 (0.06 in NorESM1; see also Table 3) and considerably improves the surface DIC and alkalinity concentrations in NorESM2 (Figs. 16a, d and 17a, d). At the same time, the spatial correlations with observations are also improved, as is the normalized RMSE (see Fig. 2a, b). We note that in a few regions, e.g. the subtropical South Pacific, the model still underestimates the observed DIC and alkalinity concentrations. We note that the PIC:POC ratio is still within the large range of other studies (0.03 to 0.25; Koeve, 2002; Sarmiento et al., 2002)

As at the surface, the spatial distributions and normalized bias of interior DIC and alkalinity in the NorESM2 improve considerably throughout the water column relative to NorESM1 (Fig. 2a, b). In both the interior Atlantic and Pacific basins, the spatial distribution and magnitude of DIC biases closely resemble that seen in the nutrients (e.g. phosphate; Figs. 16e, f, 7e, f) when multiplied by the constant stoichiometry in the model ($R_{\mathrm{C}:\mathrm{P}}=\mathrm{122}$). This suggests that the mechanisms driving the nutrient bias are also responsible for the bias in interior DIC. For instance, the positive anomaly in the equatorial Pacific between 1 and 3 km depths originates from too much biological remineralization in the model.

The simulated interior alkalinity concentrations in NorESM2 between 1 and 3 km depths are anomalously high, which can be partially attributed to excessively low interior remineralization in the Atlantic (see also AOU values Fig. 12c) or overly high CaCO₃ export production (Table 3).

5.11 Surface pCO₂ and sea–air CO₂ fluxes

The NorESM2 pCO₂ spatial correlation relative to the observations is improved in almost all seasons when compared to NorESM1, while the NRMSE is reduced for JFM and AMJ months (Fig. 2c, d). Within these months, improvements are mostly seen in the Southern Ocean, where surface mixing that is too strong in NorESM1 leads to anomalously high pCO₂. The climatological mean of contemporary surface pCO₂ in NorESM2 compares well with the observational compilations, as seen in Fig. 18a and e. Similar improvements are also seen in the CO₂ fluxes, for which NorESM2 simulates the broad spatial patterns relatively well (Fig. 18c, g). NorESM2 also produces a moderately weaker annual mean carbon flux in the northern middle to high latitudes than NorESM1. Figure 18f and h show the zonally averaged monthly surface pCO₂ and CO₂ fluxes in NorESM2. Here, the model also agrees well with the observation-based estimates (Fig. 18b, d) in terms of amplitude and temporal variability, with a distinct seasonal cycle between 20 and 45^∘ south–north. In the extratropical oceans (between 45 and 65^∘ N and south of 60^∘ S), the simulated winter (January–March and July–September in the Northern and Southern Hemisphere, respectively) pCO₂ is considerably lower than that seen in the observations. Consequently, the simulated carbon fluxes in these regions are in the opposite direction of those seen in the observations (Fig. 18d, h).

Figure 18Observed surface climatology maps of (a) the partial pressure of CO₂ and (c) sea–air CO₂ fluxes together with Hovmöller plots of the zonally averaged seasonal cycle of (b) surface pCO₂ and (d) sea–air CO₂ fluxes. Corresponding values from the NorESM2-MM simulation are shown in panels (e–g).

Figure 19Mean seasonal cycle of vertically integrated contemporary sea-to-air CO₂ fluxes in the (a) North Atlantic (north of 20^∘ N), (b) North Pacific (north of 20^∘ N), (c) tropics (20^∘ S–20^∘ N), (d) mid-latitude Southern Hemisphere (20–40^∘ S), (e) Southern Ocean (40–60^∘ S), and (f) polar Southern Ocean (south of 60^∘ S). Shown are values from NorESM1-ME (black), NorESM2-MM (red), and observational estimates (blue). Numbers within each panel depict the regional mean values averaged over all months.

Figure 20Time series of global mean annual (a) surface air temperature anomalies relative to the 1850–1879 period, (b) surface ocean and atmospheric pCO₂, (c) ocean CO₂ uptake, and (d) atmospheric CO₂ growth rate. Shown are values from the NorESM2-MM (grey) pre-industrial control and (orange) historical simulations. Green depicts results from the emissions-prescribed historical simulation with NorESM2-LM. The purple line in panel (c) represents the anthropogenic carbon uptake estimated by subtracting the natural carbon uptake from the total uptake (see also Sect. 3.2). Data for temperature and ocean carbon sinks (in blue) are estimates from HadCRUT4 and the Global Carbon Project (Morice et al., 2012; Le Quéré et al., 2018), respectively.

Download

Figure 19 depicts the mean seasonal cycle of sea–air CO₂ fluxes from NorESM models and observations, averaged over different regions. In the North Atlantic, North Pacific, and the mid-latitude Southern Hemisphere (between 20 and 40^∘ S), the seasonal phase and mean integrated annual CO₂ fluxes in NorESM2 compare better with the observations than those in NorESM1 (Fig. 19a, b, d). Nevertheless, NorESM2 simulates a stronger than observed carbon sink in the Labrador Sea and along the Kuroshio current extension during boreal winter months. We note that the Landschützer et al. (2014) data estimate higher surface pCO₂ (hence a weaker carbon sink) in these regions than the observational estimates from Takahashi et al. (2009). Despite the warm SST bias in most of the tropics, NorESM2 simulates a moderately weaker than observed CO₂ outgassing throughout the year (red line compared to blue line in Fig. 19c). A potential explanation for this is the weaker upwelling rates of DIC-rich deepwater masses in the model. In the Southern Ocean (south of 40^∘ S), the strong bias in the seasonal amplitude in NorESM1 is considerably reduced (Fig. 19e, f) due to the improvements in the biological production in this region. Nevertheless, here the seasonal phase in NorESM2 is approximately reversed compared to observations: the strongest sink is simulated in the austral winter months, whereas observational estimates indicate summer months.

5.12 Transient changes

In this subsection, we describe the global mean surface air temperature and oceanic carbon uptake over the historical period (1850–2014) as simulated by NorESM2. For this, we discuss the two transient simulations: (i) historical and (ii) esm-hist. Similar to SST (see Fig. 3), there is a warm bias in the simulated surface air temperature. For the period of 1850–1859, the NorESM2 global mean surface air temperature is 14.5 ^∘C compared to 13.7 ^∘C from the HadCRUT4 dataset (Morice et al., 2012). Nevertheless, the simulated warming between 2005–2014 and 1850–1859 is 0.7 ^∘C, comparable to that from observations of 0.8 ^∘C (Fig. 20a). The warming in esm-hist is lower at 0.6 ^∘C. Figure 20a also shows that NorESM2 (orange and green lines) simulates warming rates that are comparable to observations (blue line) between 1970 and 2010, but it shows a stronger cooling between 1950 and 1970 (historical and esm-hist).

In both simulations, the global mean surface ocean pCO₂ follows the respective atmospheric CO₂ trend over the simulation periods. Towards the end of the historical period, the oceanic pCO₂ (solid green and orange lines) diverges from, i.e. is lower than, the atmospheric counterpart (dashed green and orange lines) by approximately 20 ppm (Fig. 20b). Consistent with the lower surface ocean (than the atmosphere) pCO₂, the ocean carbon uptake continues to increase over most of the transient period in both simulations (Fig. 20c). In the historical simulation, the carbon sink stabilizes at roughly 0.7 Pg C yr⁻¹ between 1910 and 1960. This flattening of the uptake strength is associated with the slowdown in the atmospheric CO₂ growth rate. Figure 20d confirms that the atmospheric growth rate weakens during this time window (orange bars in Fig. 20d). In contrast, the carbon sink in esm-hist steadily increases during the same period, also consistent with the steady increase in the simulated prognostic atmospheric CO₂ (green bars in Fig. 20d).

From 1960 onward, the atmospheric CO₂ growth rates in both simulations increase almost linearly in time from 0.8 to more than 2.0 ppm yr⁻¹. Consequently, the ocean carbon uptake strengths increase as well, from approximately 1 Pg C yr⁻¹ in the 1960s up to roughly 2.5 Pg C yr⁻¹ in 2014 (in both historical and esm-hist). Rapid strengthening of ocean carbon uptake is also evident in estimates from the Global Carbon Project (GCP, blue line in Fig. 20; Le Quéré et al., 2018). For the 1980s and 1990s, the simulated carbon uptake is also well within the IPCC (Intergovernmental Panel for Climate Change, Denman et al., 2007) estimates (Table 4). In NorESM2, a new tracer, DIC^nat, has been implemented (see also Sect. 3.2), and it enables us to more accurately estimate the anthropogenic carbon uptake and storage in the ocean. Using this tracer, we can also estimate the net flux of anthropogenic carbon into the ocean during the historical period, as depicted by the purple line in Fig. 20c. The strong resemblance between the orange and purple lines suggests that the long-term trend in oceanic carbon sinks is associated with the increasing atmospheric CO₂ rather than changes in climate states (as expected for the historical runs in which the change in the climate state is relatively small). We note that, in the pre-industrial control simulation, the ocean is a weak net carbon sink of approximately 0.1 Pg C yr⁻¹. Table  4 also shows that the cumulative carbon uptake by the ocean for 1850–1994 and 1994–2011 is well within the observation-based estimates. The cumulative carbon uptake in NorESM2 is lower than that of NorESM1 and compares better with observations.

Table 4Annual and accumulated ocean carbon uptakes simulated by NorESM1 and NorESM2 as well as from observation-based estimates.

^* 1800–1994.

Download Print Version | Download XLSX

Figure 21Climatological δ¹³C values at (a) 500 m of depth and across the vertical sections of (b) the Atlantic and (c) Pacific from NorESM2-OC (colour shadings) and observations (contour lines). Differences between NorESM2-OC and observations are shown as colour shadings in panels (d)–(f).

5.13 Distribution of δ¹³C and Δ¹⁴C

As stated above, the carbon isotope tracers are not included in the coupled ESM configuration and only simulated in a pre-industrial stand-alone ocean configuration. We also note that although the stand-alone ocean simulation is not a direct representative of the coupled ESM simulation, its equilibrium ocean state is very comparable. For instance, the interior distributions of phosphate and DIC concentration in the stand-alone ocean simulation (Supplement Figs. S6, S7) share very similar characteristics with that produced in the coupled configuration (Figs. 7a–c, 16a–c). At the end of pre-industrial spin-up, Fig. 21a–c show that the spatial pattern of δ¹³C is in fair agreement with the gridded pre-industrial observational estimate of Eide et al. (2017). Nevertheless, the ocean component of NorESM2 simulates weaker variability than observed; i.e. it is unable to reproduce very low (negative) and high (positive) values. Therefore, the simulated δ¹³C concentrations in the relatively younger water masses, such as the North Atlantic Deep Water (NADW) and the Antarctic Intermediate Water (AAIW), have a negative bias relative to the observations (Fig. 21e). In contrast, positive biases are depicted in the older Antarctic Bottom Water (AABW) and in the deep North Pacific (Fig. 21e, f). The δ¹³C can be decomposed into biological and residual (air–sea gas exchange and circulation) components (Broecker and Maier-Reimer, 1992; Eide et al., 2017):

where δ¹³C_BIO can be estimated as a function of phosphate and DIC distribution and their global mean values:

Applying this decomposition, we can attribute the majority of the differences between model and observational estimates of δ¹³C to biases in the biological component (Fig. 22). The δ¹³C_BIO reveals that in the Southern Ocean (south of ∼55^∘ S), as well as in its exported deep waters, NorESM2-OC simulates δ¹³C that is too high due to the excessively small regenerated signal in these waters. The negative δ¹³C anomaly in the North Atlantic of about −0.5 ‰ is the combination of a δ¹³C_BIO signature that is too low in NADW (Fig. 22e) and an overly negative influence from air–sea gas exchange and circulation compared to observational estimates (the residual term δ¹³C_AS). In the Pacific the positive anomaly from the Southern Ocean persists throughout most of the basin due to excessively weak remineralization, except at intermediate depths at which remineralization is high (see also AOU; Fig. 12).

Figure 22Climatological fields of the biological component of δ¹³C (i.e. δ¹³C_BIO) at (a) 500 m of depth and across the vertical sections of (b) the Atlantic and (c) Pacific from NorESM2-OC (colour shadings) and observations (contour lines). Differences between NorESM2-OC and observations are shown as colour shadings in panels (d)–(f).

Figure 23Climatology Δ¹⁴C values at (a) 500 m of depth and across the vertical sections of (b) the Atlantic and (c) Pacific from NorESM2-OC (colour shadings) and observations (contour lines). Differences between NorESM2-OC and observations are shown as colour shadings in panels (d)–(f).

The simulated Δ¹⁴C distribution is in general agreement with the pre-industrial observational estimate by Key et al. (2004). Nevertheless, interior model–observation biases are in the range of ±50 ‰. Waters south of 60^∘ S are positively biased compared to observations (Fig. 23), indicating waters in this region that are too young and well ventilated. This positive bias of 20 ‰–50 ‰ is carried into the Atlantic and Pacific basins (Fig. 23e–f) and is equivalent to an underestimation of radiocarbon age by 200–400 years. North of the Equator, a negative bias of similar magnitude indicates water masses that are too old in the Atlantic basin below 3 km of depth as well as throughout the water column in the northern half of the Pacific basin. The top 800 m of the water column has a negative bias of 100 ‰, which we attribute to the overly negative atmospheric p¹⁴CO₂ before calibration (−36 ‰) and biases in air–sea equilibration and fractionation. Calibrated Δ¹⁴C_atm is by definition 0 ‰, δ¹³C_atm equilibrates at −7.5 ‰, and atmospheric pCO₂ is 294 ppm at the end of the simulation.