Marine biogeochemical cycling and climate-carbon cycle feedback simulated bythe NUIST Earth System Model: NESM-2.0.1

Abstract. In this study, we evaluate the performance of Nanjing University of Information Science & Technology Earth System Model, version 2.0.1 (hereafter NESM-2.0.1). We focus on model simulated historical and future oceanic CO2 uptake, and analyze the effect of global warming on model-simulated oceanic CO2 uptake. Compared with available observations and data-based estimates, NESM-2.0.1 reproduces reasonably well large-scale ocean carbon-related fields, including nutrients (phosphate, nitrite and silicate), chlorophyll, and net primary production. However, some noticeable discrepancies between model simulations and observations are found in the deep ocean and coastal regions. Model-simulated current-day oceanic CO2 uptake compares well with data-based estimate. From pre-industrial time to 2011, modeled cumulative CO2 uptake is 144 PgC, compared with data-based estimates of 155 ± 30 PgC. Diagnosed from the end of the benchmark 1 % per year CO2 increase simulations, carbon-climate feedback parameter, which represents the sensitivity of ocean CO2 uptake to climate change, is −7.1 PgC/K; Carbon-concentration feedback parameter, which represents the sensitivity of ocean CO2 uptake to increase in atmospheric CO2 is 0.81 PgC/ppm. These two feedback parameters diagnosed from model simulations are consistent with the mean value diagnosed from the CMIP5 (Coupled Model Intercomparison Project phase 5) model simulations under the same 1 % per year CO2 simulations (−7.8 PgC/K and 0.80 PgC/ppm, respectively). Our results demonstrate that NESM-2.0.1 can be used as a useful tool in the investigation of feedback interactions between the ocean carbon cycle, atmospheric CO2, and climate change.


Introduction
Global carbon cycle plays an important role in the climate system.The increase in atmospheric carbon dioxide (CO 2 ) is responsible for a large part of observed increase in global mean surface temperature (Ciais et al., 2013).From 1750 to 2016, about 645±80 PgC (1 PgC =10 15 gram carbon) of anthropogenic carbon has been emitted to the atmosphere, including 420±20 PgC from fossil fuels and industry and 225±75 PgC from land-use-change (Le Quéré et al., 2017).This CO 2 emission causes atmospheric CO 2 concentration to increase by 45% from an annual meanpreindustrial value of ~277 parts per million (ppm) Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-68Manuscript under review for journal Geosci.Model Dev. Discussion started: 2 May 2018 c Author(s) 2018.CC BY 4.0 License.(Joos and Spahni, 2008) to 404 ppm in 2016 (NOAA ESRL Global Monitoring Division, 2016).
The ocean carbon cycle is one of the main components of global carbon cycle.As a large carbon reservoir, the global ocean contains more than 50 times of carbon than that of the atmosphere (Menon et al., 2007).The global ocean also plays a key role in anthropogenic CO 2 uptake (Ballantyne et al., 2012;Wanninkhof et al., 2013).Since preindustrial time to now, about 25% of anthropogenic CO 2 (about 160±20 PgC) has been absorbed by the ocean (Le Quéré et al., 2017).The ocean carbon cycle involves many complex physical, chemical, and biological processes such as chemical reaction of water molecular with dissolved CO 2 , biological uptake of CO 2 through photosynthesis by phytoplankton in the upper ocean, and transport of inorganic and organic carbon into the ocean interior (Sarmiento and Gruber, 2006).
Increase in atmospheric CO 2 , by perturbing atmospheric radiation balance, leads to climate change.Changes in atmospheric temperature, precipitation, evaporation, and wind, induces changes in ocean physical properties such as temperature, salinity, and ocean circulation.(Levitus et al., 2000;Gregory et al. 2005;Pierce et al., 2012).These changes in ocean physical properties in turn affect the ocean carbon cycle (Sarmiento and Gruber, 2006).For example, increasing temperature can directly decrease CO 2 solubility and result in the reduction of ocean CO 2 uptake (Najjar 1992;Teng et al., 1996).Changes in wind speed can also directly influence air-sea carbon exchange by changing gas transfer velocity (Wanninkhof 1992;Wanninkhof and Trinanes 2017).Previous studies also show that global warming would lead to a weakening of global thermohaline circulation and increase ocean stratification (Gregory et al. 2005), both of which result in a slower transport of carbon from the surface to the ocean interior, reducing ocean's uptake of anthropogenic CO 2 (Cox et al., 2000;Zickfeld et al., 2008).Therefore, it is important to gain a good understanding of the potential effect of global warming on the ocean carbon cycle.Friedlingstein et al., (2006) proposed that the response of oceanic uptake of atmospheric CO 2 can be represented by the linear sum of two components: 1) carbon-concentration feedback, which refers to the response of ocean CO 2 uptake to increasing atmospheric CO 2 alone; 2) carbon-climate feedback, which refers to the response of ocean CO 2 uptake to CO 2 -induced climate change alone.Adopting the above conceptual framework, Friedlingstein et al., (2006) analyzed the carbon-concentration and carbon-climate feedback from Coupled Climate-Carbon Cycle Model Intercomparison Project (C 4 MIP) results.Since the work of Friedlingstein et al., (2006), a number of studies have analyzed the magnitudes of the carbon-concentration and carbonclimate feedback and their interactions under different CO 2 emission and concentration scenarios using CMIP5 model results (Gregory et al., 2009;Boer and Arora, 2009;Arora et al., 2013).
Given the importance of carbon cycle feedback in current and future global climate change, it is necessary to include the representation of global carbon cycle in Earth system models (Bretherton, 1985;Menon et al., 2007).In 2014, a new Earth system model, Nanjing University of Information Science and Technology Earth System Model version 1.0, was developed.NESM v1.0 reproduces realistic large-scale present-day climatic fields such as sea surface temperature (SST) and precipitation.
Large-scale climate variability and climate mode such as El Niño-Southern Oscillation (ENSO) and Madden-Julian oscillation (MJO) are also well reproduced (Cao et al., 2015).Moreover, NESM v1.0 also shows good skills in simulation of Indian monsoon precipitation (Li et al., 2017).Recently, the new version of NESM-2.0.1 is developed based on NESM v2.0.In this new version, Pelagic Interactions Scheme for Carbon and Ecosystem Studies (PISCES v2) is coupled to the ocean circulation model to represent the ocean biogeochemical processes (Aumont et al., 2015).PISCES can be used for both regional and global simulation of lower tropic levels of marine ecosystems and ocean carbon cycle (Bopp et al., 2005;Resplandy et al., 2012) The Earth System model of the Institute Pierre et Simon Lappace and Centre National de Recherche en Météorologie (IPPSL-CNRM) ESM, which uses PISCES to represent the ocean carbon cycle, contributed to CMIP5 (Séférian et al., 2013).
In this study, we present the performance of the ocean carbon cycle component in NESM-2.0.1.In Section 2, we describe the model of NESM-2.0.1 with the focus on the ocean carbon cycle component, as well as the setup of simulation experiments.In Section 3, we present the model performance.We evaluate modeled nutrients and ecosystem fields against available observations in Section 3.1.In Section 3.2, we evaluate modeled oceanic uptake of anthropogenic CO 2 during the historical period against data-based estimates.In Section 3. System Model that is designed to study interactions between different components of the Earth climate system and its response to natural and anthropogenic forcing.NESM-2.0.1 is developed based on the framework of NESM v2 with active ocean biogeochemical cycle, and the physical components of NESM-2.0.1 is same as NESM v2.The physical components of the NESM are described in detailed in Cao et al., (2015).Here we briefly introduce its main features.NESM consists of three main component models, including European Centre Hamburg Atmospheric Model (ECHAM v5.3) (Roeckner et al., 2003),   (Madec and NEMO team, 2012) and Los Alamos sea-ice model version 4.1 (CICE v4.1) (Hunke and Lipsomb, 2010).The three component models are coupled by the Ocean-Atmosphere-Sea-Ice-Soil (OASIS v3.0) Model coupling Toolkit (OASIS3-MCT) (Larson et al. 2005).The atmospheric resolution used in NESM v2 is T42L31 which has a horizontal resolution of ~ 2.8° latitude by 2.8° longitude and 31 vertical levels.The land surface model is a simple surface scheme that is implicitly coupled with the atmosphere, in which surface fluxes and temperature are computed using energy balance equation (Schulz et al., 2001).Surface albedo depends on background, snow, forest and canopy, and snow processes from canopy to soil are prognostically calculated (Roesch et al., 2001).A simple mixed-layer lake is also used.Ocean component runs with an ORCA2 global ocean configuration, which is a type of tripole grid.It is based on a 2 degree Mercator mesh and has 31 vertical levels with the thickness of ocean layer increasing from 10m at the upper ocean to 500m at 5000m depth.A local transformation is applied in the Tropics to refine the resolution up to 0.5 degree at the Equator.Sea-ice component includes four vertical ice layers and one snow layer with a multi- layer thermodynamic scheme.

