Extending the relationship between global warming and cumulative carbon emissions to multi-millennial timescales

To visualize the time-dependent multi-centennial climate response to carbon emissions, a 1000 yr simulation with the fully coupled carbon cycle–climate model GFDL ESM2M (Dunne et al 2012) has been performed. The physical state model underlying GFDL ESM2M is an updated version of the coupled model CM2.1 (Delworth et al 2006), consisting of a 1° version of the MOM4p1 ocean model (Griffies 2009) coupled to an approximately 2° atmospheric model (AM2, Anderson et al 2004). The ocean biogeochemical model is Tracers of Ocean Phytoplankton and Allometric Zooplankton code version 2 that includes the representation of the carbon cycle. The model also includes representation of the land biogeochemistry.

The GFDL ESM2M is forced with prescribed atmospheric CO₂ increasing at a rate 1% per year from preindustrial 286 to 745 ppm (year 99 of the simulation) when global annual mean surface temperature is 2 °C above preindustrial levels, after which carbon emissions stop. The diagnosed cumulative carbon emissions that are compatible with the prescribed atmospheric CO₂ concentration pathway by year 99 are 1947 Gt C (gray line in figure 1(a)). Here, the carbon emissions are stopped at 2 °C global warming, as this temperature target has been widely discussed and commonly regarded as an adequate means of avoiding dangerous climate change in policy making. An extension of the standard 1% CO₂ increase per year scenario has been chosen as this scenario is usually utilized to derive the TCRE (Matthews et al 2009). The CO₂ emissions stoppage is a rather unrealistic setup for temperature stabilization. However, intermediate complexity model simulations suggest that this scenario gives virtually identical long-term temperature changes as when using the same amount of total cumulative carbon emissions following more realistic scenarios (Eby et al 2009). Forcing from non-CO₂ greenhouse gases and other radiative forcing agents were kept at their preindustrial values throughout the simulations. Note that the GFDL ESM2M is the same model as used in the idealized carbon pulse experiment by Frölicher et al (2014). They used similar cumulative carbon emissions (1800 Gt C) and also suggest continued multi-centennial global warming after zero carbon emissions.

Zoom In Zoom Out Reset image size

To estimate the GFDL ESM2M equilibrium temperature response in year 10 000 (T_{10 000}), the ECS of the GFDL ESM2M is scaled with a factor of ln(CO_{2,10 000}/CO_2,PI)/ln(2), following the conventional form for CO₂ forcing:

$$\begin{eqnarray}&&{{\rm{T}}}_{\mathrm{10\; 000}}\left[^\circ {\rm{C}}\right]={\rm{ECS}}\left[^\circ {\rm{C}}\right]\displaystyle \frac{\mathrm{ln}\left(\displaystyle \frac{{{\rm{CO}}}_{\mathrm{2,10\; 000}}\left[{\rm{ppm}}\right]}{{{\rm{CO}}}_{2,{\rm{PI}}}\left[{\rm{ppm}}\right]}\right)}{\mathrm{ln}\;2}\end{eqnarray} \tag{ 1 }$$

The ECS of 3.3 °C is estimated from a 4000 yr 1% CO₂ to doubling experiment with the GFDL ESM2M model. After CO₂ emissions are ceased in year 99 it takes ESM2M about an additional 900 years to be in equilibrium with the contemporaneous CO₂ radiative forcing i.e. the ocean heat uptake is virtually zero (black line labeled with N in figure 1(c)). The atmospheric CO₂ in year 10 000, CO_{2,10 000} is obtained using

$$\begin{eqnarray}&&{{\rm{CO}}}_{\mathrm{2,10\; 000}}\left[{\rm{ppm}}\right]=\\ &&{{\rm{CO}}}_{2,{\rm{PI}}}\left[{\rm{ppm}}\right]+0.8\displaystyle \frac{{{\rm{emis}}}_{99}\left[{\rm{Gt}}\;{\rm{C}}\right]{{\rm{af}}}_{1000}}{2.12[{\rm{Gt}}\;{\rm{C}}\;{{\rm{ppm}}}^{-1}]}\end{eqnarray} \tag{ 2 }$$