Ocean biogeochemical component
NESM-2.0.1 employs standard PISCES v2 code to represent the ocean biogeochemical cycle in the atmosphere-ocean-sea-ice modeling system.No modification is applied to biological parameters and their formulations in the source code.There are 24 prognostic tracers in total in NESM-2.0.1, including dissolved inorganic and organic carbon, alkalinity, chlorophyll, nutrients and etc.Total Variance Dissipation (TVD) scheme (Zalesak, 1979) is used to describe the advection of these tracers.Carbonate chemistry including air-sea CO 2 exchange is formulated based on the Ocean Carbon-Cycle Model Intercomparison Project (OCMIP-2) protocol (more information can be access for free at http://ocmip5.ipsl.jussieu.fr/OCMIP/).
The inorganic and organic carbon cycle is simulated via the processes between nutrients, phytoplankton, and zooplankton.
These processes including phytoplankton growth and mortality, grazing by zooplankton, aggregation and remineralization of organic matter (Aumont et al., 2015).Growth rate of phytoplankton is limited by availability of nutrients, including phosphate, nitrate, silicate, iron, and ammonium.Also, phosphate, nitrate and carbon in phytoplankton are linked by a constant Redfield ratio.In NESM-2.0.1, the Redfield ratio of C: N: P is set to be 122:16:1 (Takahashi et al., 1985) and the -O/C ratio is set to 1.34 (Kortzinger et al., 2001).In contrast, the Fe / C, Chlorophyll / C, and Silicon / C ratio are prognostically simulated by the model.Mortality and aggregation of phytoplankton is associated with the calcifying organisms and bigenic silica.In NESM-2.0.1, 50% of the calcified organisms are associated with the shell.As a consequence, 50% of the dying calcified biomass is routed to the fast sinking particulate organic matter (POM).The degradation of the particles depends on the local temperature and while the remineralization of depends on local oxygen concentration.For POM, a simple two-compartment scheme is used.In that case, POM is modeled using two tracers corresponding to the two size classes: smaller class (1-100μm) and larger class (100-500μm), which also differ from their sinking speed.
The exchange between sediments and the ocean interior at the bottom boundary can be represented either with or without a sediment model.In NESM-2.0.1, the sediment model is inactive, basic but different parameterizations are applied depending on tracers considered.For example, buried biogenic silica is assumed to balance the external source, buried POM is determined by the organic carbon sinking rate at the bottom, and all the particulate iron would buried into the sediment once they reach the ocean bottom.The amount of unburied parts dissolve back to ocean instantaneously.Expect the interior processes, the external sources of nutrients are also described as geochemical boundary conditions in NESM-2.0.1, including atmospheric dust deposition of iron and silicon, river recharge of nutrients, dissolved carbon, and alkalinity, atmospheric deposition of nitrogen, and sediment mobilization of sedimentary iron.
In addition, light penetration parameterization is regarded as an important factor which could influence ocean thermal structure.
In NESM-2.0.1, incoming solar radiation is distributed among the upper ocean layers (391m), and bio-model penetration parameterization scheme is used to calculate the distribution of solar radiation.In that case, solar radiation penetrate rate is dependent on modeled chlorophyll concentration in each ocean layer (Lengaigne et al., 2009).

Simulations
First, NESM-2.0.1 was spun up for 1500 years with all related parameters set to pre-industrial values including orbit parameters and greenhouse gases (GHGs) concentrations (280 ppm for CO 2 , 720 ppb for CH 4 , 270 ppb for N 2 O, and 0 ppt for both CFC 11 and CFC 12 ).Averaged over the last 100 years of the spun-up simulation, globally integrated air-sea CO 2 flux is 0.03 PgC yr -1 , and the linear drift is 0.0006 PgC yr -1 per year, indicating that a quasi-equilibrium state has been reached for the global ocean carbon cycle.Meanwhile, global mean SST averaged over the last 100 years of spun-up simulation is 13.1 Celsius (℃) with the linear drift of 0.001 ℃ per year, and ocean mean temperature 3.4 ℃ with the linear drift of 0.0004 ℃ per year, indicating that ocean climate field has also reached a quasi-equilibrium state.
Using the end of the 1500-year spun simulation as the initial state of the nominal pre-industrial year of 1800, the model is further integrated to year 2100 with prescribed time-series of atmospheric concentrations of GHGs, including CO To separate the effect of atmospheric CO 2 and global warming on the ocean carbon cycle, for each of the three simulations above, three sensitivity experiments: biogeochemically coupled, radiatively coupled, and fully coupled were performed as summarized in Table 1.These types of simulations were also performed by previous studies that investigated the effect of CO 2 and CO 2 -induced climate change on the global carbon cycle (Friedlingstein et al., 2006;Arora et al., 2013;Schwinger et al., 2014).
1) Biogeochemically coupled (BC) simulations in which the code of the ocean carbon cycle sees changing atmospheric CO 2 , but the code of atmospheric radiation sees constant pre-industrial concentration of CO 2 and other GHGs as set in the spun-up simulation.In this way, the ocean carbon cycle is only affected by changing atmospheric CO 2 , but not GHG-induced climate change; 2) Radiatively coupled (RC) simulations in which the code of the ocean carbon cycle sees pre-industrial atmospheric CO 2 , but the code of atmospheric radiation sees changing concentration of atmospheric CO 2 and other GHGs.In this way, the ocean carbon cycle is only affected by GHG-induced climate change, but not the direct effect of changing atmospheric CO 2 .
3) Fully-coupled (FC) simulations in which both the codes of the ocean carbon cycle and atmospheric radiation see changing concentrations of atmospheric CO 2 and other GHGs.In this way, the ocean carbon cycle is affected by changes in both