with the unit conversion factor of 2.12 Gt C ppm⁻¹, the simulated cumulative carbon emissions in year 99, emis₉₉, of 1947 Gt C, and the preindustrial atmospheric CO₂ concentration CO_2,PI of 286 ppm. The GFDL ESM2M simulates an airborne fraction of anthropogenic carbon emissions, af₁₀₀₀, of 0.21 in year 1000. This af₁₀₀₀ is scaled by a factor of 0.8 under the assumption that 80% of the anthropogenic carbon emissions in year 1000 are still in the atmosphere in year 10 000 (Archer et al 2009).

We use existing model data from 12 ESMs that participated in CMIP5 (supplementary table 1). The temperature time series of the first 99 years (T₉₉) are obtained from available 1% increasing to 4xCO₂ (286–1144 ppm) experiments. The temperatures in year 1000, T₁₀₀₀, are estimated using:

$$\begin{eqnarray}&&{{\rm{T}}}_{1000}\left[^\circ {\rm{C}}\right]=\\ &&{\rm{ECS}}\left[^\circ {\rm{C}}\right]\displaystyle \frac{\mathrm{ln}\left(\displaystyle \frac{\displaystyle \frac{{{\rm{emis}}}_{99}\left[{\rm{Gt}}\;{\rm{C}}\right]{{\rm{af}}}_{1000}}{2.12[{\rm{Gt}}\;{\rm{C}}\;{{\rm{ppm}}}^{-1}]}+{{\rm{CO}}}_{2,{\rm{PI}}}\left[{\rm{ppm}}\right]}{{{\rm{CO}}}_{2,{\rm{PI}}}\left[{\rm{ppm}}\right]}\right)}{\mathrm{ln}\;2}\end{eqnarray} \tag{ 3 }$$

with cumulative carbon emissions in year 99, emis₉₉, from each individual model (forth row in supplementary table 1). Here, we assume that the models are in equilibrium with the contemporaneous CO₂ forcing after 1000 years, i.e. the ocean heat uptake is zero. The effective equilibrium climate sensitivities, ECS (second row in supplementary table 1), are obtained from Andrews et al (2015), who estimated the effective ECS with a linear extrapolation method using years 21–150 from CO₂ instantaneously quadrupling experiments. At present, the instantaneously quadrupling experiments are usually used to estimate the ECS, but these simulations are still not long enough to accurately estimate the ECS (Frölicher et al 2014, Paynter and Frölicher 2015). This can be exemplified with the GFDL ESM2M which when run to equilibrium has an ECS of 3.3 °C, but an estimated effective ECS of 2.6 °C when using the linear extrapolation method. However, this uncertainty is small compared to the spread in ECS across the models.

The airborne fraction of carbon emissions in year 1000, af₁₀₀₀ = 0.25, is obtained from an average of the eight EMIC simulations (sixth row in supplementary table 2) that were run for 1000 years under a 1800 Gt C delta-pulse experiment (Eby et al 2013). 1800 Gt C is close to the 1919 ± 154 Gt C emitted in the CMIP5–ESMs by year 99 (fourth row in supplementary table 1). It has been shown that the long-term atmospheric CO₂ response to a carbon pulse simulated by CMIP5–ESMs is within the uncertainty range simulated by EMICs (Joos et al 2013). To test the robustness of our results, we also use a lower (af_1000,min = 0.17) and upper bound (af_1000,max = 0.31) for the airborne fractions of anthropogenic CO₂ emissions, which is the range simulated across the eight EMICs considered in this study (sixth row in supplementary table 2). As the CMIP5–ESMs are not run over period 99–10 000, we assume that the temperature change between years 99 and 1000 and 10 000 is linear as a function of time. The resulting change between these years is best viewed as a trend rather than the actual pathway that the model will follow.

Global mean temperature change in year 10 000 is obtained from equations (1) and (2) using emis₉₉ and ECS of the individual models. Again, the temperature change between T₁₀₀₀ and T_{10 000} is estimated using linear interpolation between these two years.