Evaluation of NESM-2.0.1 simulated present-day ocean ecosystem
In this section, we compare model-simulated ocean biogeochemical fields that are directly related to the ocean ecosystem and carbon cycle, including nutrients, chlorophyll, and net primary production (NPP) against available observational and databased estimates.A detailed description of the observational data used is given in Appendix A. To have a direct comparison between NESM-2.0.1 results and observations, we interpolate NESM-2.0.1 results onto the corresponding grids of observational data using the distance-weighted averaged remapping method, as detailed in Appendix A.
Nutrients limit the growth of phytoplankton and play a vital role in the ocean biogeochemical cycles.Figure 1 compares modelsimulated annual mean spatial distributions of macronutrients, including nitrate (N), phosphate (P), and silicate (Si) during 1990s against WOA09 observations during the same period (Garcia, et al., 2010).As shown in Fig. 1, the model reproduces reasonably well the large-scale pattern of macronutrient distributions at ocean surface.Relatively high nutrient concentrations are found in the Southern Ocean, subarctic Pacific Ocean, and the mid-east Pacific Ocean where strong vertical mixing and upwelling bring nutrient-rich deep water to the surface (Whitney, 2011).Relatively low concentrations of nutrients are found in subtropical regions where the vertical mixing between surface and the deep ocean is relatively weak.Some noticeable discrepancies between model results and observations are noticed (Fig. 1).For example, modeled surface silicate concentrations are much higher than the observation in the Equatorial Pacific Ocean.Also, simulated surface phosphate concentration is lower than the observation in the Equatorial Pacific.
Figure 2 shows globally averaged vertical profiles of macro-nutrients during 1990s and WOA09 observations.In general, the concentrations of silicate, phosphate, and nitrate are more enriched in the deep ocean as a result of remineralization of organic matter in the ocean interior.NESM-2.0.1 simulates the vertical distribution of silicate reasonably well.However, the model underestimates the concentration of phosphate and nitrate between500 to 2000m depth, and overestimates their concentrations below ~2000m.The modeled-deficiencies of nutrient distributions is associated with the modeled deficiencies in the ocean dynamics or/and parameterizations of the ocean biological processes.The inclusion of natural radiocarbon ( 14 C) in the model, which is not implemented yet, would be useful in separating the effect of modeled ocean dynamics and biology (Stocker and Wright, 1996).We will further discuss this issue in the Section 4.
Figure 3 shows modeled spatial distribution of annual mean surface chlorophyll concentration during 1990s compared with SeaWifs observational data (Behrenfeld andFalkowski, 1997a, 1997b.).The model simulates reasonably well the large scale pattern of ocean surface chlorophyll concentration with high levels of chlorophyll in subarctic Pacific Ocean, North Atlantic Ocean, equatorial Pacific, and low levels of chlorophyll in subtropical ocean.NESM-2.0.1 simulates relatively higher chlorophyll concentration along the extratropical (except Arctic) coastal regions, but compared to observations, the model generally underestimates chlorophyll concentration at the tropical coastal regions, especially over the tropical Indian Ocean and Atlantic Ocean.This underestimation is probably associated with the deficiencies in modeled coastal dynamics, which is usually not represented well by the relatively coarse resolution global ocean models (Aumont et al., 2015).It is reported that the observed chlorophyll distribution is better reproduced when PISCES is coupled to a higher resolution ocean circulation model (Lee et al., 2000;Hood et al., 2003;Kone et al., 2009).In the Southern Ocean where the seawater is typically characterized by high nutrients and low chlorophyll (Lin et al., 2016), the model apparently overestimates chlorophyll concentration when compared with satellite-derived observations.Previous studies pointed out that chlorophyll concentrations derived from reflectance by standard algorithms tend to be underestimated by a factor of about 2 to 2.5, especially in intermediate concentrations regions such as the Southern Ocean (Garcia et al., 2005;Kahru and Mitchell, 2010).Therefore, the overestimate of chlorophyll concentration by NESM-2.0.1 in the Southern Ocean may partly be explained by the underestimates of satellite-derived chlorophyll.
Net primary production (NPP) of phytoplankton in NESM-2.0.1 is calculated as a function of nutrient availability, ocean temperature, and photosynthetic available radiation (PAR), which is associated with distribution of chlorophyll.Here we compare modeled NPP with data-based estimate that is calculated as a function of sea surface temperature, PAR, and satelliteobserved chlorophyll concentrations (Fig. 4).probably related to the model's bias in simulation of the Pacific cold tongue, which was shifted westward with overestimated intensity (Cao et al., 2015).NPP are underestimated in the northern North Atlantic, Arabian Sea, subarctic North Pacific and Arctic coastal regions, which is likely related to the model's underestimation of veridical mixing and coastal upwelling.
Averaged over 1990s, globally integrated ocean NPP from NESM-2.0.1 simulation is 41.5 PgC yr -1 , compared with data-based estimate of 37 to 67 PgC yr -1 .The large range of data-based estimate of global NPP is a result of different satellite observations and different algorithms for the NPP estimation (Longhurst et al., 1995;Antoine et al., 1996;Behrenfeld and Falkowski, 1997b;Behrenfeld et al., 2005).Global NPP simulated by CMIP5 models also show a wide range of values from 30.9 to 78.7 PgC yr - 1 (Bopp et al., 2013).NESM-2.0.1 simulated values of global NPP is therefore within the range of data-based estimates and current model simulations.Of the NESM-2.0.1 simulated global ocean NPP, 19% is contributed from diatoms, and 81% is contributed from nanoplankton.For comparison, the data-based estimate of NPP associated with diatoms accounts for 7% to 32% of the total NPP (Uitz et al.,2010;Hirata et al., 2011), and ocean biogeochemical models estimates that 15% to 30% global NPP is from diatoms (Aumont et al., 2003;Dutkiewicz et al., 2005;Yool and Popova, 2011).
Spatial resolution of the oceanic component of NESM-2.0.1 is relatively coarse, especially in high latitudes.It is reported that an ocean model with higher spatial resolution would produce a larger NPP as mesoscale and submesoscale processes would significantly enhance ocean biological productivity (McGillicuddy et al., 1998;Oschlies and Garçon, 1998;Lévy et al., 2001).Also, coastal regions would also be better presented in a higher resolution model.Therefore, it is possible that an ocean model with higher spatial resolution would simulate a higher NPP than what is presented here.