We also use model data from eight EMIC models to estimate the temperature response in year 99, 1000 and 10 000 using equations (1)–(3). All values are given in supplementary table 2. T₉₉, emis₉₉, af₁₀₀₀ and ECS are obtained from Eby et al (2013). Here, af₁₀₀₀ and T₉₉ are both simulated by the EMIC models, but under slightly different scenarios (4xCO₂ pulse experiments for af₁₀₀₀ and 1% CO₂ increase for T₉₉). The uncertainty in the multi-millennial temperature response for the EMICs is overall smaller than for the CMIP5–ESMs because the full range across the eight EMICs is used for the CMIP5–ESMs. Note also that the ECS estimates of the EMICs might slightly be low-biased as the ECS were obtained from relatively short 1000 yr long 2xCO₂ simulation implying that the models were not yet in equilibrium, i.e. ocean heat uptake is not zero as can be seen under http://climate.uvic.ca/EMICAR5/results/4Xd/ for a similar scenario. The use of EMIC models allows us to investigate if there is a fundamental difference in the multi-centennial warming response between CMIP5–ESMs and EMICs.

The multi-millennial climate response to carbon emissions in GFDL ESM2M can be characterized by three phases (indicated with dashed vertical lines in figure 1(b)). During phase I (years 0–99), while carbon emissions steadily increase, global mean temperature (black line in figure 1(b)) increases nearly linearly until it reaches 2 °C above preindustrial (green point in figure 1(b)). At the end of phase I, the ratio between the simulated temperature (black line in figure 1(b)) and the temperature in thermal equilibrium with the contemporaneous CO₂ radiative forcing (gray line in figure 1(b)), defined as the realized warming fraction, is 0.43.

After the phase I, carbon emissions are stopped at 2 °C warming. The atmospheric CO₂ is simulated to decrease from 747 ppm in year 99 to 476 pm in year 1000, indicating that 21% of the cumulative carbon emissions are still in the atmosphere in year 1000. The land carbon uptake diminishes after 200 years and the ocean is the only sink thereafter. In this post-emissions phase II (figure 1(b)), global mean temperature continues to increase on multi-centennial timescales by 0.5 °C (or by 25%) until peak temperature is reached at the end of the 1000 year simulation (red point in figure 1(b)). Although the realized warming fraction in year 99 of 0.43 implies that if CO₂ were held fixed the global mean temperature would increase by a further 2.6 K, the decline in CO₂ means that the actual increase by year 1000 is only a further 0.5 K.

Significant warming occurs post-emissions despite a decline in radiative forcing R that exceeds the decline in ocean heat uptake N (figure 1(c)). The difference between radiative forcing and ocean heat uptake is usually used to explain the simulated temperature response in other models (Matthews et al 2009, Solomon et al 2009). The reason for the warming in GFDL ESM2M is that decreasing ocean heat uptake warms the Earth more in the global mean than the CO₂ forcing decrease is cooling. This is because the ocean heat uptake is localized in the subpolar ocean (Winton et al 2010, Frölicher et al 2014, Rose et al 2014) and perturbations at the high latitudes are restored less strongly by radiation to space than those at lower latitudes. In addition, there is no common reason to believe that the warming influence of decreasing ocean heat uptake should counteract exactly the forcing decline due to CO₂ uptake (Winton et al 2013, Frölicher et al 2015). Recent analysis haven shown that ocean storage of anthropogenic carbon and heat have distinct patterns, which is not consist with the view of passive process acting on both anthropogenic carbon and heat (Frölicher et al 2015).

After 1000 years, the warming is equal to the equilibrium climate temperature with the contemporaneous CO₂ levels (gray line in figure 1(b)), and both atmospheric CO₂ and global mean surface temperature start to slowly decrease. In this phase III, the ocean waters pH that has decreased in phase I and II will be restored by reactions between oceanic dissolved CO₂ and calcium carbonate of seafloor sediments, partly replenishing the buffer capacity of the ocean and further drawing down atmospheric CO₂ (Ciais et al 2013). Assuming a decrease in the airborne fraction of anthropogenic CO₂ in year 1000 by 20% to year 10 000 (see methods), the global mean temperature is estimated to be still 2.1 °C above preindustrial levels after 10 000 years (blue point in figure 1(b)) and still 0.1 °C above the temperature at the point of stopping emissions in year 99. Note that in the next several hundred thousands of years, CO₂ would be restored to its preindustrial levels through reactions with silicate rocks, a very slow process of CO₂ reaction with calcium silicate and other mineral of igneous rocks (Sundquist 1990). As a result, global mean temperature would slowly decrease to preindustrial levels.

Under this idealized global warming scenario, the GFDL ESM2M suggests the peak temperature change in response to cumulative carbon emissions 25% larger than the temperature change when emissions were stopped, and the multi-millennial temperature change after 10 000 years still 5% larger than the temperature change when emissions were stopped. The fact that peak temperature is larger than the temperature at the point when stopping emissions has profound implications for allowable carbon emissions to remain below a certain temperature target. This will be discussed in section 3.4.