Evaluation of NESM-2.0.1 simulated ocean carbon cycle with observations
Here we compare NESM-2.0.1 simulated air-sea CO 2 flux and oceanic uptake of atmospheric CO 2 during the historical period against available observations.
We compare model-simulated air-sea CO 2 flux against observational-based estimate for the reference year of 2000 (Takahashi et al., 2009).Air-sea CO 2 flux in NESM-2.0.1 is calculated following OCMIP protocol that is defined as: Where   is the CO 2 gas transfer coefficient calculated as a function of wind speed,   is the solubility of CO 2 at the sea surface,   and   are CO 2 partial pressure at the overlying atmosphere and sea surface, respectively.Defined in this way, a positive value of F represents CO 2 flux out of the ocean, i.e. from ocean to the atmosphere.
As shown in Fig. 5, NESM-2.0.1 reproduces well observed large-scale pattern of air-sea CO 2 flux with CO 2 outgassing in the equatorial oceans and CO 2 uptake in the mid-to-high latitude oceans.Strong CO 2 uptake is observed in the North Atlantic Ocean where sea surface temperature is low and formation of deep water is active, both of which act to absorb CO 2 from the atmosphere.It appears that the model overestimates both the amount and extent of CO 2 outgassing in the equatorial Pacific Ocean.According to Eq. ( 1), air-sea CO 2 flux is calculated from the combined effects of three factors: surface wind speed, CO 2 solubility and partial pressure difference between atmosphere and ocean (dpCO 2 ).As expected, the global pattern of airsea CO 2 flux is mainly determined by simulated difference in CO 2 pressure (dpCO 2 ).Simulated bias in wind speed would also contribute to the bias in the air-sea CO 2 flux.As pointed out by Cao et al., (2015), wind speed in the equatorial Pacific Ocean is overestimated relative to observations, which would also contribute to the overestimated air-sea CO 2 flux in this region.
Now we compare model-simulated ocean storage of anthropogenic CO 2 during the 1990s (i.e.cumulative ocean CO 2 uptake from pre-industrial time as a result of ocean's uptake of anthropogenic CO 2 emissions) with data-based estimates from GLODAP (Key et al., 2004) in terms of both latitude-depth distribution (Fig. 6a, c) and vertically integrated column inventory (Fig. 6b, d).NESM-2.0.1 captures the observed distribution of anthropogenic CO 2 in the ocean.As pointed out by Sabine et al., (2004), deep penetrations of anthropogenic CO 2 are typically associated with convergence zones at temperate latitudes and high latitude oceans where vertical mixing is strong.For both data-based estimates and model simulations, substantial amount of anthropogenic CO 2 has penetrated down to the ocean interior as deep as 1000 m with two penetration tongues near 40°N and 40°S.In the North Atlantic where deep water forms, substantial amount of anthropogenic CO 2 has penetrated down to more than 2000 m of the ocean.NESM-2.0.1 also captures the observed large-scale pattern of vertically integrated column inventory of anthropogenic CO 2 with largest storage in the North Atlantic Ocean.However, the model appears to overestimate the carbon storage in the Pacific Ocean at the northern hemisphere.
Globally integrated anthropogenic CO 2 from NESM-2.0.1 simulation and observation for the 1990s is 120 PgC and 106 ± 17 PgC (Sabine et al., 2004), respectively.We also calculated the present-day anthropogenic CO 2 budget over different periods (1980s, 1990s, 2000s, 2002-2011, and from pre-industrial to 2011) and compared NESM-2.0.1 simulated results against the data-based estimate provided by IPCC AR5 (Table 2).The model-simulated ocean uptake of anthropogenic CO 2 is slightly lower than that from IPCC AR5, but within the estimated uncertainty range.For example, from the pre-industrial time to year 2011, NESM-2.0.1 simulated cumulative oceanic CO 2 uptake is 144 PgC, compared with data-based estimates of 155 ± 30 PgC.
Figure 7 compares spatial pattern of NESM-2.0.1 simulated carbon-related fields with corresponding observations using Taylor diagrams (Taylor, 2001).As shown in Fig. 7, model-simulated statistic pattern of surface nitrite and phosphate compares well with observations with correlation coefficients r>0.9 and normalized standardized deviation (SD) close to 1.0.However, simulated spatial pattern of surface silicate shows larger deviations from the observation.Simulated spatial pattern of chlorophyll and NPP compares poorly with observations with a correlation of 0.45 and 0.40, respectively.Simulated spatial pattern of air-sea CO 2 flux shows good agreement with observations with a correlation coefficient of 0.75 and a normalized SD close to 1.It is noted that chlorophyll, NPP and air-sea CO 2 flux is not directly observed but diagnosed from observationalbased data, and there is considerable uncertainty associated with the data.
In short, NESM-2.0.1 reproduces reasonably well the large-scale features of marine biogeochemical fields, including nutrients, chlorophyll, NPP, air-sea CO 2 flux, and anthropogenic carbon storage, indicating that the model can be used as a useful tool to investigate the response of the ocean carbon cycle to changes in atmospheric CO 2 and global climate.

Response of the ocean CO2 uptake to atmospheric CO2 and global warming
In this section, we first present NESM-2.0.1 simulated physical climate change and oceanic CO 2 uptake under the prescribed atmospheric CO 2 concentration pathway of RCP 8.5 scenario.We then present NESM-2.0.1 simulated oceanic CO 2 uptake in the benchmark simulations with 1% per year increase in atmospheric CO 2 .

NESM-2.0.1 simulated physical climate change under RCP 8.5
In this section, we investigate NESM-2.0.1 simulated oceanic CO 2 uptake in response to increasing atmospheric CO 2 concentration and CO 2 -induced climate change.Increasing atmospheric CO 2 directly affects air-sea CO 2 flux and thus oceanic CO 2 uptake.Meanwhile, CO 2 -induced warming also affects the ocean carbon cycle via changes in climate fields such as temperature, ocean stratification, and ocean circulation.Previous studies proposed that the response of the ocean carbon uptake can be approximated as a linear sum of two parts: 1) the response to changes in atmospheric CO 2 concentration (carbonconcentration feedback), and 2) the response to changes in surface temperature (carbon-climate feedback) (Friedlingstein et al., 2006;Arora et al., 2013).
To investigate carbon-concentration feedback and carbon-climate feedback and the nonlinear interactions between the two, we use the NESM-2.0.1 to perform three sensitivity simulations: biogeochemically-coupled (BC), radiatively-coupled (RC) and fully-coupled (FC).The detailed setup of these experiments were described in Section 2 and listed in Table 1. is within the range of CMIP5 model results of 3.7±0.7K(Collins and Knutti, 2013;Knutti and Sedláček, 2013).Associated with increasing atmospheric temperature, global ocean also becomes warmer in FC and RC simulations, reducing CO 2 solubility in seawater and acting to reduce oceanic CO 2 uptake.
MLD and AMOC show much stronger interannual fluctuation than SAT, but both of them show long-term trend of decrease.
The reduction of mixed layer depth, which is associated with a relatively faster warming of the surface ocean and a slower response of the deep ocean, indicates a more stabilized upper ocean with global warming (Held et al., 2010).The pre-industrial value of the NESM-2.0.1 simulated AMOC index is 24 Sv (1Sv =10 6 m 3 s -1 ), compared with the value of 14 to 31 Sv from CMIP5 models (Weaver et al., 2012).The modeled annual mean of AMOC transport at 26°N averaged from 2000 to 2004 is 25 Sv while the observation record from RAPID/MOCHA (Rapid Climate Change programme / Meridional Ocean Circulation and Heatflux Array) is 17.5 ± 3.8 Sv during 2004 to 2011 (Rhein et al., 2013).A substantial weakening of AMOC intensity in the RC and FC simulations is observed in 21th century under the RCP8.5 scenario, which is associated with ocean surface warming and increased freshwater input at the North Atlantic (Gregory et al., 2005).By 2100, the simulated intensity of AMOC declines to about half of its pre-industrial value.The simulated 54% weakening of AMOC by the end of this century is at the higher end of what are simulated by CMIP5 models that range from 15% to 60% under the RCP 8.5 scenario (Cheng et al., 2013).

NESM-2.0.1 simulated oceanic CO2 uptake under RCP 8.5
Changes in the ocean carbon cycle are tightly associated with the change in physical climate processes (Doney et al., 2004).
In the FC simulation, weakening of ocean vertical mixing, as indicated by the reduced mixed layer depth and weakening of AMOC, will reduce the vertical exchange of CO 2 between the upper ocean and the ocean interior, and thus suppress oceanic CO 2 uptake.A warmer surface ocean would reduce CO 2 solubility, also suppressing oceanic CO 2 uptake.
Figure 9 shows time evolution of model-simulated oceanic CO 2 uptake for the simulations of BC (biogeochemical coupled), RC (radiative coupled), FC (fully coupled), and the linear sum of BC and RC.In BC simulation, only increasing atmospheric CO 2 affects the ocean carbon cycle.By year 2100, the modeled global ocean has absorbed a total of 604 PgC anthropogenic CO 2 from the atmosphere.In RC simulation, constant atmospheric CO 2 is seen by the ocean carbon cycle, and the atmosphere radiation sees increasing atmospheric CO 2 concentration.As discussed above, increasing sea surface temperature, increasing ocean stratification, and reduced AMOC all act to decrease CO 2 uptake.As a result, CO 2 -induced warming alone causes the ocean to release CO 2 to the atmosphere (negative CO 2 uptake).By year 2100, modeled cumulative CO 2 uptake is -37.6 PgC.
In FC simulation, oceanic CO 2 uptake is affected by both of the increase in atmospheric CO 2 and CO 2 -induced global warming.
By the end of 21th century, the NESM-2.0.1 simulated cumulative oceanic CO 2 uptake in the FC run is 516 PgC, which is close to the middle value of ~500 PgC from CMIP5 models under the same RCP 8.5 scenario (Jones et al., 2013).
As seen from Figure 9, the sum of the simulated oceanic CO 2 uptake from the BC and RC simulations (566 PgC) is larger than that from the FC run (516 PgC), indicating that the biogeochemical effect of atmospheric CO 2 (carbon-concentration feedback) and radiative effect of atmospheric CO 2 (carbon-climate feedback) on oceanic CO 2 uptake is not exactly additive.This nonlinearity has also been found in previous studies (Boer and Arora, 2009;Gregory et al., 2009;Schwinger et al., 2014).The NESM-2.0.1 simulated nonlinearity (discrepancy between the sum of biogeochemical and radiative effect on ocean carbon uptake and the total carbon uptake, i.e.BC+RC-FC) is 50.4PgC by the end of 21th century.This nonlinearity is about 9.8% of the total ocean uptake, and the magnitude of the nonlinearity (50.4 PgC) is comparable to the magnitude of the radiative effect on ocean carbon uptake (-37.6 PgC).
To better understand oceanic CO 2 uptake in response to changing atmospheric CO 2 and CO 2 -induced climate change, Figure 10 shows the spatial distribution of anthropogenic air-sea CO 2 flux at the end of 21th century (averaged over year 2091 to 2100) under the RCP8.5 scenario for FC, RC, and BC simulations, respectively.Also shown is the difference between FC and the sum of RC and BC.Positive values represent CO 2 flux out of ocean, and negative values represent CO 2 flux into the ocean.
In the BC simulation, the ocean absorbs atmospheric CO 2 in most regions except for a few scattered grid points of Pacific Ocean at middle latitude with slightly CO 2 outgassing.The strongest CO 2 uptake is observed in the North Atlantic and Southern Ocean.Results from RC simulation show CO 2 outgassing in large parts of the global ocean as a result of CO 2 -induced warming that reduces CO 2 solubility and the rate of CO 2 transport from the surface to deep ocean.In the Arctic Ocean, warming induces a net uptake of CO 2 as a result of reduced sea ice extend under global warming, which allows more open seawater to absorb atmospheric CO 2 .FC simulation show the combined effect of increasing atmospheric CO 2 and CO 2 -induced warming on oceanic CO 2 uptake (Fig 10c).Positive oceanic CO 2 uptake is observed in most regions, indicating the dominant role of CO 2 biogeochemical effect in the total oceanic carbon uptake.Similar to the BC simulation, strongest CO 2 uptake is observed in the North Atlantic and Southern Ocean.Somewhat enhanced outgassing is observed in part of subtropical Pacific Ocean where the mixed layer depth is shallow and the water becomes more stratified.Oceanic CO 2 uptake in the subtropical Pacific Ocean is small in BC simulation, and it appears that the radiative effect of atmospheric CO 2 is dominant in these regions.
Figure 10d shows the spatial distribution of the difference in air-sea CO 2 flux between FC simulation and the sum of BC and RC simulation.The difference represents the nonlinearity between carbon-climate feedback and carbon-concentration feedback.In NESM, relatively large nonlinearity is observed in the northern North Atlantic Ocean and Southern Ocean (especially the southern South Atlantic), which is consistent with the finding of previous studies (Zickfeld et al., 2011;Schwinger et al., 2014).The nonlinearity can be explained by the interaction between elevated CO 2 concentration and climate change.For example, compared to pre-industrial CO 2 background, under elevated CO 2 concentrations (RCP 8.5 here), reduced North Atlantic Deep Water (NADW) has a much larger impact on reducing CO 2 transport from the surface to deep ocean.
Therefore, anthropogenic CO 2 uptake in the North Atlantic in the fully coupled simulation (FC) is smaller than the linear sum of the biogeochemically-coupled (BC) and radiatively-coupled (RC) simulations, Also, in the FC simulation, compared with the RC simulation, more ocean carbon is subject to increasing ocean temperature and therefore ocean warming would have a larger effect on CO 2 outgassing compared to the warming effect on top of pre-industrial CO 2 state.Thus, the sum of RC and BC simulation underestimates CO 2 outgassing from the ocean in the FC simulation.