Next, we will compare the GFDL ESM2M results with estimates from CMIP5–ESMs and EMICs to test the robustness of the GFDL ESM2M results and to test if CMIP5–ESMs and EMICs show a fundamentally different response. The particular carbon emissions scenario used for GFDL ESM2M has only been run over phase I by other CMIP5–ESMs and by EMICs, but the global temperature response in year 1000 and year 10 000 can be estimated using preexisting simulations (see methods).

During phase I, all models simulate a nearly linear increase in global mean temperature, consistent with earlier studies (Matthews et al 2009, Gillett et al 2013). At the end of phase I, the CMIP5–ESMs and EMICs as a class simulate similar temperature changes (2.5 ± 0.5 °C in CMIP5–ESMs; 2.5 ± 0.6 °C in EMICs; figures 2(a) and (b). Interestingly, the simulated realized warming fraction is significantly lower in CMIP5–ESMs (0.58 ± 0.06; fifth row in supplementary table 1) than in EMICs (0.70 ± 0.08; fifth row in supplementary table 2). In other words, the relative influence of the ocean in reducing the response to CO₂ radiative forcing is much higher in CMIP5–ESMs than in EMICs. The simulated temperature response and the realized warming fraction is lower in GFDL ESM2M than the mean of the CMIP5–ESMs and EMICs in phase I (figures 2(a), (b), figure 1(b), supplementary table 1 and 2).

Zoom In Zoom Out Reset image size

During phase II, the estimated temperature response between the CMIP5–ESMs (including GFDL ESM2M) and the EMICs diverge. Assuming an airborne fraction of anthropogenic CO₂ of 25% at the end of phase II in year 1000 (see methods), the CMIP5–ESMs are estimated to simulate an increase in temperature of 0.2 ± 0.5 °C over the period 99–1000, whereas the EMICs as a class are estimated to simulate a decrease in temperature of −0.6 ± 0.2 °C. Note however that the uncertainties in the estimated temperature response for the CMIP5–ESMs is large, especially when additionally taking into account the uncertainties originating from the estimated airborne fraction (rows seven and eight in supplementary table 1). Here, it is assumed for all models that the climate system is in equilibrium with the contemporaneous CO₂ radiative forcing in year 1000. Extending the equilibration time would not significantly change the temperature response, as the airborne fraction of atmospheric CO₂ only slowly decreases after year 1000. Reducing the equilibration time, however, would lead to an additional increase in temperature after zero emissions as the airborne fraction is higher before year 1000.

During phase III (period 1000–10 000), the temperature is estimated to decrease by −0.3 °C ± 0.1 °C in the CMIP5–ESMs. The temperature in year 10 000 is only 0.1 ± 0.5 °C below the temperature in year 99 in CMIP5–ESMs when carbon emissions are ceased. The EMIC models simulate an additional temperature decrease of −0.2 ± 0.1 °C from year 1000 to 10 000.

The estimated multi-centennial decrease in temperature by the EMICs is consistent with the study by Zickfeld et al (2013) who analyzed preindustrial CO₂ emissions commitment in a few future emissions scenarios with the same EMIC models. They show that even when the positive non-CO₂ forcing is kept constant after zero CO₂ emission, the temperature decreases on multi-centennial timescales. The few existing multi-centennial-scale modeling simulations with comprehensive ESMs show contrasting responses. The GFDL ESM2M for example, simulates ongoing warming after a CO₂ emissions stoppage, not simulated by both NCAR CSM1.4-carbon (Frölicher and Joos 2010) and CanESM1 (Gillett et al 2011). The estimated multi-millennial warming responses by year 10 000 of 2.8 ± 0.8 °C per 1919 ± 154 Pg C emissions for CMIP5–ESMs and of 2.1 ± 0.6 °C per 1658 ± 90 Pg C emissions for EMICs are consistent with previous theoretical predictions of 1.6–3.8 °C per 2000 Pg C (Williams et al 2012), and qualitatively in line with single intermediate complexity model simulations that were run out for the entire 10 000 years (Eby et al 2009).

Variations in the long-term warming response between models are caused by both variations in the airborne fraction of anthropogenic carbon emissions and variations in the realized warming fraction at time when CO₂ emissions cease. Figure 3 shows the relationship between the multi-centennial temperature response (T₁₀₀₀–T₉₉) after zero carbon emissions and the realized warming fraction at time of zeroing carbon emissions (i.e. T₉₉/ECS₉₉) for all models investigated in this study. There is a nearly linear relationship between the multi-centennial temperature response and the realized warming fraction, especially in the CMIP5–ESMs. Models with low realized warming fraction show a positive multi-centennial temperature change, whereas models with high realized warming fraction show negative multi-centennial temperature change. The CMIP5–ESMs have an overall lower realized warming fraction than the EMICs and therefore show a positive temperature response post CO₂ emissions. According to the CMIP5–ESMs, models with a realized warming fraction lower than 0.61 (0.47 in EMICs) are estimated to simulate ongoing warming after zero carbon emissions. One might argue that the CMIP5–ESMs results are highly uncertain as these models are not run until year 1000 and therefore the airborne fraction of the anthropogenic CO₂ emission is not known. However, the EMICs, for which the airborne fractions in year 1000 are known, show a similar relationship between realized warming fraction and long-term climate change (blue dashed line for CMIP5–ESMs and red dashed line for EMICs in figure 3), although with lower confidence.

Zoom In Zoom Out Reset image size

The fact that peak temperature is larger than the temperature at the point when stopping emissions in GFDL ESM2M and CMIP5–ESMs has profound implications for allowable carbon emissions to remain below a policy-driven temperature target, such as the 2 °C global warming target. As illustratively shown for GFDL ESM2M in figure 4(a), the relationship between cumulative carbon emissions and global mean temperature change shows a close to linear relationship when carbon emissions increase (phase I). This near constancy of TCRE is consistent with earlier modeling results (Mathews et al 2009, Gillet et al 2013). The simulated increase in temperature of 0.5 °C in GFDL ESM2M (from green to red point in figure 4(a)) without a change in cumulative carbon emissions during phase II and the estimated temperature decrease in phase III (from red point to blue point in figure 4(a)) cause changes in the coefficient relating global warming to cumulative carbon emissions. The increase in temperature after zero emissions in GFDL ESM2M suggest that estimates of allowable carbon emissions required to remain below the 2 °C global warming target are 20% lower (1558 Gt C instead of 1947 Gt C) than previously thought on a 1000 years timescale, and 5% lower (1854 Gt C instead of 1947 Gt C) on a 10 000 years timescale (figure 4(b)).

Zoom In Zoom Out Reset image size

We introduce two new metrics that better characterize the time-dependent temperature response to cumulative carbon emissions, the ECRE and the MCRE. The ECRE is calculated by using the temperature change at the time when ocean heat uptake is zero, i.e. in year 1000 when the system is in equilibrium with its contemporaneous CO₂ forcing, and then by dividing this by the cumulative CO₂ emissions up to this time. Usually, ECS is defined as the change in global mean temperature that results when the climate system attains a new equilibrium with the radiative forcing change resulting from a doubling of atmospheric CO₂. In our particular simulation, atmospheric CO₂ will further decrease. Therefore, we also introduce the MCRE, which takes into account the decrease in atmospheric CO₂. The MCRE is calculated identically to the ECRE, but defined on multi-millennial timescales (at year 10 000) when the exchange of carbon between the atmosphere and global ocean ceases (i.e. the ocean and land biosphere are in equilibrium with the contemporaneous atmospheric CO₂), and only weathering with silicate rocks can restore atmospheric CO₂ levels to preindustrial.

Applying the TCRE, ECRE and MCRE to a 2 °C global warming target in CMIP5–ESMs, TCRE is 1388 ± 388 Gt C, ECRE is 1328 ± 389 Gt C and MCRE is 1465 ± 427 Gt C (figure 4(b)). Thus, MCRE is 6% larger than TCRE, but ECRE is 4% smaller than TCRE in CMIP5–ESMs. This implies that less carbon emissions are required to hit the 2 °C target when using ECRE than when using TCRE. This is not the case for EMICs, for which peak temperature response is attained when carbon emissions are ceased and therefore ECRE (1590 ± 482 Gt C) and MCRE (1760 ± 536 Gt C) are generally larger than TCRE (1189 ± 224 Gt C; figure 4(b)).