carbon-concentration and carbon-climate feedback diagnosed from NESM-2.0.1
We further investigate the carbon-concentration and carbon-climate feedback.Friedlingstein (2006) proposed that cumulative oceanic CO 2 uptake can be decomposed approximately using the linear sum of carbon-concentration feedback and carbon-climate feedback: Where ′ is the air-sea CO 2 flux change and ∫ ′ biogeochemically-coupled simulations and radiatively-coupled simulations, and diagnosed these two coefficients using CMIP5 model outputs.In the biogeochemically-coupled simulations where the ocean carbon uptake is only affected by changing atmospheric CO 2 , Eq. ( 2) can be simplified as: where F' represent air-sea CO 2 flux change in the biogeochemically coupled simulation.In the radiatively-coupled simulations where the ocean carbon uptake is only affected by changing climate, Eq. ( 2) can be simplified as: where F' represent air-sea CO 2 flux change in the radiatively coupled simulation.
In this study, we estimate the strength of carbon-concentration and carbon-climate feedback following Arora et al., (2013) using equations ( 3) and (4).Figure 11 shows the change in ocean carbon storage as a function of change in and atmospheric CO 2 (Fig. 11a) and global annual mean surface temperature (Fig. 11b), respectively.The derived evolution of carbonconcentration feedback parameter  as a function of atmospheric CO 2 concentration and carbon-climate feedback parameter  as a function of change in temperature is shown in Fig. 11c and 11d, respectively.
As shown in Fig. 11, in the BC simulation modeled ocean storage of anthropogenic CO 2 scales roughly linearly with atmospheric CO 2 , and in the RC simulation modeled ocean storage of anthropogenic CO 2 scales roughly linearly with changes in global mean surface temperature.Increasing atmospheric CO 2 alone increases oceanic CO 2 uptake whereas increasing temperature alone decreases CO 2 uptake.Therefore, carbon-climate parameter  is negative while carbonconcentration parameter  is positive (Fig. 11).After an initial adjustment period, a decreasing trend of  is observed, which reaches to the value of -6.2 PgC per K at the end of 21th century.The decreasing trend of  (more negative) indicates that with enhanced warming, one degree of surface temperature increase would induce more CO 2 outgassing form the ocean (Fig 11d).The parameter  initially increases with atmospheric CO 2 , and then slightly decreases (Fig. 11c) with the value of 0.83 PgC/ppm at the year 2100.The decreasing trend of  is consistent with the slowdown of the increasing trend of CO 2 flux at the end of 21th century.Similar trends of carbon-climate feedback parameter and carbon-concentration feedback parameter are also found in the previous research (Arora et al., 2013).The increased sensitivity of CO 2 outgassing to increasing temperature and decreased sensitivity of CO 2 uptake to increasing atmospheric CO 2 indicate that for a given amount of CO 2 emission, a larger fraction of anthropogenic CO 2 would remain in the atmosphere.Arora et al., (2013) analyzed carbon-concentration and carbon-climate feedback parameters from CMIP5 models using the benchmark simulations in which atmospheric CO 2 is assumed to increase at a rate of 1% per year for 140 years to reach 4× CO 2 .To have a direct comparison with CMIP5 results, we performed another set of simulations under the same CO 2 concentration pathway.

carbon-concentration and carbon-climate feedback parameters from 1% per year CO2 simulations
Figure 12 shows time evolution of model-simulated oceanic CO 2 uptake for the simulations of BC, RC, FC, and the linear sum of BC and RC in 1% per year CO 2 increase experiments.To some extent, Figure 12 shows similar results with time evolution of model-simulated oceanic CO 2 uptake under RCP 8.5 scenario as shown in Figure 9.By year 140, largest positive oceanic CO 2 uptake (685 PgC) is observed in BC simulation while relative small negative oceanic uptake (-36 PgC) is seen in RC simulation, and the sum of the oceanic CO 2 uptake from the BC and RC simulations (649 PgC) is larger than that from the FC run (592 PgC).It is reported that the nonlinearity between biogeochemical and radiative feedback accounts for 3.6% -10.6% of the total ocean carbon uptake for CMIP5 models in 1% per year CO 2 increase experiments (Schwinger et al., 2014).For comparison, at the end of 1% per year CO 2 increase simulations, NESM-2.0.1 simulated nonlinearity is 9.6% (57 PgC) of the total oceanic CO 2 uptake, which is comparable with the relative magnitude of the nonlinearity estimated by CMIP5 models.
The comparison of feedback parameters is provided in Fig. 13.At the end of 1% increasing CO 2 simulation, the diagnosed value of β from CMIP5 models ranges from 0.69 to 0.91 PgC/ppm with a multi-model mean value of 0.

Conclusion and discussion
In this study, we evaluated the performance of the NUIST Earth System Model (NESM) in simulating the present-day ocean carbon cycle.We also investigated the model-simulated oceanic CO 2 uptake in terms of the oceanic CO 2 uptake in response to increasing atmospheric CO 2 alone, e.g., carbon-concentration feedback, and the oceanic CO 2 uptake in response to CO 2induced climate change alone, e.g., carbon-climate feedback.
The model simulates reasonably well the large-scale patterns of surface nutrient concentration including nitrate, phosphate, and silicate.The model also does a reasonable job in simulating large-scale distribution of chlorophyll and ocean net primary production (NPP).The global integrated ocean NPP simulated by NESM-2.0.1 is 42 PgC yr -1 , while observation-based estimate ranges from 37 to 67 PgC yr -1 and ocean NPP simulated from CMIP5 models ranges from 31 to 79 PgC yr -1 .The NESM-2.0.1 simulated cumulative anthropogenic CO 2 uptake from the pre-industrial time to year 2011 is 144 PgC, compared well with data-based estimates of 155 ± 30 PgC (Ciais et al., 2013).
As proposed by Friedlingstein et al. (2006), the response of oceanic CO 2 uptake can be decomposed by the sum of two components, the response to increasing atmospheric CO 2 alone (carbon-concentration feedback) and the response to CO 2induced climate change alone (carbon-climate feedback): In the simulation where atmospheric CO 2 increases by 1% per year, by year 140, the NESM-2.0.1 simulated carbon-concentration feedback parameter is 0.81 PgC/ppm and carbon-climate feedback parameter is -7.1 PgC/K, indicating that increasing atmospheric CO 2 alone increases oceanic CO 2 uptake while global warming alone decreases oceanic CO 2 uptake.These estimated feedback parameters are in general agreement with those estimated by CMIP5 models that yield carbon-concentration feedback parameters ranging from 0.69 to 0.91 PgC/ppm and carbon-climate feedback parameters ranging from -2.4 to -12.1 PgC/K.These results indicate that the NESM-2.0.1 simulated response of oceanic CO 2 uptake to increasing atmospheric CO 2 and climate change are in general agreement with the state-ofthe-art Earth system models.
It is also found from NESM-2.0.1 simulations that the sum of oceanic CO 2 uptake in response to changes in atmospheric CO 2 alone and climate changes alone is somewhat larger than the oceanic CO 2 uptake in response to the combined effect of atmospheric CO 2 concentration and CO 2 -induced climate change.The difference between the total oceanic CO 2 uptake and the linear sum of carbon-concentration feedback and carbon-climate feedback indicates the nonlinearity between these two feedbacks.In the NESM-2.0.1 simulation with 1% per year increase in atmospheric CO 2 , by year 140, the nonlinearity is 9.6% of the total oceanic CO 2 uptake.For comparison, nonlinearity between the carbon-concentration and carbon-climate feedback from CMIP5 models accounts for 3.6% to 10.6% of the total oceanic CO 2 uptake.While the results presented here show some success of the NESM-2.0.1 model in the simulation of processes involving the ocean ecosystem and ocean CO 2 uptake, the model simulation shows a number of caveats that requires better understanding of the source of the deficiencies and calls for future improvements.In particular, it is not well understood (a) why in the equatorial Pacific the simulated surface silicate concentrations are much higher than the observation, but the surface phosphate concentration is lower than the observation (Fig. 1), (b) why the model captures well the vertical distribution of silicate, but not the nitrate and phosphate concentrations (Fig. 2), and (c) why the model simulates chlorophyll concentration well along the extratropical coastal regions, but not along the tropical coastal regions (Fig. 3).The simulated annual mean distribution of vertical integrated NPP has notorious errors (Fig. 4).Both the simulated NPP and air-sea CO 2 flux have large biases in the equatorial Pacific (Fig. 4 and 5), which are likely associated with the model's bias in simulated SST and upwelling fields.
These simulated biases are likely to be associated with the climatological bias in simulation of the Pacific cold tongue.The simulated cold tongue shifted westward compared with observations (Cao et al., 2015), implying that the upwelling also shifted westward (Jin 1996;Li and Xie 2014), which may lead to deficiencies in nutrients distribution, NPP and air-sea CO2 flux.
This suggests that overcoming of the mean climatological cold tongue bias should be a major target in further improvement of the oceanic carbon cycle simulation.
In the future, the NESM, together with its ocean carbon cycle component, will be continuously developed and improved.To better evaluate the NESM-2.0.1 simulated ocean dynamics and the ocean carbon cycle, the simulation of natural and bomb 14 C will be implemented because their distributions in the ocean is a good indicator of the strength of ocean mixing and deep ocean circulation (Levin and Vago, 2000;Matsumoto 2007;Skinner et al, 2017).Also, the model appears to underestimate the concentration of chlorophyll and NPP in the tropical coastal regions as a result of relatively low model resolution, which cannot capture well coastal processes.The development of a higher-resolution (~10km) NESM-2.0.1 is planned in the near future, which could better capture the mesoscale processes and coastal dynamics (Griffies et al. 2015).In the future, the NESM-2.0.1 will be used to conduct a series of studies to better understand interactions between the ocean carbon cycle, atmospheric CO 2 , and climate change.
3, we analyze modeled carbon-concentration feedback and carbon-climate feedback under different CO 2 concentration scenarios and compare our results with CMIP5 model results.Conclusion and discussions are given in Section 4. Nanjing University of Information Science and Technology (NUIST) Earth System Model (NESM) is a comprehensive Earth Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-68Manuscript under review for journal Geosci.Model Dev. Discussion started: 2 May 2018 c Author(s) 2018.CC BY 4.0 License.
The large-scale distribution of relatively high NPP in the central-eastern Pacific, subarctic Pacific, and North Atlantic is well reproduced by NESM-2.0.1.Major discrepancies are seen in the central-eastern Pacific, northern North Atlantic and Arctic coastal regions.The overestimated NPP in the equatorial central-eastern Pacific is Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-68Manuscript under review for journal Geosci.Model Dev. Discussion started: 2 May 2018 c Author(s) 2018.CC BY 4.0 License.

Figure 8
Figure 8 shows NESM-2.0.1 simulated changes (relative to pre-industrial level) in global annual mean surface air temperature (SAT), mixed layer depth (MLD), and the intensity of Atlantic meridional overturning circulation (AMOC) at 30°N from 1900 to 2100 under RCP8.5 scenario.As expected, substantial change of SAT, MLD and AMOC intensity is observed in FC and RC simulation where increasing atmospheric CO 2 is allowed to affect atmospheric radiative transfer.No significant long-term trend is observed in BC simulation where increasing atmospheric CO 2 is only allowed to affect the ocean carbon cycle.NESM-2.0.1 simulated annual mean SAT anomalies over the period of 2080 to 2100 (relative the period of 1986-2005) is 4.0K, which Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-68Manuscript under review for journal Geosci.Model Dev. Discussion started: 2 May 2018 c Author(s) 2018.CC BY 4.0 License.

𝑡 0 𝑑𝑡
represent the cumulated ocean carbon uptake. represents the sensitivity of ocean carbon storage to climate change and ∆ is the change in global mean surface air temperature. represents the sensitivity of ocean carbon storage to the change in atmospheric CO 2 concentration and ∆  is the atmospheric CO 2 concentration change.Arora et al., (2013) obtained these two parameters following equation (2) with different sensitive experiments, i.e.
Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-68Manuscript under review for journal Geosci.Model Dev. Discussion started: 2 May 2018 c Author(s) 2018.CC BY 4.0 License.change.The mean value of carbon-climate feedback parameter from CMIP5 models is -7.8 PgC/K.For comparison, at the end of simulation, the NESM-2.0.1 simulated value of carbon-climate parameter  with the scenario of 1% per year increase in atmospheric CO 2 is -7.1 PgC/K.These results indicate that NESM-2.0.1 simulated response of oceanic CO 2 uptake to atmospheric CO 2 and climate change is in general agreement with ensemble means of the CMIP5 models.

Figure 5 :
Figure 5: NESM-2.0.1 simulated air-sea CO2 flux (g C m -2 yr -1 ) at year 2000 against observational data.Left panels are geographical distribution and right panel is zonal mean pattern.Positive values represent CO2 flux out of ocean, and negative values represent CO2 flux into the ocean

Figure 7 :
Figure 7: Taylor diagram comparing statistic pattern of annual mean carbon-related fields between NESM-2.0.1 simulations and corresponding observations, including surface nitrate, phosphate, silicate, alkalinity, chlorophyll concentration, verticallyintegrated net primary production, and air-sea CO2 flux.All fields are normalized by the standard deviation of corresponding

Figure 9 .
Figure 9.The NESM-2.0.1 simulated (a) annual mean oceanic CO2 uptake and (b) cumulative oceanic CO2 uptake for the simulations of BC (biogeochemical coupled), RC (radiative coupled), FC (fully coupled), and the linear sum of BC and RC under the RCP 8.5 scenario.

Figure 10 .
Figure 10.Spatial distribution of anthropogenic air-sea CO2 flux at the end of 21th century (averaged over year 2091 to 2100) under RCP8.5 scenario for BC (biogeochemically coupled), RC (radiatively coupled), and FC (fully coupled) simulations, respectively.Also shown is the difference between FC simulation and the sum of RC and BC simulations (FC-RC-BC).Positive values represent CO2 flux out of ocean, and negative values represent CO2 flux into the ocean.

Figure 11 :Figure 12 :
Figure 11: The NESM-2.0.1 simulated cumulated ocean uptake against (a) the atmospheric CO2 in the BC experiment and (b) the global mean surface air temperature change in the RC experiment.Also shown is time evolution of diagnosed carbonconcentration feedback parameter as a function of atmospheric CO2 (c) and carbon-climate feedback parameter as a function of global mean surface air temperature change (d).
PgC/K at the end of simulation.The larger model-dependence of carbon-climate feedback parameter is expected given the large uncertainty in model-simulated climate change and dependence of carbon cycle processes on climate