CO 2 Increase Experiments Using the CESM: Relationship to Climate Sensitivity and Comparison of CESM1 to CESM2

Abstract

We examine the response of the Community Earth System Model Versions 1 and 2 (CESM1 and CESM2) to abrupt quadrupling of atmospheric CO₂ concentrations (4xCO2) and to 1% annually increasing CO₂ concentrations (1%CO2). Different estimates of equilibrium climate sensitivity (ECS) for CESM1 and CESM2 are presented. All estimates show that the sensitivity of CESM2 has increased by 1.5 K or more over that of CESM1. At the same time the transient climate response (TCR) of CESM1 and CESM2 derived from 1%CO2 experiments has not changed significantly—2.1 K in CESM1 and 2.0 K in CESM2. Increased initial forcing as well as stronger shortwave radiation feedbacks are responsible for the increase in ECS seen in CESM2. A decomposition of regional radiation feedbacks and their contribution to global feedbacks shows that the Southern Ocean plays a key role in the overall behavior of 4xCO2 experiments, accounting for about 50% of the total shortwave feedback in both CESM1 and CESM2. The Southern Ocean is also responsible for around half of the increase in shortwave feedback between CESM1 and CESM2, with a comparable contribution arising over tropical ocean. Experiments using a thermodynamic slab‐ocean model (SOM) yield estimates of ECS that are in remarkable agreement with those from fully coupled Earth system model (ESM) experiments for the same level of CO₂ increase. Finally, we show that the similarity of TCR in CESM1 and CESM2 masks significant regional differences in warming that occur in the 1%CO2 experiments for each model.

Introduction

The coupled climate system responds in complicated ways to anthropogenic changes in greenhouse gas concentrations, aerosol emissions, and land use, among other factors. To investigate climate model response to these forcings, two idealized configurations were introduced in the Coupled Model Intercomparison Project Phase 5 (CMIP5; Taylor et al., 2012): (1) the abrupt 4xCO2 increase experiment and (2) the 1%CO2 increase experiment. For both experiments, a fully coupled atmosphere‐ocean general circulation model (AOGCM) or Earth system model (ESM) is run to a spun‐up state using estimated preindustrial (year ≈ 1850) greenhouse gas concentrations, aerosol emissions, land use, and other climate forcings. The spun‐up preindustrial control run (piCTL) is then subjected to an abrupt quadrupling of atmospheric CO₂, or to 1% annually increasing CO₂, while holding all other forcings at preindustrial levels. Both experiments are part of the initial Diagnostic, Evaluation and Characterization of Klima (DECK) requirements for participation in Phase 6 of the CMIP (CMIP6; Eyring et al., 2016).

Equilibrium climate sensitivity (ECS) is defined as the equilibrium warming that would occur under a doubling of CO₂ (Charney et al., 1979). The abrupt 4xCO2 increase scenario was introduced to evaluate model climate sensitivity. The CMIP 4xCO2 experimental design calls for 150 yr of simulation, although the long oceanic time scales in the climate system imply that coupled simulations may require ∼1,000 years to reach a new equilibrium (e.g., Danabasoglu & Gent, 2009; Rugenstein et al., 2019). ECS has been estimated from 4xCO2 experiments using linear regression to global‐mean top‐of‐atmosphere (TOA) (or top‐of‐model [TOM]) radiative imbalance $\bar{\mathcal{N}}$ and global‐mean warming $\mathrm{\Delta }\bar{T}$ (Gregory et al., 2004). The linear fit to $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ is extrapolated to $\bar{\mathcal{N}}={\bar{\mathcal{N}}}_{\text{pi}}$, the spun‐up preindustrial control value, to estimate an equilibrium warming $\mathrm{\Delta }{\bar{T}}_{eq}$, which is divided by 2 (under the assumption of linearity) to estimate ECS. We will refer to the ECS estimate derived in this way as the inferred ECS or iECS. The iECS approach was applied to 150‐yr 4xCO2 AOGCM/ESM simulations to derive the published ECS values for CMIP5 (Flato et al., 2014). Recent analysis by Zelinka et al. (2020) shows that the ensemble mean iECS has increased from 3.3 K in CMIP5 models to 3.9 K in CMIP6 with many CMIP6 models possessing an iECS greater than 4.5 K.

Another approach to estimating ECS was proposed by Danabasoglu and Gent (2009), using a thermodynamic slab‐ocean model (SOM) rather than a full dynamical ocean in abrupt CO₂ increase experiments to eliminate the long time scales produced by the slow deep‐ocean responses to warming. The SOM experiments equilibrate in decades rather than centuries, yielding a SOM‐based estimate of ECS (ECS‐SOM).

Both the iECS and ECS‐SOM approaches to estimating the true ECS of a coupled model have shortcomings. A drawback of the iECS approach is that $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ may be a nonlinear function of $\mathrm{\Delta }\bar{T}$, leading to iECS values that depend on the number of years in the regression analysis (e.g., Andrews et al., 2012; Williams et al., 2008). The radiative response to abrupt CO₂ increase is also known to be nonlinear (e.g., Etminan et al., 2016). In the case of ECS‐SOM, it is unclear whether details in the construction of a SOM configuration can affect the resulting ECS (e.g., Danabasoglu & Gent, 2009; Senior & Mitchell, 2000; Stouffer & Manabe, 1999; Williams et al., 2008).

Gettelman, Hannay, et al. (2019) examined the evolution of ECS‐SOM in CESM and utilized prescribed SST configurations with perturbed SSTs and aerosols to explore potential causes for the increased ECS‐SOM in CESM2 (5.3–5.5 K) compared to that of CESM1 (4.0–4.2 K). Using prescribed SST+4 K experiments (Cess et al., 1990) they identified possible impacts on model sensitivity from warm phase and mixed‐phase microphysics including the removal of a negative phase feedback over the Southern Ocean in moving from CESM1 to CESM2. In this study we will examine the climate response in fully coupled 4xCO2 and 1%CO2 CMIP simulations using CESM1 and CESM2. As in Gettelman, Hannay, et al. (2019) we find that Southern Ocean clouds are an important component of the increased ECS of CESM2, but we find less compelling evidence for the role of phase feedbacks. In addition, we show an important and complex set of tropical cloud changes that have also contributed to raising the ECS of CESM2. Finally, we relate results from the 4xCO2 experiments to the 1%CO2 scenario and transient climate response (TCR; Taylor et al., 2012) using the step model approach of Gregory et al. (2015). We note finally that the bulk of the analysis here is carried out using standard CMIP experiments and output.

Basic Results

Table 1 summarizes values of ECS‐SOM and iECS for Versions 1 and 2 of the Community Earth System Model (CESM). All estimates of ECS have increased substantially in Version 2 of CESM (CESM2; Danabasoglu et al., 2020). In addition to the increases in ECS‐SOM already noted, iECS has also increased in CESM2. Values of iECS based on 150 yr of simulation results have increased by almost 2 K in CESM2 compared to CESM1 (Column 3 of Table 1). Values of iECS based on 800 yr of simulation results (Column 4), which are substantially higher than those based on 150 yr, have increased by 2.3 K.

1 TableMeasures of Climate Response Discussed in This Analysis

Equilibrium climate	ECS‐SOM(4x)	Inferred ECS		Transient
sensitivity based on	based on	(iECS) based	iECS	climate
2xCO2‐SOM experiments	4xCO2‐SOM	on 150‐yr	based on	response
(ECS‐SOM)	experiments	regression	800 yr	(TCR)
CESM1
4.0, 4.2	4.2 (0.03)	3.4 (0.04)	4.2 (0.05)	2.1 (0.07)
CESM2
5.5 (0.03)	6.5 (0.07)	5.3 (0.22)	6.5 (0.07)	2.0 (0.04)

Note. All values in degrees kelvin (K). SOM‐based numbers for CESM1 equilibrium climate sensitivity (ECS‐SOM; first row, first column) are taken from Gettelman et al. (2012). Two numbers given are for 1° and 2° horizontal resolutions, respectively. All other numbers were calculated for this study. Details of the calculations are given in Appendix A. Second column shows ECS‐SOM(4x) based on SOM runs subject to a 4xCO2 increase (section 5). Inferred ECS (iECS; third and fourth columns) is derived from linear regression analysis of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ from abrupt CO₂ increase experiments (Gregory et al., 2004). Transient climate response (TCR; fifth column) is derived from experiments subject to a 1% annual CO₂ increase (section 6). Standard errors, where available, are shown in parentheses.

Figure 1 illustrates key features of 4xCO2 experiments using CESM1 and CESM2. Figure 1a shows global‐mean TOM radiative imbalance $\bar{\mathcal{N}}$ as a function of global‐mean surface temperature $\mathrm{\Delta }\bar{T}$ for CESM1 (black) and CESM2 (red). The equilibrium temperature of the respective piCTL simulation has been subtracted from $\bar{T}$ to give $\mathrm{\Delta }\bar{T}$. Although the fully coupled 4xCO2 runs shown in Figure 1a are over 800 yr in length, they have not equilibrated. Also, we see that $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ for both CESM1 and CESM2 exhibits nonlinearity (e.g., Andrews et al., 2012), that is, a change in the slope of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ with warming. The presence of such nonlinearity has been attributed to rapid nonlinear low‐cloud SST feedbacks (Williams et al., 2008) and multiple time scales of deep‐ocean heat uptake (e.g., Held et al., 2010; Li et al., 2013; Senior & Mitchell, 2000).

1 Figure. (a) Annual‐mean, global top‐of‐model radiation imbalance N‾ as a function of annual‐mean, global‐mean surface temperature change ΔT‾ for abrupt 4xCO2 experiments CESM1‐4xCO2 (black) and CESM2‐4xCO2 (red). Dashed lines show linear fits to N‾(ΔT‾) for Years 100–800. Two points are indicated on each N‾(ΔT‾) relationship: Values of linear fits at Year 100 and diagnosed inflection points. (b) Inferred equilibrium climate sensitivities (iECS) from linear regressions: Horizontal axis gives number of years used in the regression. Long curves extending to 800 yr and beyond show iECS derived for CESM1‐4xCO2 (black) and CESM2‐4xCO2 (red) from linear regressions of N‾(ΔT‾). Shorter red curves show iECS derived from a 2xCO2‐SOM experiment with CESM2 (CESM2‐2xCO2‐SOM, Table 2) and from a 4xCO2‐SOM experiment with CESM2 (CESM2‐4xCO2‐SOM). Short black indicates iECS derived from CESM1b‐4xCO2‐SOM. Black and red triangles on right vertical axis show values of ECS‐SOM for CESM1 (4.0 and 4.2 K) and CESM2 (5.5 K). (c) Global‐mean surface temperature T‾ as a function of time for CESM1‐4xCO2 (black) and CESM2‐4xCO2 (red). (d) As (c) except focusing on first 200 yr of experiments. Gray line shows results for CESM1b‐4xCO2.

Figure 1b shows ECS inferred from linear regressions (iECS) of $\bar{\mathcal{N}}$ versus $\mathrm{\Delta }\bar{T}$ as a function of years used in the regression. CESM1 exhibits a long initial period (∼150 yr) during which iECS is relatively constant near 3.5 K, or even weakly decreasing, before increasing to values slightly over 4 K by Year 800. In CESM2, however, iECS increases rapidly from Year 20 onward and quickly exceeds the published ECS of 5.4 K (Gettelman, Hannay, et al. 2019) between Years 150 and 200. The iECS for CESM1 derived from the full 150 yr of the prescribed 4xCO2 experiment is around 3.4 K, well below the value derived from SOM runs or from longer periods of the 4xCO2 run. In CESM2, the iECS in Year 150 is around 5.5 K but approaches 6.5 K as more years are used in the regression. In section 5 we will show that the iECS at long times agrees with ECS‐SOM with 4xCO2 forcing for both CESM1 and CESM2.

Figures 1c and 1d show time series of $\bar{T}$ for CESM1 and CESM2, again with interesting differences between the two models. In CESM1 an extended pause (hiatus) in warming sets in after a short initial period of rapid warming. The hiatus lasts for around 100 yr, after which gradual warming resumes. Warming in CESM2 has no such hiatus; rates of warming decrease consistently over the integration. The warming hiatus in CESM1 appears to be the ultimate cause of the local minimum in iECS around Year 100 (Figure 1b).

TCR, defined as the global‐mean warming averaged over Years 60–80 in the 1%CO2 experiment with respect to piCTL, is also shown in Table 1. Interestingly, TCR values have changed little between CESM1 (2.1 K) and CESM2 (2.0 K) despite the large changes in ECS.

Outline

In the remainder of this paper we will address three topics: (1) Global and regional radiation feedbacks and their role in increasing climate sensitivity between CESM1 and CESM2; (2) The relationship between SOM‐based estimates of ECS and those from fully coupled ESM runs using a dynamic ocean; and (3) The behavior of the 1%CO2 configurations of CESM1 and CESM2 and its relation to TCR. We find that the increased climate sensitivity of CESM2 arises from both stronger shortwave radiation feedbacks with surface temperature T_s and from increased initial forcing ${\bar{\mathcal{N}}}_0$. The strengthened shortwave feedbacks in CESM2 originate primarily in low‐cloud shortwave feedbacks over the Southern Ocean and also from a complex pattern of tropical cloud feedbacks. We find that SOM‐based estimates of ECS agree with those based on full ESM simulations, despite differences in regional warming patterns. We will also see that 1%CO2 experiments for CESM1 and CESM2 differ more than is implied by the similar values of TCR. In particular, TCR does not capture significant regional variations between the models.

The paper is organized as follows: Section 2 briefly describes CESM1, CESM2, the CESM Slab‐Ocean Model, and the experimental setups used in this study. A notable feature of this study is a comparison of fully coupled 4xCO2 ESM integrations with 4xCO2‐SOM integrations. Section 3 details the model variables examined and describes analysis methods, including a consistent treatment of regional versus global feedback parameters. Section 4 describes results from the fully coupled 4xCO2 experiments, including analysis of longwave and shortwave radiative responses (section 4.1), regional decomposition of feedbacks (section 4.3), and an analysis of physical processes contributing to feedbacks (section 4.4). Section 5 describes SOM‐based abrupt CO₂ increase experiments and compares them with full ESM results. Section 6 examines results from 1%CO2 experiments using CESM1 and CESM2. Finally, section 7 summarizes our results and discusses implications of the various measures of climate response.

Models and Experimental Design

CESM2 and CESM1

The CESM Version 2 (CESM2; Danabasoglu et al., 2020) was developed over 5 yr for participation in CMIP6 (Eyring et al., 2016). This development was finished in December 2018, and CMIP6 DECK simulations with CESM2 are now complete. Its predecessor model, CESM1 (Hurrell et al., 2013), has been extensively documented. The versions of CESM1 examined here are those used in the Last‐millenium ensemble project (LME Otto‐Bliesner et al., 2016) and the CESM Large Ensemble project (LENS; Kay et al., 2015). The only differences between these versions are the atmospheric horizontal resolution, 2° for LME and 1° for LENS, as well as some retuning of low‐cloud fraction. Results of the preindustrial and twentieth century historical simulations using the LME version of CESM1 were contributed to the CMIP5 archive as “CESM1 (CAM5.1, FV2)”.

CESM2 incorporated major changes to several component models, including atmosphere, land, and ocean. A new interactive model of the Greenland Ice Sheet (Lipscomb et al., 2019) was also introduced. (Ice sheet elevation and extent were held fixed, however, in the simulations analyzed here.) In addition to component development, emissions data sets and other forcing data sets were substantially revised for CMIP6 (Hoesly et al., 2018).

The CESM2 atmosphere component differs substantially from that in CESM1. Every physics parameterization, except for the rapid radiative transfer model for GCM applications (RRTMG; Iacono et al., 2008), was replaced or modified (Richard Neale, personal communication). The major physics changes relevant to cloud and turbulence processes are the replacement of shallow convection, boundary layer turbulence, and cloud macrophysics schemes in CESM1 with the Cloud Layers Unified by Binormals (CLUBB; Bogenschutz et al., 2013) scheme and an update of cloud microphysics from the Morrison‐Gettelman scheme Version 1 (MG1; Morrison & Gettelman, 2008) to MG2 (Gettelman et al., 2015). In addition modifications were made to implementation of the Zhang‐Mcfarlane deep convection scheme (ZM; Zhang & McFarlane, 1995) that resulted in generally shallower convective updrafts and cloud.

CLUBB is a turbulence and shallow‐convection scheme based on higher‐order closure, employing 10 higher‐order moments of subgrid vertical velocity w′, temperature T′, and total moisture qt′. CLUBB also produces large‐scale cloud fraction and partitions between condensed and vapor phase water. MG2 is a sophisticated two‐moment cloud microphysics scheme that explicitly models the interactions between clouds and aerosols. MG2 extends MG1 by including prognostic equations for rain and snow in addition to cloud ice and liquid. MG2 also includes changes to the treatment of mixed‐phase ice nucleation that have led to increased amounts of supercooled liquid in mixed‐phase clouds.

Updates to ocean, land, land ice, and s ea ice components in CESM2 are discussed by Danabasoglu et al. (2020) and references therein.

Experimental Design

Abrupt 4xCO2 and 1%CO2 increase experiments are branched from spun‐up, fully coupled preindustrial control (piCTL) experiments in which all forcing (e.g., aerosol emissions, greenhouse gases, and land use) is fixed at estimated 1850 levels. A CESM piCTL run is considered equilibrated if TOM radiative imbalance $|\bar{\mathcal{N}}|<0.1$ W m⁻² in a 20‐yr mean. The CESM1 and CESM2 piCTL experiments used to initialize the CO₂ increase experiments are each over 1,150 yr in length. The 4xCO2 and 1%CO2 scenarios were branched off in Year 1,000 of the CESM1 piCTL experiment and in Year 501 of the CESM2 piCTL.

In the 4xCO2 scenario, atmospheric CO₂ is abruptly quadrupled after branching, and the climate is allowed to evolve freely. The typical evolution of such runs is illustrated in Figure 1. In 1%CO2 experiments, an annually compounding increase in atmospheric CO₂ is imposed after branching, with other forcing fixed to piCTL values. For the CESM2 experiments discussed here, radiatively active species other than CO₂, notably ozone, are specified from piCTL experiments using the high‐top Whole Atmosphere Community Climate Model (WACCM; Gettelman, Mills, et al. 2019) with fully interactive chemistry. This procedure is discussed in detail by Danabasoglu et al. (2020). Impacts of this procedure on the evolution of CO₂ increase scenarios using CESM are under investigation, but will not be discussed here.

Table 2 summarizes the experiments discussed in this paper. We examine results from the 4xCO2 experiment performed for CMIP6 (CESM2‐4xCO2) as well as two 4xCO2 experiments using CESM1: CESM1‐4xCO2, performed with the LME version at 2° horizontal resolution, and CESM1b‐4xCO2, performed with the LENS version at 1° horizontal resolution. As noted in the table, the CESM1‐4xCO2 and CESM2‐4xCO2 experiments are significantly longer that the 150 yr requested in the CMIP protocol. As seen in Figure 1, equilibration of 4xCO2 experiments may take ∼1,000 yr or longer. We also examine results from the CESM2 1%CO2 run performed for CMIP6 (CESM2‐1%CO2 and from a CESM1‐1%CO2 run performed with the LME version of CESM1.

2 TableCESM Experiments Discussed in This Analysis and Their Shorthand Designations

Designation	Model version	Horizontal resolution	Setup	Length (years)
Fully coupled, Earth System Model (ESM) runs
CESM1‐piCTL	CESM1	2°	Preindustrial control	100
CESM1‐4xCO2	CESM1(LME)	2°	Abrupt 4xCO2 increase	800
CESM1b‐4xCO2	CESM1(LENS)	1°	Abrupt 4xCO2 increase	200
CESM2‐4xCO2	CESM2.1	1°	Abrupt 4xCO2 increase	1,000
CESM1‐1%CO2	CESM1(LME)	2°	1% annual CO2 increase	190
CESM2‐1%CO2	CESM2.1	1°	1% annual CO2 increase	150
Slab‐ocean model (SOM) runs
CESM1b‐4xCO2‐SOM	CESM1(LENS)	1°	Abrupt 4xCO2 increase	30
CESM2‐4xCO2‐SOM	CESM2.0	1°	Abrupt 4xCO2 increase	100
CESM2‐2xCO2‐SOM	CESM2.0	1°	Abrupt 2xCO2 incease	100

SOM Experiments

We also conducted abrupt CO₂ increase experiments using the CESM SOM. The CESM‐SOM configuration relies on ocean parameters derived from equilibrated, preindustrial control simulations and is designed to reproduce the climate of the fully coupled ESM configuration (Bitz et al., 2012). The parameters used by the SOM are 2‐D annual‐mean estimates of ocean mixed layer depths along with 2‐D monthly heat flux anomalies to the deep ocean. These parameters are used to drive an interactive thermodynamic slab that is forced from above by atmospheric fluxes. By construction, the global‐mean deep‐ocean heat flux is identically 0. ECS estimates for CESM and predecessors using 2xCO2‐SOM simulations have been reported in several studies (e.g., Bitz et al., 2012; Danabasoglu & Gent, 2009; Gettelman et al., 2012; Gettelman, Hannay, et al., 2019). Here we will examine both 4xCO2 and 2xCO2‐SOM experiments with CESM to quantify nonlinearity in ECS estimates and to enable direct comparison with fully coupled experiments.

In the following, we append “SOM” to any experiments using the slab‐ocean configuration. Experiments using fully coupled CESM do not normally have a descriptive suffix, for example, “CESM2‐4xCO2”. If clarity is a concern, the latter are designated as “ESM” experiments.

Model Output and Analysis Methods

The analyses presented here use monthly and annually averaged output from CESM, including radiative fluxes, cloud condensates, and surface temperature. We use TOM radiation fluxes rather than estimated TOA fluxes, and surface temperature T_s rather than 2‐m air temperatures T_2m. The results are not sensitive to the TOM versus TOA distinction or the T_s versus T_2m distinction (Table A1 and Appendix A). Throughout this analysis T will always refer to surface temperature T_s.

Net TOM shortwave and longwave fluxes are denoted by $\mathcal{S}$ and $\mathcal{L}$, respectively. The TOM radiative imbalance $\mathcal{N}$, already introduced in Figure 1, is simply [Image Omitted. See PDF]

We follow the usual atmospheric convention of defining upward longwave radiation flux and downward shortwave flux as positive.

CESM atmospheric model output also includes shortwave and longwave cloud radiative effect (CRE) ${\mathcal{S}}_{cld}$ and ${\mathcal{L}}_{cld}$, as well as TOM clear sky fluxes ${\mathcal{S}}_{clr}$ and ${\mathcal{L}}_{clr}$. These are calculated directly in the CAM radiation scheme in each grid column and time step and are approximately related to all‐sky fluxes by [Image Omitted. See PDF] [Image Omitted. See PDF] where a small residual (∼0.05 W m⁻²) exists due to the definition of CRE at TOA instead of TOM. CESM follows the usual sign conventions for CRE: Negative ${\mathcal{S}}_{cld}$ indicates reflection of shortwave radiation by clouds, and positive ${\mathcal{L}}_{cld}$ indicates downward longwave radiation from clouds.

We also examine simulated total cloud amount c from CESM. This is calculated using the random overlap assumption across three cloud macrolayers bounded by the surface, 700, 400, and 50 hPa.

This means the net clear fraction across adjacent macrolayers is parameterized as (1 − c_k)(1 − c_k − 1). Within each cloud macrolayer a fraction is calculated using maximum‐random cloud overlap, meaning that when there are clouds in adjacent layers they are maximally overlapped, that is, the net fraction is the greater of the two. Finally, we will examine liquid and ice cloud condensate paths (LWP and IWP, g m⁻²). An estimate of in‐cloud condensate paths is calculated by dividing monthly grid means of LWP and IWP by the cloud amount c, that is, [Image Omitted. See PDF] [Image Omitted. See PDF]

Regional and Global Feedback Parameters

Studies of climate sensitivity focus on feedback relationships of the form [Image Omitted. See PDF] where X is a flux or other quantity of interest, T is surface temperature, and λ_X is a feedback parameter (slope) that linearly relates changes in X and T. X and T may represent regional or global‐mean quantities (e.g., Armour et al., 2013). Below, we will establish quantitative relationships between regional feedbacks and global feedbacks. We will be primarily interested in feedbacks between radiative fluxes and temperatures.

The global mean of X can be written as a sum of regional means over N regions, [Image Omitted. See PDF] where X_k is the mean of X in region k, T_k is the regional‐mean surfce temperature, and a_k is the areal fraction of region k. Global means will be denoted by $\overline{()}$ throughout this analysis.

The regional means X_k on the RHS of Equation 5 may depend on variables other than the regional surface temperature, including surface temperatures in other regions, or other meteorological variables such as vertical velocity or stability. We will assess the functional relationships between regional quantities and regional surface temperature T_k by examining scatterplots of X_k versus T_k. If the points in a scatterplot fall close to a curve, over a range of T_k, we assume we are justified in assuming X_k ≈ X_k(T_k).

The global feedback parameter ${\bar{\lambda }}_X$ between $\bar{X}$ and $\bar{T}$ can then be estimated from a sum of regional feedbacks according to [Image Omitted. See PDF]

We approximate the derivatives on the RHS of Equation 6 with slope parameters from linear regressions of X_k versus T_k and of T_k versus $\bar{T}$. The linear regression slope of X_k versus T_k is simply the regional feedback parameter for X in region k and will be denoted λ_X; k. The linear regression slope of T_k versus $\bar{T}$ is the regional warming rate divided by the global rate. This is the amplification factor for regional warming and will be denoted by A_k. With these approximations, we rewrite Equation 6: [Image Omitted. See PDF]

The global feedback parameter ${\bar{\lambda }}_X$ has thus been written as a weighted sum of local feedbacks λ_X; k. The validity of regional decomposition can be tested by comparing the sum in Equation 7 with an independent regression using global‐mean quantities. This will be shown in section 4.3.

CESM1‐4xCO2 has large interannual variability compared CESM2‐4xCO2 (e.g., Figure 1d), likely related to strong ENSO. This is associated with correlated subdecadal variations in $\mathcal{S}$ and T that have small but significant effects on linear regression estimates of ${\lambda }_{\mathcal{S}}$. For the analysis of long‐term regional feedbacks we apply a decadal average to model results. Decadal averaging has negligible impacts on the analysis of CESM2‐4xCO2 results. Its impacts in the analysis of CESM1‐4xCO2 are largely restricted to calcuation of shortwave feedbacks in the tropics, and will be discussed further in section 4.

Approximate Partial Radiative Perturbations

We will examine cloud contributions to shortwave radiative forcing using the approximate partial radiative perturbation approach (APRP; Taylor et al., 2007). APRP constructs an analog to the full shortwave radiation calculation in an atmospheric model using monthly fields of clear‐sky and all‐sky shortwave fluxes at TOM and at the surface, as well as monthly total cloud amounts. The result is a reconstructed planetary albedo $\mathcal{A}$ that depends on seven parameters [Image Omitted. See PDF] where c again is total cloud amount; α_clr and α_oc are clear‐sky and overcast surface albedos; μ_clr and μ_cld are clear‐sky and cloudy‐sky absorption coefficients; and γ_clr and γ_cld are clear‐sky and cloudy‐sky scattering coefficients. The albedo and net all‐sky TOM shortwave flux $\mathcal{S}$ are related by [Image Omitted. See PDF] where ${\mathcal{S}}^{\downarrow }$ is the incoming shortwave radiation at TOM. The APRP method provides estimates of the albedos, and absorption and scattering coefficients as well as an analytical expression for $\mathcal{A}$ that can be used to calculate partial derivatives and quantify the impact of different processes on shortwave radiation in the atmosphere. Given the importance of high‐latitude responses in warming climates (e.g., Kay et al., 2014), it is particularly important to distinguish the roles of surface and cloud processes in the overall feedback.

Rapid and Long‐Term Time Scales

Several studies (e.g., Held et al., 2010) have noted the existence of multiple time scales in the adjustment of the coupled climate system to abrupt perturbations. The behavior of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ shown in Figure 1a suggests the existence of at least two phases in the evolution of CESM after an abrupt quadrupling of CO₂. There is an initial phase with rapid warming and steep negative slope in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$, followed by a slower adjustment with nearly constant but shallower negative slope in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$, that persists until the end of both 4xCO2 experiments. The time evolution of $\bar{T}$ in CESM1 includes a long pause in warming from Years 20 to 100 (Figures 1c and 1d). During this pause, there is little evolution of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$, with values of $\mathrm{\Delta }\bar{T}$ and $\bar{\mathcal{N}}$ fluctuating around 5 K and 2 W m⁻², respectively. Then warming in CESM1 resumes, and $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ is approximately linear with a slope of about −0.6 W m⁻² K⁻¹. Based on this behavior, we identify Years 1–20 as representative of the rapid initial adjustment of both 4xCO2 experiments. Year 100 is approximately when the slope in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ for CESM1‐4xCO2 changes and the transition in the slope of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ occurs earlier in CESM2‐4xCO2. For simplicity, we choose Years 100–800 to describe the long‐term behavior of both experiments.

We use linear regressions of $\bar{\mathcal{N}}$, $\bar{\mathcal{S}}$, and $\bar{\mathcal{L}}$ versus $\bar{T}$ over Years 1–20 of the 4xCO2 experiments, extrapolated to their corresponding piCTL equilibrium $\bar{T}$ values to estimate initial radiative forcing ${\bar{\mathcal{N}}}_0$ and ultrarapid longwave and shortwave adjustments $\mathrm{\Delta }{\mathcal{L}}_0$ and $\mathrm{\Delta }{\mathcal{S}}_0$, which are given in Table 3.

3 TableInitial Radiative Imbalance ${\bar{\mathcal{N}}}_0$ and Rapid Initial Adjustments to Longwave ( $\mathrm{\Delta }{\bar{\mathcal{L}}}_0$) and Shortwave Fluxes ( $\mathrm{\Delta }{\bar{\mathcal{S}}}_0$) in 4xCO2 Experiments

${\bar{\mathcal{N}}}_0$ (W m−2)	$\mathrm{\Delta }{\bar{\mathcal{L}}}_0$ (W m−2)	$\mathrm{\Delta }{\bar{\mathcal{S}}}_0$ (W m−2)
CESM1‐4xCO2
7.4	−7.6	−0.2
CESM2‐4xCO2
8.6	−7.6	1.0

Note. Numbers are diagnosed from linear fits to $\bar{\mathcal{N}}$, $\bar{\mathcal{L}}$, and $\bar{\mathcal{S}}$ during Years 1–20 of CESM1‐4xCO2 and CESM2‐4xCO2. Regression parameters are used to extrapolate $\bar{\mathcal{N}}$, $\bar{\mathcal{L}}$, and $\bar{\mathcal{S}}$ to the equilibrium $\bar{T}$ from the corresponding piCTL experiment (or equivalently to $\mathrm{\Delta }\bar{T}=0$).

Results From 4xCO2 Experiments

Here we will examine results from the extended 4xCO2 experiments, focusing on processes that contribute to the increased climate sensitivity of CESM2 compared to that of CESM1. As described in Appendix A, iECS is derived from linear fits to $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$. [Image Omitted. See PDF] where ${\bar{\mathcal{N}}}_I$ and ${\bar{\lambda }}_{\mathcal{N}}$ are the intercept and slope of the linear fit, and the factor of 0.5 scales 4xCO2 results to a 2xCO2 scenario assuming linearity (see Appendix A). In physical terms, ${\bar{\lambda }}_{\mathcal{N}}$ is the net radiation feedback with respect to T and ${\bar{\mathcal{N}}}_I$ is an estimate of the initial radiative forcing (which is equal to ${\bar{\mathcal{N}}}_0$ defined previously, for a regression over Years 1–20).

Nonlinearity in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ means that the linear fit parameters ${\bar{\lambda }}_{\mathcal{N}}$ and ${\bar{\mathcal{N}}}_I$ (slope and intercept) will change with the number and range of years used in the regression. Nevertheless, Equation 10 is a useful starting point to examine factors controlling climate sensitivity. We see that sensitivity increases both as ${\bar{\mathcal{N}}}_I$ increases, and as the magnitude of ${\bar{\lambda }}_{\mathcal{N}}$ decreases.

Shortwave and Longwave Contributions to Feedback and Initial Forcing

Figure 2 shows net shortwave and longwave TOM radiation fluxes, $\bar{\mathcal{S}}$ and $\bar{\mathcal{L}}$, as functions of $\bar{T}$ for CESM1‐4xCO2 (black) and CESM2‐4xCO2 (red). Figure 2 also shows equilibrium conditions for the piCTL experiments, in which $\bar{\mathcal{S}}$ and $\bar{\mathcal{L}}$ are within 0.1 W m⁻² of each other. Tables 3 and 4 give values of ${\bar{\mathcal{N}}}_0$, $\mathrm{\Delta }{\mathcal{L}}_0$, and $\mathrm{\Delta }{\mathcal{S}}_0$ as well as feedback parameters (slopes) ${\bar{\lambda }}_{\mathcal{N}}$, ${\bar{\lambda }}_{\mathcal{S}}$, and ${\bar{\lambda }}_{\mathcal{L}}$.

2 Figure. Annual‐mean, global top‐of‐atmosphere net shortwave S‾ and longwave L‾ radiative fluxes as functions of annual‐mean, global‐mean surface temperature T‾ for CESM1 (black) and CESM2 (red). Filled circles show annual‐mean S‾ for 4xCO2 experiments, and filled triangles show L‾. Large circles with error bars (2σ) show equilibrated multiyear means of S‾ and L‾ as functions of T‾ from the corresponding preindustrial control runs (piCTLs) for each model. Note that in the piCTLs, multiyear means of S‾ and L‾ are within 0.1 W m−2 of each other. Long dashes show extrapolations of linear regression fits to S‾ for Years 100–800 for CESM1‐4xCO2 extrapolation (black dashed line) and CESM2‐4xCO2 (red dashed line). Dotted lines show linear fits for Years 1–20. Slopes λ‾S of these lines are given in Table 4.

4 TableGlobal Feedback Parameters for Shortwave Flux ${\bar{\lambda }}_{\mathcal{S}}$, Longwave Flux ${\bar{\lambda }}_{\mathcal{L}}$, and Net Radiative Imbalance ${\bar{\lambda }}_{\mathcal{N}}$ for CESM1‐4xCO2 and CESM2‐4xCO2

Years	${\bar{\lambda }}_{\mathcal{S}}$ (W m−2 K−1)	${\bar{\lambda }}_{\mathcal{L}}$ (W m−2 K−1)	${\bar{\lambda }}_{\mathcal{N}}$ (W m−2 K−1)
CESM1‐4xCO2
1–20	0.99 (0.08)	2.05 (0.04)	−1.06 (0.09)
100–800	1.23 (0.02)	1.82 (0.01)	−0.59 (0.02)
100–800 (decadal averaging)	1.32 (0.02)	1.81 (0.02)	−0.49 (0.02)
CESM2‐4xCO2
1–20	0.87 (0.06)	2.01 (0.03)	−1.15 (0.07)
100–800	1.50 (0.01)	1.86 (0.01)	−0.36 (0.01)

Note. Note that since $\bar{\mathcal{N}}=\bar{\mathcal{S}}-\bar{\mathcal{L}}$, the fourth column is simply the difference of the second and third columns. Standard errors for the regression slopes are shown in parentheses. Results for regressions using decadally averaged quantities are shown for CESM1‐4xCO2. Decadal averaging has no effect on CESM2‐4xCO2 results.

When CO₂ is quadrupled, $\bar{\mathcal{L}}$ decreases rapidly by about 7.6 W m⁻² in both CESM1‐4xCO2 and CESM2‐4xCO2, while $\bar{\mathcal{S}}$ adjusts by +1 W m⁻² in CESM2‐4xCO2 and around −0.2 W m⁻² in CESM1‐4xCO2. This yields a larger net initial forcing ${\bar{\mathcal{N}}}_0$ of 8.6 W m⁻² in CESM2‐4xCO2 than 7.4 W m⁻² in CESM1‐4xCO2 (Table 3). So increased initial forcing, arising from a larger shortwave adjustment, is one component of the increased sensitivity of CESM2.

The overall behavior of $\bar{\mathcal{L}}(\bar{T})$ in Figure 2 is quite similar in CESM1‐4xCO2 and CESM2‐4xCO2, despite a small offset of about 2 W m⁻². We have already seen that in both experiments there is an initial adjustment in $\bar{\mathcal{L}}$ of around −7.6 W m⁻². Table 4 shows that the longwave feedback parameters ${\bar{\lambda }}_{\mathcal{L}}$ are also similar; initially around 2 W m⁻² K⁻¹ and becoming slightly smaller during Years 100–800, 1.82 W m⁻² K⁻¹ for CESM1‐4xCO2 and 1.86 W m⁻² K⁻¹ for CESM2‐4xCO2.

The long‐term value of ${\bar{\lambda }}_{\mathcal{S}}$ for CESM2‐4xCO2 is 1.50 W m⁻² K⁻¹, higher than in CESM1‐4xCO2 (1.23 W m⁻² K⁻¹). This produces the increased sensitivity in CESM2 by reducing the magnitude of long‐term ${\bar{\lambda }}_{\mathcal{N}}$ ( $={\bar{\lambda }}_{\mathcal{S}}-{\bar{\lambda }}_{\mathcal{L}}$) from −0.59 W m⁻² in CESM1‐4xCO2 to −0.36 W m⁻² in CESM2‐4xCO2 (Table 4), overwhelming the small increase in ${\bar{\lambda }}_{\mathcal{L}}$ from CESM1 to CESM2. Thus, both factors that can lead to increased iECS in CESM2, ${\bar{\mathcal{N}}}_0$ and ${\bar{\lambda }}_{\mathcal{N}}$, are modified through the shortwave component $\bar{\mathcal{S}}$. The stronger nonlinearities in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ for CESM2 also emerge from $\bar{\mathcal{S}}$.

We estimate the impact on ECS of the 1.2 W m⁻² increase in ${\bar{\mathcal{N}}}_0$ between CESM1 and CESM2 using the Year 100–800 linear fits shown in Figure 1a. The linear fit values of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ and $\mathrm{\Delta }\bar{T}$ at Year 100 are indicated in the figure. For CESM2‐4xCO2 we have $\mathrm{\Delta }\bar{T}(100)=6.58$ K and ${\bar{\mathcal{N}}}_{lin}(100)=2.55$ W m⁻². Using a slope  ${\bar{\lambda }}_{\mathcal{N}}=-0.36$ W m⁻² K⁻¹ (Table 4), we calculate an equilibrium warming of 6.58 $+\frac{2.55}{0.36}\approx$13.7 K, that is, the x intercept of the red dashed line in Figure 1a. Lowering ${\bar{\mathcal{N}}}_{lin}(100)$ by 1.2 to 1.35 W m⁻² would yield an adjusted equilibrium warming of 6.58 $+\frac{1.35}{0.36}\approx$10.3 K, corresponding to a climate sensitivity of 5.15 K. So with ${\bar{\lambda }}_{\mathcal{N}}$ as given in Table 4, reducing ${\bar{\mathcal{N}}}_0$ for CESM2‐4xCO2 to its value in CESM1‐4xCO2 gives a substantial reduction in ECS but would still yield a sensitivity larger than 5 K.

For comparison, we calculate the ECS that CESM2 would have if the long‐term, net radiative feedback in CESM2‐4xCO2 had the same value as in CESM1‐4xCO2, that is, −0.59 W m⁻² K⁻¹ instead of −0.36 W m⁻² K⁻¹. From Figure 1a, we see a slope change in $\bar{\mathcal{N}}$ near $\mathrm{\Delta }\bar{T}=5$ K for both CESM1‐4xCO2 and CESM2‐4xCO2. The value of the linear regression fit to $\bar{\mathcal{N}}$ at $\mathrm{\Delta }\bar{T}=5$ K for CESM2‐4xCO2 is 3.1 W m⁻². If the slope of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ in CESM2‐4xCO2 were steepened to −0.59 W m⁻² K⁻¹ at this point, there would be additional warming of about $\frac{3.1}{0.59}\approx$5.3 K, yielding a total warming of 10.3 K, again corresponding to an ECS of around 5.15 K.

We have seen that increased initial shortwave radiative forcing and increased shortwave radiation feedbacks play comparable roles in the greater sensitivity of CESM2‐4xCO2 relative to CESM1‐4xCO2. An important question which we cannot address here is how these two components of the sensitivity would change in an abrupt 2xCO2 ESM experiment. However, experiments with the CESM2‐SOM configuration (section 5) suggest that feedback strength ${\bar{\lambda }}_{\mathcal{N}}$ in 2xCO2 and 4xCO2 experiments is similar, while there is nonlinearity in ${\bar{\mathcal{N}}}_0$.

Impact of Subdecadal Variability

Table 4 shows that decadal averaging has a small but appreciable impact on regression estimates of shortwave feedback in CESM1‐4xCO2. We believe this impact arises because subdecadal variations in $\mathcal{S}$ and T are negatively correlated over large areas of the tropical ocean in CESM1‐4xCO2 (not shown). The origin of these correlated variations is not completely understood but is likely related to strong ENSO in the LME version of CESM1 (Otto‐Bliesner et al., 2016; Stevenson et al., 2016). It is worth emphasizing that the difference between the estimates of ${\bar{\lambda }}_{\mathcal{S}}$ using decadal and annual averages is not a reflection of statistical uncertainty in either estimate.

We will not address high‐frequency variability further in this study. However, it is clear that this variability could have impacts on calculations of iECS from 4xCO2 experiments in some models.

Global Distribution of Feedbacks

Figure 3 shows maps of long‐term linear regression slopes of quantities involved in shortwave radiative feedback for Years 100–800 in CESM1‐4xCO2 and CESM2‐4xCO2. The annual‐mean fields of $\mathcal{S}$ and T have been smoothed in time with a running 10‐yr window, and in space with an 8° rectangular latitude‐longitude window, before performing the linear regression.

3 Figure. Slopes from linear regressions over Years 100–800 of CESM1‐4xCO2 (a, c, e, g) and CESM2‐4xCO2 (b, d, f, h) as functions of latitude and longitude: (a, b) A(x, y), local warming amplification factor from regression of local temperature versus global‐mean temperature T‾; (c, d) λS(x,y), local shortwave feedback from regression of shortwave radiation S versus temperature; (e, f) slope of local shortwave flux versus global‐mean temperature T‾; (g, h) product of A(x, y) and λS(x,y).

Figures 3a and 3b show regression slopes of T(x, y) versus $\bar{T}$. This is a local amplification factor for warming, which we denote by A(x, y) and is the gridpoint analog of A_k in Equation 7. Both CESM1‐4xCO2 and CESM2‐4xCO2 exhibit polar amplification in both northern and southern high latitudes, although relative warming in the Arctic is much stronger in CESM1. This is likely related to differences in sea ice, as will be shown below. With the exception of the Arctic in CESM1‐4xCO2, warming in both models is generally stronger in the Southern Hemisphere (SH) than in the north. Both models show weak warming A(x, y) < 0.5 in the northwest Atlantic, accompanied by similarly weak warming in the northwest Pacific in CESM1‐4xCO2. An El Niño‐like warming pattern is present in the equatorial and southeastern Pacific.

Figures 3c and 3d show regression slopes of $\mathcal{S}(x,y)$ versus local T(x, y). This is the local feedback between shortwave radiation and surface temperature, which we denote by ${\lambda }_{\mathcal{S}}(x,y)$ and is the gridpoint analog of ${\lambda }_{\mathcal{S};k}$ in Equation 7. Despite the substantial changes in boundary layer and cloud physics parameterizations between CESM1 and CESM2, there are rough similarities in ${\lambda }_{\mathcal{S}}(x,y)$, particularly where low clouds are likely to control the shortwave response. Positive slopes with values between 3 and 5 W m⁻² K⁻¹ are evident in the midlatitude storm tracks (NH and SH) and stratus/stratocumulus regions of both models. This suggests the presence of positive low‐cloud SW feedbacks (i.e., thinner low clouds with higher T) of comparable magnitudes in both models. Shortwave feedbacks over the Southern Ocean stormtracks, however, are stronger in CESM2‐4xCO2 by about 1 W m⁻² K⁻¹. Also, CESM2‐4xCO2 has a large ${\lambda }_{\mathcal{S}}(x,y)$ in the deep convective region over the western tropical Pacific, whereas this strong positive feedback (>5 W m⁻² K⁻¹) is absent in CESM1.

Figures 3e and 3f show regression slopes of $\mathcal{S}(x,y)$ versus $\bar{T}$ in CESM1‐4xCO2 and CESM2‐4xCO2. Although the direct physical meaning of this regression quantity is unclear, this quantity is of interest since simple area integrals give the global feedback ${\bar{\lambda }}_{\mathcal{S}}$ (Andrews et al., 2015). Figures 3g and 3h show ${\lambda }_{\mathcal{S}}(x,y)\times A(x,y)$. This quantity should be close to the regression slopes of $\mathcal{S}$ versus $\bar{T}$ shown in Figures 3e and 3f, and this is in fact the case. The agreement between Figures 3e and 3f and Figures 3g and 3h argues that regional feedbacks on decadal time scales and ∼8° spatial scales can be accurately decomposed according to Equations 6 and 7.

In addition, comparison of Figures 3e and 3f or Figures 3g and 3h with Figures 3c and 3d highlights the role of regional warming in modulating the regional contributions to the global shortwave feedback. In particular, the relatively strong warming of the Southern Ocean amplifies its contribution to the global shortwave feedback, while weaker warming in the tropics reduces the contribution relative to its size in Figures 3c and 3d. Nevertheless, as will be shown in the next section the large size of the tropical ocean regions results in a large contribution to the global feedback, comparable to that of the Southern Ocean.

Figure 4 shows the distribution of shortwave, longwave and net radiation feedbacks in CESM2 and their evolution from CESM1. There is substantial compensation between shortwave and longwave feedbacks in the tropics (Figures 4a and 4c). In particular we see that the large positive shortwave feedbacks in the western tropical Pacific in CESM2 are accompanied by similarly large longwave feedbacks. Similar compensation exists for the negative shortwave and longwave feedbacks in the central tropical Pacific.

4 Figure. Slopes from linear regressions over Years 100–800 of CESM2‐4xCO2 (a, c, e) and differences from CESM1‐4xCO2 (b, d, f) as functions of latitude and longitude: (a) λS(x,y), local shortwave feedback; (b) ΔλL(x,y), difference from shortwave feedback in CESM1; (c, d) same as (a) and (b) except for longwave feedbacks; (e, f) same as (a) and (b) except for net radiation feedbacks.

Net radiation feedbacks in CESM2 are shown in Figure 4e. Not surprisingly over middle‐ and high‐latitude oceans the net feedback follows the shortwave feedback, with weak longwave compensation. This simply reflects the dominant role of low‐cloud feedbacks in these regions. In the tropics and subtropics a more complicated picture exists. In the northwestern tropical Pacific, longwave feedback exceeds shortwave feedback leading to more negative net feedback. By contrast in the eastern Pacific, positive net radiation feedback is present. In section 4.4.3 we will show that this is a consequence of overlapping middle‐ and high‐cloud feedbacks related to the pattern of upper‐level warming in the model, and that may ulitmately tied to changes in the implementation of the deep convection scheme in CESM2.

Figure 4f shows the change in net radiation feedback $\mathrm{\Delta }{\lambda }_{\mathcal{N}}(x,y)$ between CESM1 and CESM2. This quantity has a positive global mean reflecting the reduced negative global‐mean feedback in CESM2 (Table 4), which along with increased initial forcing has caused the increased ECS in CESM2. The pattern of $\mathrm{\Delta }{\lambda }_{\mathcal{N}}(x,y)$ is not as well defined as those of its shortwave and longwave components in Figures 4b and 4d. However, the uncompensated increase in shortwave feedback over the Southern Ocean is reflected in $\mathrm{\Delta }{\lambda }_{\mathcal{N}}(x,y)$. In the tropics we see generally positive values of $\mathrm{\Delta }{\lambda }_{\mathcal{N}}(x,y)$ except over the west Pacific warm pool and, interestingly, over subtropical low‐cloud regions indicating that these regions are acting to reduce ECS in CESM2 with respect to CESM1. Another region with negative $\mathrm{\Delta }{\lambda }_{\mathcal{N}}(x,y)$ is the Arctic, driven by reduced shortwave feedbacks there. In the following sections we will quantify regional contributions to radiation feedbacks and attempt to tie them to physical processes in CESM1 and CESM2.

Regional Feedbacks and Their Contribution to Global Climate Sensitivity

Figure 5 shows regions that have been selected to examine regional radiation feedbacks: (a) Arctic Ocean; (b) North Atlantic and North Pacific north of 30°N (NAtlPac); (c) tropical oceans between 30°S and 30°N (Trop_Ocn); (d) midlatitude Southern Ocean between 30°S and 60°S (SHml_Ocn); (e) high‐latitude Southern Ocean south of 60°S (SHhl_Ocn); (f) land north of 30°N (NH_Land); (g) tropical land between 30°S and 30°N (Trop_Land); (h) land south of 30°S (SH_Land); and (i) global. The fractional global area of each region is shown in the panels. The North Atlantic/North Pacific and midlatitude Southern Ocean regions (Figures 5b and 5d) are chosen to characterize generally ice‐free midlatitude oceans, while Arctic and high‐latitude Southern Ocean regions (Figures 5a and 5e) characterize high‐latitude oceans in which sea ice feedbacks may play a role.

5 Figure. Regions used for feedback analyses: (a) Arctic Ocean; (b) North Atlantic and North Pacific north of 30°N (NAtlPac); (c) ocean between 30°S and 30°N (Trop_Ocn); (d) midlatitude Southern Ocean between 30°S and 60°S (SHml_Ocn); (e) high‐latitude Southern Ocean south of 60°S (SHhl_Ocn); (f) land north of 30°N (NH_Land); (g) land between 30°S and 30°N (Trop_Land); (h) land south of 30°S (SH_Land); and (i) global. Approximate fractional area of regions are given in each panel.

Despite the obvious longitude dependence of tropical feedbacks in Figures 3 and 4, we have chosen to keep a single tropical ocean region to simplify the presentation of the regional analysis and also because the tropics represent a reasonably distinct self‐contained dynamical regime. The zonal structure of tropical feedbacks will be described in more detail in section 4.4.3.

Figure 6 shows time series of T in the analysis regions. After a rapid initial warming, there is a pause in warming, or even cooling, for about 100 yr in the Arctic, North Atlantic/North Pacific and northern land regions (Figures 6a, 6b, and 6f) in both CESM1‐4xCO2 and CESM2‐4xCO2; however, this feature is stronger in CESM1. In CESM2, rapid warming in the tropics (Figures 6c and 6g) and SH (Figures 6d, 6e, and 6h) overwhelms the effect of northern middle to high latitudes in the global mean (Figure 6i). In CESM1, the northern ocean cooling is strong enough to produce the noticeable hiatus or pause in global warming from around Year 20 to Year 150 seen here (Figure 6i) and in Figures 1c and 1d. Notably, the corresponding regional time series in CESM1b‐4xCO2 (not shown) and global time series (shown in Figure 1d, gray line) are nearly identical to those from CESM1‐4xCO2, despite different atmosphere resolution and ocean initialization. This consistency suggests that the NH land/ocean behavior shown in Figure 6 is a robust response of CESM1 to 4xCO2 forcing scenarios, not a result of internal variability. The complex response of northern high latitudes in the 4xCO2 scenario is of great interest, but will not be explored in this study. The figure also highlights the greater subdecadal, interannual variability in CESM1, which is particularly evident in the tropics (Figures 6c and 6g).

6 Figure. Regional‐mean time series of surface temperature T for regions in Figure 5. Black shows CESM1‐4xCO2, and red shows CESM2‐4xCO2. Solid lines show annual means subjected to a running 10‐yr mean. Symbols show annual means.

Figure 7 shows scatterplots of decadally averaged annual‐mean ${\mathcal{S}}_k$ versus T_k in CESM1‐4xCO2 and CESM2‐4xCO2 for the regions in Figure 5. The figure shows that reasonable functional relationships exist between decadally averaged ${\mathcal{S}}_k$ and T_k in all regions, that is, there is only small scatter about a curve fitted through the points in the scatterplot. Similar results are obtained for longwave radiation (not shown). The figure highlights the regional variations in ${\mathcal{S}}_k(T)$ as well as the large absolute differences between shortwave fluxes in CESM1 and CESM2. Regional‐mean differences of over 10 W m⁻² are present, with $\mathcal{S}$ in CESM1 generally lower (stronger shortwave CRE) than in CESM2 in the tropics, and $\mathcal{S}$ in CESM1 higher than in CESM2 in midlatitudes. The behavior of ${\mathcal{S}}_k$ in tropical ocean (Figure 7c) is especially noteworthy showing clearly stronger feedback in CESM2 (consistent with the patterns in Figures 3c and 3d), even though ${\mathcal{S}}_k$ is higher, which indicates thinner clouds. This is noteworthy because feedbacks are often assumed to increase as clouds become thicker.

7 Figure. Regional‐mean, net shortwave radiation Sk as a function of mean surface temperature Tk in CESM1‐4xCO2 (black circles) and CESM2‐4xCO2 (red circles) for regions in Figure 5. Larger circles show decadal averages for entire 4xCO2 simulations. Smaller circles show annual means for Years 1–20.

Regional Linear Regression Analyses

To quantify the contributions of the regions in Figures 5a–5h to global feedbacks between radiative fluxes and T, we perform linear regressions of ${\mathcal{S}}_k$, ${\mathcal{L}}_k$ and ${\mathcal{N}}_k$ versus T_k to determine regional feedback parameters ${\lambda }_{\mathcal{S},L,N;k}$, as well as regressions of T_k versus $\bar{T}$ to determine regional warming amplification factors A_k. These regression parameters are then used in Equation 7. We perform regressions over two periods: Years 1–20, to characterize the initial adjustment; and Years 100–800, to characterize the long‐term slow adjustment. As indicated in section 3.1, model results for Years 100–800 are decadally averaged before linear regression is performed. The subdecadal variability present in the tropics of CESM1 can be expected to affect the regressions for Years 1–20. We note this possibility but will not attempt to address it further in this analysis.

Figure 8 shows regional contributions to radiation feedbacks calculated using Equation 7 for shortwave, longwave and net TOM radiation fluxes $\mathcal{S}$, $\mathcal{L}$, and $\mathcal{N}$. The figure quantifies how much each analysis region contributes to the total global feedback parameters shown in Table 4. The bars in Positions 1–8 of each panel show the complete summands $a_k{\lambda }_{\mathcal{S},L,N;k}$ used in Equation 7 for the regions indicated. CESM1‐4xCO2 is shown by the black bars, and CESM2‐4xCO2 by the red bars, with differences shown in green. The bars in Position 9 show the direct sum over the eight regions, while Position 10 shows independent regressions of global means $\bar{\mathcal{S}}$, $\bar{\mathcal{L}}$, and $\bar{\mathcal{N}}$ versus $\bar{T}$. The close agreement between the direct sums in Position 9 and the independent regression estimates in Position 10 validates the regional decomposition in Equation 7.

8 Figure. Regional contributions to shortwave, longwave, and net radiation feedback parameters λ‾S, λ‾L, and λ‾N computed using Equation 7. Left panels(a, c, e) show results for the early phase of the 4xCO2 runs (Years 1–20), and right panels (d, b, f) show results for the later “slow‐adjustment” phase (Years 100–800). (a, b) Complete regional shortwave contributions akAkλS;k. (c, d) Complete regional longwave contributions akAkλL;k. (e, f) Complete regional contributions to net radiation feedback λN;k Black bars indicate results for CESM1‐4xCO2, red bars indicate CESM2‐4xCO2, and green bars show differences between CESM2 and CESM1. Each panel shows 10 pairs of bars. Positions 1–8 show quantities for the regions shown in Figure 5. The bars in Position 9 show direct sums over the eight terms shown to the left, while Position 10 shows independent regressions of global means S‾, L‾, and N‾ versus T‾.

The nonlinear time behavior in radiation feedbacks can be visually evaluated by comparing the early regression period (Years 1–20, Figures 8a, 8c, and 8e) with the later period (Years 100–800, Figures 8b, 8d, and 8f). The largest regional contributions to the nonlinearity in shortwave feedback are from tropical and midlatitude Southern Oceans (Figures 8a and 8b, Positions 3 and 4), accounting for almost all of the later‐period increase in global shortwave feedback in both models. In contrast, contributions to shortwave feedback from middle‐ and high‐latitude Northern Hemisphere (NH) and tropical land (Positions 6 and 7) decrease by 0.1–0.3 W m⁻² K⁻¹ between the early and later periods. Regional longwave feedbacks in both models also change with time. Notable decreases in later‐period longwave feedbacks appear over tropical oceans (Position 3) and tropical and northern land (Positions 6 and 7) driving the 0.15‐0.2 W m⁻² K⁻¹ decreases in global longwave feedbacks shown in Table 4. Interestingly, longwave feedbacks over the Southern Ocean (Position 4) in both models increase with time. This is primarily a clear‐sky effect (not shown), but the ultimate cause of this behavior is not yet understood.

The increased ECS in CESM2 is driven by the increased initial forcing examined in section 4.1 and by changes to the long‐term (later period) radiative feedbacks shown in Figures 8b, 8d, and 8f. The evolution of net radiation feedbacks from CESM1 to CESM2 is summarized by the green bars in Figure 8f. The fact that the global net radiation feedbacks (Positions 9 and 10) have become more positive (weaker) in CESM2 simply restates the content of Table 4 (Column 4) in graphical form. The regional contributions to the global net radiation feedbacks can be read directly from the bars in Positions 1–8. Tropical oceans and Southern Ocean play the largest role in reducing the strength of net radiation feedbacks and thereby increasing the ECS of CESM2. Examining Figures 8b and 8d we see that the important regional changes to net radiative feedback arise from changes to shortwave feedbacks. In fact, the long‐term longwave feedback in CESM2 is slightly stronger in most regions, which would tend to reduce ECS.

The evolution of Arctic feedbacks (Position 1) is notable, because it suggests a possible buffering effect on ECS from rapid sea ice loss in CESM2‐4xCO2. In the following sections we will examine the behavior of feedbacks in the key regions identified by Figure 8; Southern Ocean, tropical oceans and Arctic, and also examine the role of cloud and surface processes in determining these feedbacks.

Cloud and Surface Processes

Figure 9 shows the regional breakdown of radiation feedbacks into all‐sky, cloudy (CRE) and clear‐sky components for CESM1‐4xCO2 and CESM2‐4xCO2 for Years 100–800 of the experiments. We focus on the slow adjustment because these feedbacks are ultimately responsible for determining the model climate sensitivity. Our initial analysis looks at CESM outputs of total (all‐sky) longwave and shortwave TOM radiation and longwave and shortwave cloud radiative forcing, which are then used to diagnose clear‐sky fluxes according to Equation 2. This gives a first impression of the role of cloud feedbacks. Shortwave cloud feedbacks are then further analyzed using the APRP approach.

9 Figure. Decomposition of radiation feedbacks for Years 100–800 in CESM1‐4xCO2 (a–c), CESM2‐4xCO2 (d–f), and differences (g–i) into all‐sky (black, red and green bars), cloud radiative effect (CRE, gray bars) and clear‐sky (blue bars) components by region as in Figure 8. First column (a, d, g) shows total regional contributions to global shortwave feedbacks. Second column (b, e, h) shows total regional contributions to global longwave feedbacks. The longwave CRE contribution has been multiplied by −1 so that bars for clear‐sky and CRE feedbacks are additive in the same sense as in the shortwave. Third column (c, f, i) shows contributions to net TOM radiation feedbacks. More negative values of net TOM radiation feedback correspond to reduced climate sensitivity. Thus, positive green bars in in Panel (i) indicate a regional contribution to increased climate sensitivity in CESM2.

In the shortwave (Figure 9a, 9d, and 9g) the large increase in feedback between CESM1 and CESM2 arises from the cloudy component (gray bars), with approximately equal contributions from tropical oceans and Southern Ocean (Figure 9g, Positions 3 and 4). In CESM1, clear‐sky shortwave feedbacks (blue bars) are large in the high‐latitude ocean regions (Arctic, Position 1, and high‐latitude Southern Ocean, Position 5), and over NH land, while in CESM2, clear‐sky feedbacks are noticeable only over middle‐ to high‐latitude land regions. Positive high‐latitude clear‐sky feedbacks over high‐latitude oceans produce a global positive clear‐sky shortwave feedback in CESM1 that is actually larger than the cloudy feedback. The positive clear‐sky feedbacks are accompanied and partially compensated by negative shortwave cloud feedbacks. The net shortwave feedback in these regions nevertheless remains positive in CESM1‐4xCO2 as highly reflective snow and ice surfaces disappear and are replaced by somewhat less reflective clouds (e.g., Frey et al., 2018).

Longwave feedbacks (Figures 9b, 9e, and 9h) have changed less in the evolution from CESM1 to CESM2. This is clearly seen by comparing the difference plots for shortwave and longwave feedbacks (Figures 9g and 9h). Clear‐sky longwave feedback is much larger than longwave CRE feedback in both models. Nevertheless, clear‐sky and CRE feedback both make comparable contributions to the small differences in longwave feedback between CESM1 and CESM2.

Regional contributions to the net radiation feedback are shown in Figures 9c, 9f, and 9i. Figure 9i, in particular, is a useful summary of the net radiation feedback changes that have occurred between CESM1 and CESM2. Changes to the net radiation feedbacks are clearly driven by changes in shortwave feedbacks (Figure  9g). Furthermore, all changes leading to increased climate sensitivity in CESM2 (positive sign in Figure 9i) arise in CRE feedbacks (gray bars). In high‐latitude ocean regions, increased CRE feedback in CESM2 is opposed by clear‐sky feedback (blue bars). Finally, it is worth noting that the increased tropical ocean shortwave feedback in CESM2 is not compensated by longwave feedbacks (Figures 9h and 9i).

The results in the section and the previous section have highlighted the role of Southern Ocean and tropical ocean shortwave cloud feedbacks in changing the net radiation feedback and increasing the ECS of CESM2 compared to that in CESM1. These results are roughly consistent with SST+4 K results discussed by Gettelman, Hannay, et al. (2019), although the contribution of tropical shortwave feedbacks to increasing CESM2's ECS appears to be larger in the 4xCO2 experiments examined here. In addition the regional results shown here point to interesting differences in Arctic and high‐latitude Southern Ocean feedbacks between CESM1‐4xCO2 and CESM2‐4xCO2.

Sea Ice Evolution

Figure 10 shows sea ice concentrations and surface albedo (calculated from model shortwave fluxes at the surface) in the Arctic and high‐latitude Southern Oceans in CESM1‐4xCO2 and CESM2‐4xCO2. Sea ice concentrations decrease rapidly in CESM2‐4xCO2 with little sea ice remaining in either high‐latitude ocean region after Year 200. The effective surface albedo in these regions is then essentially constant, explaining the lack of long‐term clear‐sky shortwave feedback in CESM2‐4xCO2. Sea ice and surface albedo in CESM1‐4xCO2 decrease much more slowly, especially in the Arctic, and remain at appreciable levels throughout the 800 yr of the experiment. This explains the presence of the large, long‐term, clear‐sky shortwave feedbacks seen for CESM1 in Figure 9.

10 Figure. (a) Annual‐mean sea ice fraction as a function of time for Arctic and high‐latitude Southern Oceans in CESM1‐4xCO2 (black) and CESM2‐4xCO2 (red). (b) As in (a) except for surface albedos as functions of time. Dashed lines show fraction and surface albedo in the Arctic Ocean (Figure 5a), and solid lines show fraction and surface albedo in the high‐latitude Southern Ocean (Figure 5e).

Figures 11a and 11b show regional‐mean cloud condensates as functions of surface temperature in the Arctic and high‐latitude Southern Oceans. As sea ice decreases in CESM1 (Figure 10), cloud condensate amounts increase with T throughout the experiment, contributing to the negative shortwave CRE feedback obtained for these regions in CESM1 (Figure 9a). In CESM2 we see an initial increase in condensate amounts in high‐latitude oceans, but during Years 100–800 condensate amounts become nearly constant, consistent with the lack of long‐term SW CRE feedbacks over high‐latitude oceans in Figure 9d.

11 Figure. Regional‐mean, in‐cloud condensate paths (IWP∗ and LWP∗, Equation 3) in g m−2 as functions of regional‐mean Tk in CESM1‐4xCO2 (black) and CESM2‐4xCO2 (red): (a) Arctic Ocean; (b) high‐latitude Southern Ocean; and (c) midlatitude Southern Ocean. Circles show cloud liquid water path LWP∗. Triangles show cloud ice water path IWP∗. Gray lines show linear fits over Years 100–800. Regression slopes λLWP∗;k and λIWP∗;k for these fits are given in upper right corner of each panel.

The sea ice behavior shown has an interesting implication for ECS in CESM2. The rapid loss of sea ice in CESM2 removes a large positive shortwave feedback from the climate system early in 4xCO2 experiment. If preindustrial sea ice in CESM2 were thicker it is possible it would continue to produce strong shortwave feedbacks longer into the 4xCO2 simulation, leading to even higher ECS in CESM2.

APRP Analysis

We use the APRP approach of Taylor et al. (2007) to further decompose shortwave radiation feedbacks into components related to specific physical processes. Figure 12 compares shortwave CRE feedbacks with respect to T(x, y), that is, ${\lambda }_{{\mathcal{S}}_{cld}}(x,y)$ over Years 100–800 with the quantities [Image Omitted. See PDF] [Image Omitted. See PDF] where $\mathcal{A}$, and γ_cld are APRP reconstructions of the planetary albedo and cloud scattering (Equation 8); c is total cloud amount used in the APRP calculation; and ${\mathcal{S}}^{\downarrow }$ is the incoming solar radiation at TOM. Partial derivatives are evaluated using the analytical expressions for $\mathcal{A}$ in Taylor et al. (2007) (their Equations 7, 13, 14, and 15) employing the Year 100–800 average values for all parameters in the evaluation. The feedback parameters λ_c and ${\lambda }_{{\gamma }_{cld}}$ are determined from linear regressions of c and γ_cld versus T(x, y) over Years 100–800.

12 Figure. Cloud‐related shortwave feedbacks as functions of latitude and longitude over Years 100–800 of CESM1‐4xCO2 and CESM2‐4xCO2: (a, b) Linear regression slopes for shortwave CRE Scld versus T, that is, λScld. (c, d) Cloud scattering contribution Λγcld (Equation 11b) to shortwave feedback. (e, f) Cloud amount contribution Λc (Equation 11a) to shortwave feedback. (g, h) Sum of Λγcld and Λc. Left column (a, c, e, g) shows results for CESM1‐4xCO2, and right column (b, d, f, h) shows results for CESM2‐4xCO2.

The quantities Λ_c and ${\mathrm{\Lambda }}_{{\gamma }_{cld}}$ are the dominant cloud‐related contributions to the shortwave feedback. Comparing Figures 12a and 12b with Figures 12g and 12h we see that the sum of Λ_c and ${\mathrm{\Lambda }}_{{\gamma }_{cld}}$ is very close to the shortwave CRE feedback (and to the all‐sky shortwave feedbacks in Figures 3c and 3d away from high latitudes). The individual components represent separate feedbacks associated with cloud scattering properties ( ${\mathrm{\Lambda }}_{{\gamma }_{cld}}$, Figures 12c and 12d) and cloud amount (Λ_c Figures 12e and 12f). Away from the tropics, these two components of the feedback have comparable magnitudes (1 to 2 W m⁻² K⁻¹) in both models. The cloud amount feedback is slightly more positive in CESM2 (Figure 12f) than in CESM1 (Figure 12e). In particular, Λ_c over the midlatitude Southern Ocean is similar in CESM1 and CESM2.

However, pronounced differences between CESM1 and CESM2 appear in ${\mathrm{\Lambda }}_{{\gamma }_{cld}}$, the cloud scattering component of the shortwave feedback (Figures 12c and 12d). Strong scattering feedbacks ∼4 W m⁻² K⁻¹ are present in CESM2 in the tropics, which are the main contribution to the stronger overall tropical ocean shortwave feedback noted in Figures 8 and 9 for CESM2. Over the midlatitude Southern Ocean we also see larger values of ${\mathrm{\Lambda }}_{{\gamma }_{cld}}$ in CESM2 which produce most of the increase in overall shortwave feedback there compared to CESM1.

The main conclusion of Figure 12 is that cloud scattering feedback explains more of the increased shortwave feedback in CESM2 than cloud amount feedback. Frey and Kay (2018) found similar increases in scattering feedback and climate sensitivity in CESM1 when they perturbed the model microphysics to increase the amount of supercooled liquid present in clouds. They discuss the possible role of phase feedbacks in suppressing Southern Ocean shortwave feedbacks in the default CESM1, that is, as ice cloud is replaced by more reflective liquid in a warming climate, cloud albedo increases. Gettelman, Hannay, et al. (2019) also cite reduced phase feedbacks as a possible contributor to increased ECS in CESM2.

Figure 11c shows average in‐cloud liquid and ice phase condensate paths (IWP^∗ and LWP^∗, Equation 3) over the midlatitude Southern Ocean. There is strong long‐term decrease of LWP^∗ with T in CESM2 compared to that in CESM1, coupled with a weak increase in IWP^∗. In CESM1, both LWP^∗ and IWP^∗ decrease with T in the long term, although a clear initial bump in LWP^∗ occurs. For Years 100–800, ${\lambda }_{LW{P}^{\ast }}=-1.67$ W m⁻² K⁻¹ in CESM1‐4xCO2, more than double ${\lambda }_{IW{P}^{\ast }}=-0.72$ W m⁻² K⁻¹, while in CESM2‐4xCO2 ${\lambda }_{LW{P}^{\ast }}=-4.47$ W m⁻² K⁻¹. We cannot quantify how much of the increased SW feedback in CESM2‐4xCO2 is due simply to the stronger loss of total condensate with T, and how much is due to the presence of negative phase feedback in CESM1‐4xCO2. However, the large decrease in CESM2 Southern Ocean total condensate compared to that of total condensate in CESM1 would suggest that total condensate loss is the dominant process. Quantitative analysis is left for a future study.

In the tropics we have an incomplete understanding of how changes to model physics changes between CESM1 and CESM2 have contributed to the evolving feedbacks in the model. Below we present more detail on the evolution of tropical responses to warming between CESM1 and CESM2 that we hope may point the way to a more physical understanding of these changes.

Tropical Feedbacks

Changing tropical radiation feedbacks make an important contribution to the increased ECS of CESM2. Figure 8f shows that tropical oceans contribute about half of the reduction in net feedback strength that occurs between CESM1 and CESM2. This is an area integral over a complex pattern of local feedbacks (Figure 4f) encompassing both deep convective and shallow coastal stratocumulus clouds. Generally speaking tropical stratocumulus in CESM2 are associated with stronger negative ${\lambda }_{\mathcal{N}}$ (i.e., tend to reduce ECS) compared to CESM1. This is evident in Figure 4f off the western coasts of North and South America, as well as off the Namibian coast and northwestern Australia.

However, net radiation feedbacks in tropical convective regions change in complicated ways between CESM1 and CESM2. Over the western Pacific warm pool there are dramatic increases in both longwave and shortwave feedbacks in CESM2 (Figures 4b and 4d) leading to a negative net feedback change $\mathrm{\Delta }{\lambda }_{\mathcal{N}}$. By contrast, over the western Atlantic we see a large reduction in longwave feedback in CESM2 which drives a large $\mathrm{\Delta }{\lambda }_{\mathcal{N}}>2$ W m⁻² K⁻¹ in that region. Over the remainder of the tropics including the eastern ITCZs in the Pacific, and over the Indian Ocean, $\mathrm{\Delta }{\lambda }_{\mathcal{N}}$ has values between 0.5 and 2.0 W m⁻² K⁻¹. A further breakdown of tropical ocean radiation feedbacks into contributions from the Indian, Pacific and Atlantic oceans (not shown) reveals that the effect of the tropical Atlantic on net radiation feedback is comparable to that of the Pacific.

Figure 13 shows cross tropical cross sections of trends from Year 100 to Year 800 in temperature, cloud fraction and cloud condensates from CESM1‐4xCO2 and CESM2‐4xCO2. Trends in net radiation are also shown. All trends are normalized by the surface temperature so are directly comparable to linear regression results (e.g., Figure 4). In Figures 13a and 13b we see the signature of enhanced upper‐level warming (e.g., Santer et al., 2016) in both CESM1‐4xCO2 and CESM2‐4xCO2. Upper‐level warming has become stronger in CESM2 over the western central Pacific (longitudes 120–200°), while over the tropical Atlantic (longitudes 280–340°) it is somewhat stronger in CESM1. There is clear alignment of the strongest upper‐level warming with reductions in cloud fraction (Figures 13c and 13d) from 500 to 100 hPa.

13 Figure. Evolution of tropical cloud response to warming between CESM1 and CESM2: (a) Warming amplification as a function of longitude and pressure in CESM1 averaged over ocean between 15°S and 15°N. Calculated as the change in mean T(x, p) for Years 101–200 and Years 701–800, year divided by the change in Ts(x); (b) same as (a) except for CESM2; (c, d) same as (a, b) except for cloud fraction; (e, f) same as (a, b) except for total cloud condensate; and (g, h) trends in N(x,y) as functions of longitude and latitude. Lines in (g) and (h) show latitudes 15°S and 15°N. Arrows in (e) indicate large positive low‐level cloud condensate trends in CESM1.

Changes in middle‐ and upper‐tropospheric cloud condensate (Figures 13e and 13f) are also clearly related to tropospheric warming, but there are significant differences in the vertical structure of these changes between CESM1 and CESM2. CESM2 (Figure 13f) shows large condensate reductions in a layer between 500 and 400 hPa with weaker increases in layer between 400 and 300 hPa. In CESM1 (Figure 13e) the strongest condensate loss appears between 300 and 200 hPa. Changes to the implementation of the ZM deep convection scheme in CESM2 have resulted in reduced convective cloud top heights. However, whether these changes have resulted in the changes in Figures 13a–13f is not clear.

Net radiation $\mathcal{N}$ trends are shown in Figures 13g and 13h. There is clear alignment in both models between regions where the trend in net radiation ${\lambda }_{\mathcal{N}}$ is negative and the location of the strongest reductions in middle‐ and upper‐level tropospheric cloud quantities. This indicates that, locally, reductions to longwave CRE are overcoming reductions to the shortwave CRE induced by thinner middle‐ and upper‐level clouds. An added complication appears in CESM1 where two regions of pronounced low‐cloud condensate increase appear at longitude ∼110° and longitude ∼300°, indicated by the arrows in Figure 13e. These strengthening low clouds enhance the negative trend in $\mathcal{N}$ (Figure 13g). In CESM2 the weaker positive trends in ${\lambda }_{\mathcal{N}}$ elsewhere in the tropics appear to result from a combination of decreasing low clouds and increasing midtropospheric clouds, for example, longitudes 200–300° in CESM2 (Figure 13h). We surmise that this midlevel cloud increase not able to overcome reductions in shortwave CRE caused by the underlying low‐cloud loss, while at the same increasing longwave CRE, leading to the positive trend in net tropical radiation feedbacks over much of the eastern Tropical Pacific and Indian Oceans in CESM2 (Figure 8f). By contrast, in CESM1 midlevel cloud condensate increases are lower (500 hPa) and are overlain (300 hPa) by a layer of stronger condensate decrease.

Summary

Section 4.1 examined longwave and shortwave contributions to climate sensitivity in CESM1 and CESM2. The increased iECS in CESM2 was seen to originate predominantly from changes to shortwave radiation. As found for many CMIP6 models in Zelinka et al. (2020), increased iECS in CESM2 arises from both increased initial forcing ${\mathcal{N}}_0$ and weaker net radiation feedback ${\bar{\lambda }}_{\mathcal{N}}$. In CESM2 both of these factors arise through changes to the shortwave component

Section 4.2 showed global distributions of radiation feedbacks in 4xCO2 simulations and their evolution from CESM1 to CESM2. Regional contributions to global feedbacks were analyzed in section 4.3. We showed that the increased ECS in CESM2 is primarily driven by stronger shortwave feedbacks in the Southern Ocean and over tropical oceans. The Southern Ocean results are consistent with those for CESM shown by Gettelman, Hannay, et al. (2019) and also discussed by Zelinka et al. (2020) for the CMIP5 and CMIP6 ensembles.

In section 4.4 we examined the processes leading to the evolution of regional radiation feedbacks in CESM. Arctic and Antarctic sea ice disappears more quickly in CESM2 than in CESM1 leading to reduced shortwave feedbacks over high‐latitude ocean. We applied the APRP approach (Taylor et al., 2007) to examine shortwave feedbacks and found that scattering effects, that is, cloud condensate amounts and possibly phase, were more importatnt than cloud fraction feedbacks in determining the evolution of shortwave feedbacks from CESM1 to CESM2 in the tropics and Southern Ocean.

We also examined the evolution of tropical cross sections of cloud and temperature tendencies in CESM 4xCO2 runs (section 4.4.3). Tropical longwave and shortwave feedbacks in CESM exhibit complicated zonal variations, that when integrated over the tropics result in weaker (less negative) net radiation feedback in the tropics in CESM2. We attribute this to stronger increases in midlevel cloudiness in CESM2 away from the west Pacific warm pool. It is still unclear how the changing tropical cloud response is related to changes in the models' ZM deep scheme.

Comparison With Slab‐Ocean Experiments and Relation to ECS‐SOM

Experiments with a thermodynamic SOM have been proposed as a way of reducing the computation required to derive estimates of ECS (Bitz et al., 2012; Danabasoglu & Gent, 2009). SOM experiments approach radiative equilibrium within several decades compared to the several hundred years required for ESM simulations with a dynamic ocean. This approach has been used by several investigators to estimate ECS for various versions of CESM (Bitz et al., 2012; Gettelman et al., 2012; Gettelman, Hannay, et al. 2019).

Figure 14 shows $\bar{\mathcal{N}}$ versus $\mathrm{\Delta }\bar{T}$ for CESM1‐4xCO2 and CESM2‐4xCO2 in SOM and ESM configurations. The sparse density of points for the SOM‐4xCO2 runs (gray circles) is a consequence of their rapid equilibration. Nevertheless, there is remarkable overlap between $\bar{\mathcal{N}}$ versus $\mathrm{\Delta }\bar{T}$ in the SOM and ESM experiments despite the vastly different time scales with which warming occurs.

14 Figure. Top panels (a, b) show N‾, net annual‐mean global radiative imbalance at TOM, as a function of global‐mean surface temperature change ΔT‾ for fully coupled (ESM) and slab‐ocean (SOM) abrupt CO2 increase experiments: (a) CESM1. Gray circles show CESM1b‐4xCO2‐SOM, and black circles show CESM1‐4xCO2 (ESM); (b) CESM2. Gray circles show CESM2‐4xCO2‐SOM, red circles show CESM2‐4xCO2 (ESM), and gray triangles show CESM2‐2xCO2‐SOM. Larger gray circles in (a) and (b) show years in the SOM 4xCO2 experiments where ΔT‾ overlaps with that in the Year 100–800 range of the corresponding ESM experiments, that is, Years 5–15 of CESM1‐4xCO2‐SOM and Years 10–30 of CESM2‐4xCO2‐SOM. Bottom panels (c, d) show sea ice fraction as a function of regional‐mean surface temperature: (c) High‐latitude Southern Ocean; and (d) Arctic. Sea ice fraction in Years 100–800 in CESM1‐4xCO2 (ESM) and CESM2‐4xCO2 (ESM) is shown, along with Years 5–15 for CESM1b‐4xCO2‐SOM and Years 10–30 for CESM2‐4xCO2‐SOM.

The 4xCO2‐SOM experiments attain radiative equilibrium $\mathcal{N}\to$ 0 and surpass the warming realized in the corresponding ESM runs. Close inspection of the CESM2‐4xCO2‐SOM results in Figure 14b shows that for $\mathrm{\Delta }\bar{T}>11.5$ K there is an increase in the feedback strength ${\bar{\lambda }}_{\mathcal{N}}$, leading to a smaller equilibrium warming (∼12.6 K) than the equilibrium warming of 13.7 K predicted by extrapolating the slow‐adjustment behavior of CESM2‐4xCO2 (ESM) as discussed in section 4.1, and suggesting that feedbacks may change in the ESM even after 1,000 yr (Rugenstein et al., 2019). This warming in CESM2‐4xCO2‐SOM is leads to ECS‐SOM(4x) values of around 6.5 K (Table 1, see also Figure A1). Figure 14 also shows results for a 2xCO2‐SOM experiment with CESM2 (gray triangles). The scatter of points is large compared to the warming signal, but the overall shape of the $\bar{\mathcal{N}}$ versus $\mathrm{\Delta }\bar{T}$ relationship in the 2xCO2 experiment resembles that in the 4xCO2 experiments. Interestingly, values of ${\bar{\lambda }}_{\mathcal{N}}$ calculated over Years 1–100 of the 2xCO2‐SOM and 4xCO2‐SOM experiments are similar, −0.42 and −0.39 W m⁻² K⁻¹, respectively, suggesting that radiation feedbacks in CESM2 are not highly nonlinear with respect to CO₂. We note that sea ice evolution in the 2xCO2 and 4xCO2‐SOM experiments is different, with Arctic sea ice fraction in CESM2‐2xCO2‐SOM remaining close to 0.3 though the end of the experiment, while it disappears completely from CESM2‐4xCO2‐SOM.

We calculate an “ECS‐SOM(4x)” from these 4xCO2‐SOM runs as is done in the standard 2xCO2 setup to determine ECS‐SOM, except that we divide the equilibrium warming from the 4xCO2‐SOM runs by 2 (Appendix A). Table 1 gives ECS‐SOM(4x) values from our experiments compared to values of ECS based on 2xCO2‐SOM experiments. Minimal nonlinearity exists in ECS‐SOM estimates for CESM1, but moderate nonlinearity is present in CESM2, with ECS‐SOM(4x) about 1.15 times higher than ECS‐SOM based on 2xCO2 experiments. Plots of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ in Figure 14b suggest nonlinearity in initial forcings for CESM2‐4xCO2‐SOM and CESM2‐2xCO2‐SOM, with ${\mathcal{N}}_0\sim$ 8 W m⁻² in 4xCO2 compared with ∼3.5 W m⁻² in 2xCO2. Assuming similar long‐term slopes for $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$, these changes in ${\bar{\mathcal{N}}}_0$ would account for the nonlinearities in ECS‐SOM(4x) noted in Table 1.

To identify roughly comparable periods of long‐term adjustment in the 4xCO2‐SOM and ESM experiments we match $\mathrm{\Delta }\bar{T}$ in the SOM to the values obtained in Years 100–800 in the corresponding ESM experiments. These points are shown on the plots of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ (Figure 14) by larger symbols. For CESM1‐4xCO2‐SOM we identify Years 5–15 as the equivalent long‐term adjustment period, while for CESM2‐4xCO2‐SOM we identify Years 10–30. We recognize that this equivalence may miss important regional differences. Figures 14c and 14d show sea ice fraction in the high‐latitude Southern and Arctic Oceans as functions of T for Years 100–800 in CESM1‐4xCO2(ESM) and CESM2‐4xCO2(ESM), and the equivalent periods in CESM1‐4xCO2‐SOM and CESM2‐4xCO2‐SOM. Sea ice fractions in the SOM runs are significantly higher than in the ESM at similar values of regional‐mean T, particularly in the Arctic. Higher Arctic temperatures are reached in CESM2‐4xCO2‐SOM than in CESM2‐4xCO2(ESM) (Figure 14d). Table 5 shows radiation feedback parameters for CESM2‐4xCO2‐SOM and ESM experiments.

5 TableGlobal Feedback Parameters for Shortwave Flux ${\bar{\lambda }}_{\mathcal{S}}$, Longwave Flux ${\bar{\lambda }}_{\mathcal{L}}$, and Net Radiative Imbalance ${\bar{\lambda }}_{\mathcal{N}}$ for CESM2‐4xCO2 and CESM2‐4xCO2‐SOM

Years	${\bar{\lambda }}_{\mathcal{S}}$ (W m−2 K−1)	${\bar{\lambda }}_{\mathcal{L}}$ (W m−2 K−1)	${\bar{\lambda }}_{\mathcal{N}}$ (W m−2 K−1)
CESM2‐4xCO2
1–20	0.87 (0.06)	2.01 (0.03)	−1.15 (0.07)
100–800	1.50 (0.01)	1.86 (0.01)	−0.36 (0.01)
CESM2‐4xCO2‐SOM
1–5	0.79 (0.10)	2.11 (0.04)	−1.32 (0.10)
10–30	1.48 (0.03)	1.75 (0.02)	−0.28 (0.04)

Note. Note that since $\bar{\mathcal{N}}=\bar{\mathcal{S}}-\bar{\mathcal{L}}$, the fourth column is simply the difference of the second and third columns. Standard errors for the regression slopes are shown in parentheses.

Figure 15 shows maps of regression coefficients of T versus $\bar{T}$, that is, A(x, y); $\mathcal{S}$ versus T, that is, ${\lambda }_{\mathcal{S}}(x,y)$; and $\mathcal{L}$ versus T, that is, ${\lambda }_{\mathcal{L}}(x,y)$ for CESM2‐4xCO2‐SOM and CESM2‐4xCO2 (ESM). The regressions for CESM2‐4xCO2 (ESM) are performed over Years 100–800 and the corresponding period (Years 10–30) in the SOM experiment. The warming amplification factor A(x, y) shows large differences between SOM and ESM experiments, which are likely tied to lack of dynamically evolving heat transport in the SOM. The SOM (Figure 15a) exhibits a more hemispherically symmetrical distribution, with both northern and southern high latitudes having broad areas with A(x, y) > 1.75. In contrast, the ESM has values of A(x, y,) around 1.25 or below in northern high latitudes, but exceeding 2.5 over much of the Antarctic. Broad areas of the tropics and northern midlatitudes also warm less in the ESM, while southern midlatitudes warm more. The role of sea ice (Figures 14c and 14d) in the different pattern of polar amplification in the SOM and ESM is not yet understood.

15 Figure. Slopes from linear regressions as functions of latitude and longitude for CESM2‐4xCO2‐SOM (a, c, e) and CESM2‐4xCO2 (b, d, f): (a, b) A(x, y), the local warming amplification factor from regression of local temperature versus global‐mean temperature T‾; (c, d) λS(x,y), the local shortwave feedback from regression of shortwave radiation S versus local temperature; and (e, f) λL(x,y), the local longwave feedback from regression of shortwave radiation S versus local temperature. Regressions are performed over Years 10–30 for CESM2‐4xCO2‐SOM and Years 100–800 for CESM2‐4xCO2 (ESM).

Radiation flux feedbacks ${\lambda }_{\mathcal{S}}(x,y)$ and ${\lambda }_{\mathcal{L}}(x,y)$ shown in Figures 15c–15f are remarkably similar in the SOM and ESM experiments. Feedbacks across the tropical Pacific are somewhat more zonal in the ESM for both shortwave and longwave radiation. This is especially evident in shortwave feedbacks over the tropical eastern Pacific where strong positive feedbacks (>5 W m⁻² K⁻¹) appear in the SOM but not in the ESM. This could reflect stronger eastward shifting of tropical Pacific convection in the ESM, consistent with precipitation differences between the SOM and the ESM (not shown). Increasing high‐cloud associated with this shift masks the decrease in low cloud over the eastern Pacific associated with local warming, as reflected in the distribution of Λ_c shown in Figure 12f.

Overall, the close agreement between the final equilibrium global temperatures under 4xCO2 and in the behavior of $\bar{\mathcal{N}}(\bar{T})$ in SOM and ESM configurations is striking. It is perhaps even more striking that this agreement occurs despite significant regional and hemispheric differences in warming. It is tempting to seek an explanation based on energetic considerations. However, a convincing explanation has not yet been found by the authors. In any event, it appears that CESM's SOM configuration is capable of predicting the total global‐mean warming produced in much longer ESM 4xCO2 simulations. It is not clear whether the skillful performance of the CESM SOM is critically dependent on its design. SOM versions of other ESMs may not behave in a similar way (e.g., Boer & Yu, 2003).

Comparison With 1%CO2 Experiments and Relation to TCR

The TCR (Taylor et al., 2012) is determined from fully coupled ESM experiments in which atmospheric CO₂ concentrations are increased by 1% annually beginning from an equilibrated preindustrial (piCTL) simulation. The TCR is defined as the average warming in Years 61–80 of the 1%CO2 experiment, that is, when CO₂ concentrations are about 2 times the piCTL value (see Appendix A for details of our calculation). Figure 16 shows $\mathrm{\Delta }\bar{T}{t}_{1\%}$ (Equation A1) as a function of time for CESM1 and CESM2 1%CO2 experiments. The two curves are close to each other through Year 100, and TCR values determined from these curves are also very close, 2.1 K for CESM1 and 2.0 K for CESM2 (Table 1; Meehl, Arblaster, et al. 2020). Based on the standard errors for the TCR estimates in Table 1 we conclude that TCR in CESM1 and CESM2 is not significantly different. Nevertheless, we will see that many other aspects of the 1%CO2 experiments for CESM1 and CESM2 exhibit what appear to be large and significant differences.

16 Figure. Warming ΔT‾1% (Appendix A) as a function of time for 1%CO2 experiments using CESM1 (gray) and CESM2 (blue). Arrows positioned at Years 60–80 indicate transient climate sensitivity (TCR) values for CESM1 (2.1 K) and CESM2 (2.0 K). TCR is defined as the mean of ΔT‾1% over Years 60–80 in the 1%CO2 scenario. Dashed lines show results of step model calculations. Cyan line labeled CESM2** is constructed by replacing North Atlantic, North Pacific, and Arctic T time series in CESM2‐1%CO2 with those in CESM1‐1%CO2.

Figure 16 shows that after Year 100, the $\mathrm{\Delta }\bar{T}{t}_{1\%}$ values in CESM1‐1%CO2 and CESM2‐1% begin to diverge. Consistent with its higher sensitivity, CESM2 warms more rapidlythan CESM1 from Years 100 to 150. The linear trends over Years 100–150 are 0.41(0.02) K decade⁻¹ for CESM1 and 0.52(0.01) K decade⁻¹ for CESM2. Figure 17 shows regional time series of ocean surface temperature T_k. These exhibit dramatic differences between CESM1‐1%CO2 (gray) and CESM2‐1%CO2. (blue). Tropical ocean warming (Figure 17c) is more pronounced in CESM2 than in CESM1 throughout the 1% experiments, and temperatures in the midlatitude Southern Ocean (Figure 17d), while initially lower in CESM2‐1%CO2 than in CESM1‐1%CO2, also increase more rapidly in CESM2 throughout the 1%CO2 experiments. The behavior of T in these regions is consistent with that seen in the 4xCO2 experiments (red and black curves), that is, in both regions CESM2 warms more rapidly in both 4xCO2 and 1%CO2 scenarios.

17 Figure. Regional annual‐mean surface temperature Tk as a function of time for analysis regions in Figure 5; CESM1‐1%CO2 (gray), CESM2‐1%CO2 (blue), CESM1‐4xCO2 (black), and CESM2‐4xCO2 (red).

There is an interesting reversal of this consistency in northern ocean (Figures 17a and 17b) and land (Figure  17f) regions. In these regions, CESM1‐1%CO2 warms more strongly than CESM2‐1%CO2, albeit starting from cooler initial conditions. In the North Atlantic/Pacific region (Figure 17b) CESM1‐1%CO2 is briefly almost 1 K warmer than CESM2‐1%CO2 around Year 110. This is a marked contrast with the behavior of the 4xCO2 experiments, in which northern oceans are much warmer, and also warm more rapidly, in CESM2 than in CESM1. These regional differences clearly have implications for the interpretation of abrupt CO₂ increase experiments with respect to 1%CO2 experiments, and will be explored in future studies.

The step model approach (Gregory et al., 2015) establishes a connection between the warming ΔT_4x(t) in abrupt 4xCO2 experiments discussed earlier and the time behavior of ΔT_1%(t) in the 1%CO2 scenario. The step model as applied here is simply the sequential sum of increments fΔT_4x(t − t_i) over the 140 yr t_i since the beginning of the 1%CO2 run, where $f=\frac{1}{140}$ and the increments are nonzero only when (t − t_i) > 0. Results from the step model applied to CESM1 and CESM2 4xCO2 warming are shown in Figure 16(dashed lines). There appears to be rough agreement between the 1%CO2 and step model time series of ΔT, although late in the simulation the step model warms less rapidly than the 1%CO2 run. The CESM2 step model warms more rapidly than CESM2‐1%CO2 for all years. TCR values derived from step model time series are 2.1 K for CESM1 and 2.5 K for CESM2, indicating a more nonlinear response to increased CO₂ in CESM2.

The cyan curve in Figure 16 labeled CESM2^∗∗ shows a synthetic curve ΔT_1%(t) that was constructed by replacing the CESM2‐1%CO2 regional time series of ΔT_s; k for the Arctic and N Atlantic/N Pacific regions (Figures 16a and 16b) with those from CESM1‐1%CO2. This provides an estimate of the global impact of the reduced warming seen for these regions in CESM2‐1%CO2. The CESM2^∗∗ time series does in fact warm more rapidly than that of CESM2‐1%CO2 beginning around 60. TCR for CESM2^∗∗ is around 2.2 K. This is evidence that Arctic and other high‐latitude processes in CESM2 are responsible for the low TCR obtained for this model. Whether these processes in CESM are principally atmospheric or are instead driven by ocean heat transport is not clear at the moment.

The contrasting behavior of ECS, iECS and TCR in CESM1 and CESM2 is interesting. Clearly, these two versions of CESM do not suggest a linear relationship between TCR and ECS as identified by Flato et al. (2014). The similarity between TCR in CESM1 and CESM2 may be largely spurious, masking large and significant differences in regional warming. The existence of strong North Atlantic cooling in CESM2‐1%CO2 compared with CESM1‐1%CO2 contrasts sharply with the behavior of 4xCO2 runs and suggests an important difference in ocean heat transport in 1%CO2 versus 4xCO2 scenarios. The relationship between TCR and ECS in CMIP ensembles has been investigated by Meehl, Senior, et al. (2020).

This North Atlantic cooling may be responsible for a delayed response of Greenland temperatures and surface mass balance (SMB) in CESM2‐1%CO2 experiments. Figure 18 shows T trends for the North Atlantic and Greenland during 150 yr of increasing CO₂. The North Atlantic warms for 40 yr, after which temperatures are flat or slightly decreasing until around Year 90, and then turn sharply upward. Similarly, Greenland temperatures are flat during Years 40–90 before increasing steeply. Sellevold and Vizcaino (2020) have analyzed Greenland Ice Sheet SMB changes, which are driven mainly by increased surface melting associated with warmer temperatures, in a 150‐year CESM2‐1%CO2 experiment. They found that SMB decreases modestly, by 2.5 ± 0.4 Gt yr⁻², during Years 1–90, and much more quickly, by 15.9 ± 1.1 Gt yr⁻², after Year 90. Thus, the Greenland SMB and resulting sea level contribution in 1%CO2 experiments appear to be closely linked to North Atlantic temperatures and ocean heat transport.

18 Figure. Regional‐mean surface temperature Tk as a function of time in the CESM2‐1%CO2 experiment: (a) North Atlantic; (b) Greenland. The respective regions are shown in the panel insets.

Gregory et al. (2015) examined the role of ocean heat uptake in the 1%CO2 scenario. They found increasing rates of warming in CMIP5 1%CO2 experiments with time which they attribute to declining ocean heat uptake. This leads to the nonlinearity in ΔT_1%(t) seen in Figure 16. A convenient measure of this nonlinearity is the ratio $\frac{{⟨\mathrm{\Delta }{T}_{1\%}⟩}_{140}}{\text{TCR}}$ where ⟨ΔT_1%⟩₁₄₀ is the warming around the time of CO₂ quadrupling in the 1%CO2 scenario (year  $=140$). For the CMIP5 ensemble mean this ratio is around 2.4. We derive ⟨ΔT_1%⟩₁₄₀ of 4.9 K for CESM1‐1%CO2 and 5.1 K for ESM2‐1%CO2 (Appendix A) giving ratios of 2.3 and 2.6, respectively. These measures of warming may be better indications of expected conditions at the end of the 21st Century for various realistic scenarios (Gregory et al., 2015), and they appear to better capture differences between CESM1 and CESM2. Nevertheless these measures still mask the profound regional differences in warming evident in Figures 17 and 18.

Summary and Discussion

This study examined abrupt CO₂ and 1%CO₂ increase simulations using two versions of the CESM, CESM1 and CESM2. We used results from extended (800+ yr) 4xCO2 experiments using fully coupled, ESM configurations with a dynamic ocean to investigate the origins of CESM2's substantially increased ECS compared to CESM1. Table 1 showed several estimates of ECS for CESM1 and CESM2. Values of inferred ECS (iECS) from linear regression of net TOM radiative imbalance as a function of global‐mean temperature, $\bar{\mathcal{N}}$ versus $\mathrm{\Delta }\bar{T}$, for 4xCO2 experiments (Gregory et al., 2004) depend strongly on the number of years in the regression. In all cases, however, CESM2's iECS is 1 to 2 K higher than that of CESM1 (Figure 1b), with values of up to 6.5 K for iECS derived from 800 yr of CESM2‐4xCO2.

Contributions to the increased sensitivity of CESM2 from initial forcing and from radiation feedbacks were examined in section 4.1. We found an increase in initial forcing ${\bar{\mathcal{N}}}_0$ in CESM2 of around 1.2 W m⁻² compared to CESM1‐4xCO2 (Table 3), which appears to originate in rapid initial adjustments of shortwave fluxes (Table  3). A simple calculation showed that the increased initial forcing contributes as much as half of the increased sensitivity diagnosed from CESM2‐4xCO2. However, in CESM2 SOM experiments using 2xCO2 and 4xCO2 forcing (section 5) we found that ${\bar{\mathcal{N}}}_0$ responds nonlinearly to CO₂ increase while radiation feedbacks in CESM2‐2xCO2‐SOM and CESM2‐4xCO2‐SOM remain constant.

Longwave and shortwave contributions to the net radiation feedbacks in CESM1 and CESM2 were separated. We found that global longwave feedbacks in CESM1 and CESM2 are similar, while shortwave feedbacks in the two models are substantially different (Figure 2). Positive shortwave feedback in Years 100–800 of the 4xCO2 simulations is significantly higher in CESM2 (1.50 W m⁻² K⁻¹, Table 4) than in CESM1 (1.23 and 1.32 W m⁻² K⁻¹). The increased shortwave feedback in CESM2 is responsible for reducing the strength of the net radiation feedback ${\bar{\lambda }}_{\mathcal{N}}$ (Equation 10), which in turn increases climate senstivity. In addition, shortwave feedbacks are responsible for the highly nonlinear behavior of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ observed in CESM2‐4xCO2.

In sections 4.3 and 4.4, we analyzed regional contributions to the global radiation feedbacks using the decomposition in Equation 7. The largest contributions to changing long‐term (Years 100–800) radiation feedbacks in th models arises from midlatitude Southern Ocean low clouds, and from complex tropical ocean feedbacks involving, low, midlevel and high clouds, with about half of the change in net global feedback arising in each of these region (Figure 8f).

The APRP technique (Taylor et al., 2007) was employed (section 4.4.2)to analyze the contribution of different cloud processes to shortwave feedbacks. APRP showed that the increased feedbacks in CESM2 are related to increased cloud scattering feedback (Figure 12). We examined the evolution of cloud condensate phase in high latitudes and midlatitudes (Figure 11). CESM2 is characterized by a much larger proportion of liquid‐phase clouds. Over the midlatitude Southern Ocean we found dramatically enhanced feedback for liquid condensate in CESM2 (−4.5 g m⁻² K⁻¹) compared to CESM1 (−1.7 g m⁻² K⁻¹) but also somewhat stronger feedback for ice condensate in CESM1 than in CESM2 (−0.7 g m⁻² K⁻¹ vs. 0.2 g m⁻² K⁻¹). Thus, increased scattering feedback over the Southern Ocean in CESM2 could result from stronger condensate amount feedback, or from reduced negative cloud phase feedback (e.g., Frey & Kay, 2018). Without further analysis we cannot quantify the role of these two feedback processes. Our results are also consistent with analyses by Gettelman, Hannay, et al. (2019) who found increased Southern Ocean radiation feedbacks in CAM6 versus CAM5 in SST+4 K experiments, which they attribute to changes in the treatment of ice nucleation in the two models.

In section 5 we compared results from SOM runs with those from the the fully coupled ESM configurations of CESM1 and CESM2. ECS estimated from SOM runs (ECS‐SOM) has been proposed as a way to reduce the resources required to calculate ECS (e.g., Danabasoglu & Gent, 2009). ECS‐SOM using 2xCO2 forcing has increased from about 4 K in CESM1 to 5.4 K in CESM2 (Table 1). We found that ECS‐SOM(4x) derived from SOM runs subject to 4xCO2 increase agrees remarkably well with iECS derived from long ESM simulations. In addition there is also remarkable similarity in the evolution of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ between ECS and SOM 4xCO2 experiments (Figure 14). These similarities occur despite the presence of significant regional differences in warming (Figures 15a and 15b).

In contrast to ECS the TCR has not changed between CESM1 and CESM2 (Table 1). TCR is defined as the warming present around Year 70 in experiments subject to a 1% annual increase in CO₂, that is, around the time of CO₂ doubling. In section 6 we examined the evolution of 1%CO2 CESM1 and CESM2 experiments. While TCR has not changed between CESM1 and CESM2, there are large regional differences in warming between CESM1‐1%CO2 and CESM2‐1%CO2. Tropical and midlatitude Southern Oceans warm more rapidly in CESM2‐2%CO2 than in CESM1‐1%CO2, consistent with the higher ECS of CESM2 (Figure  17). However, the Arctic and North Atlantic/North Pacific in CESM1‐1%CO2 and CESM2‐1%CO2 behave very differently from what would be expected from their behavior in the 4xCO2 configuration. Northern oceans in CESM2‐1%CO2 warm more slowly than in CESM1‐1%CO2. The North Atlantic in CESM2‐1%CO2 shows a dramatic multidecadal cooling from Years 40 to 80 (Figure 18a). We applied the step model analysis (Gregory et al., 2015) to relate ECS and TCR changes in CESM1 and CESM2 and found that the TCR in CESM2 is lower than that predicted from the step model using 4xCO2 results (Figure 16).

This study explored the evolution of a single modeling system in response to increased CO₂ forcing. We hope this analysis will help in the design of multimodel studies that compare ECS and TCR across the CMIP5 and CMIP6 ensembles. Our study again points to the importance of shortwave CREs in determining model climate sensitivity and suggests a key role for ice phase and mixed‐phase microphysics both in high‐latitude low clouds. We also identified a complex pattern of tropical cloud changes in the CMIP5 and CMIP6 versions of our model that also contributed to increasing ECS in our CMIP6 model. This may encourage other modeling groups to examine their tropical simulations and potential impacts on ECS. Our study also suggests that model TCR may miss significant regional responses to increasing CO₂, especially in high latitudes. Both 4xCO2 and 1%CO2 experiments may yield insight into coupled model behavior in more realistic forcing scenarios.

A AppendixCalculation of ECS and TCR

Calculations of ECS and TCR are subject to uncertainties due both to internal variability in model simulations and to details in calculation procedures, such as the specification of preindustrial reference temperatures, detrending techniques, etc. Here we describe how the numbers in Table 1 were derived and examine sensitivities to details in the calculations.

Infrerred ECS (iECS) and TCR are derived from 4xCO2 and 1%CO2 simulations and their respective preindustrial control (piCTL) simulations. We denote the year in which 4xCO2 and 1%CO2 simulations branch from their piCTL by Y_b. The duration of the experiments beyond Y_b is denoted by ΔY_exp. According to the CMIP protocols ΔY_exp is 150 yr for both the 1%CO2 experiment the 4xCO2 experiment. The piCTLs for CESM1 and CESM2 also run through the period Y_b to Y_b + 150. Linear fits to the global‐mean surface temperature $\bar{T}$ from the piCTLs during this period are performed, which we denote by $T_l^{\ast }(t)$.

To calculate TCR we first subtract $T_l^{\ast }(t)$ from $\bar{T}{t}_{1\%}(t)$, the time series of global‐mean surface temperature for the corresponding 1%CO2 experiment: [Image Omitted. See PDF]

TCR is then the average of $\mathrm{\Delta }\bar{T}{t}_{1\%}(t)$ over Years 61–80 of the 1%CO2 experiment. This procedure follows that in the ESMValTool (Righi et al., 2020) except that surface temperature is used instead of 2‐m temperature. This approach gives TCR values of 2.1 K (0.07 K) for CESM1 and 2.0 K (0.04 K) for CESM2 where the standard errors are shown in parentheses. Standard errors are calculated using bootstrapping with replacement. Bootstrapping is applied to the linear fit $T_l^{\ast }$ as well as to the 20‐year mean of $\mathrm{\Delta }\bar{T}{t}_{1\%}(t)$.

A second average of the warming over Years 131–150, ${⟨\mathrm{\Delta }\bar{T}{t}_{1\%}⟩}_{140}$, is also calculated to characterize the warming attained in the 1%CO2 scenario when CO₂ values have approximately quadrupled, that is, around Year 140 (Gregory et al., 2015). The procedure is identical to that used for the TCR calculate except for the averaging period used. We obtain ${⟨\mathrm{\Delta }\bar{T}{t}_{1\%}⟩}_{140}$ values of 4.9 K (0.08 K) for CESM1‐1%CO2 and 5. K (0.08 K) for CESM2‐1%CO2.

To calculate iECS, a linear fit to 150 years of $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ from the 4xCO2 experiment is performed. Here $\mathrm{\Delta }\bar{T}$ is defined as the difference of $\bar{T}$ from the 4xCO2 experiment with respect to the average of $\bar{T}$ from the piCTL over Years Y_b to Y_b + 150. The linear fit to $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ may be expressed as [Image Omitted. See PDF] where ${\bar{\lambda }}_{\mathcal{N}}$ and ${\bar{\mathcal{N}}}_I$ are the slope and intercept of the linear fit. Note that elsewhere in the text we use ${\bar{\mathcal{N}}}_0$ to refer to the intercept for a linear fit to $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ over Years 1–20. This particular interval is used to estimate the initial radiative forcing in the 4xCO2 simulations. In the absence of nonlinearity in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ there would be no significant difference between these quantities.

Strictly, the Gregory method is applied to the change in TOA net radiation in the 4xCO2 experiment from the preindustrial values as a function of the change in global‐mean surface air temperature. However, previous studies of climate sensitivity with CESM have used surface skin temperature (e.g., Danabasoglu & Gent, 2009; Gettelman et al., 2012). TOA radiation in a low‐top (∼2 hPa) atmosphere model such as CAM is a highly parameterized quantity that has changed significantly between CESM1 and CESM2. Therefore, we have performed the bulk of our diagnostics using TOM radiation and surface temperature. However, in Table A1) we show that these details in how iECS is calculated result in only small differences of order ∼0.1 K (Table A1).

A1 TableFour Estimates of iECS—Equilibrium Climate Sensitivity Derived Using the Gregory Method—Using Different Radiation and Temperature Quantities

CESM1	CESM2
iECS (K)
$\mathcal{N}$ versus ΔTs
3.36 (0.04)	5.25 (0.22)
$\mathrm{\Delta }\mathcal{N}$ versus ΔTs
3.39 (0.04)	5.21 (0.22)
$\mathrm{\Delta }{\mathcal{N}}_{TOA}$ versus ΔTs
3.36 (0.04)	5.10 (0.20)
$\mathrm{\Delta }{\mathcal{N}}_{TOA}$ versus ΔT2m
3.36 (0.04)	5.15 (0.20)
piCTL radiation (W m−2)
$\mathcal{N}$
−0.07	0.05
${\mathcal{N}}_{TOA}$
2.04	0.66

Note. First row shows results using raw TOM radiative balance $\mathcal{N}$ and surface skin temperature warming ΔT_s. Second row shows the impact of removing the small preindustrial biases from $\mathcal{N}$. Third row shows iECS calculated using the change in TOA radiative balance from the control run $\mathrm{\Delta }{\mathcal{N}}_{TOA}$ and ΔT_s. Fourth row shows iECS calculated from $\mathrm{\Delta }{\mathcal{N}}_{TOA}$ and the warming of the 2‐m air temperature ΔT_2m. Fifth and sixth rows show mean radiative imbalances at TOM and TOA, respectively, from Years Y_b to Y_b+150 in the piCTL experiments.

Equation A2 is inverted for ${\bar{\mathcal{N}}}_l=0$ to give an equilibrium $\mathrm{\Delta }\bar{T}$, which is divided by 2 in 4xCO2 experiments to give the expression for iECS in Equation 10. This approach gives iECS values of 3.4 K (0.04 K) for CESM1 and 5.3 K (0.22 K) for CESM2

The calculation of iECS(800) based on 800 yr of 4xCO2 differs from the conventional iECS only in how the piCTL $\bar{T}$ reference is defined. Since the piCTL simulations did not extend for 800 yr past Y_b we use an average of the linear fit ${\bar{T}}_{;l}(t)$ extrapolated through Year Y_b+800 to define $\mathrm{\Delta }\bar{T}$. Using this method, we derive values of iECS(800) of 4.2 K (0.05 K) for CESM1 and 6.5 K (0.07 K) for CESM2.

Again, our approach for estimating iECS from 4xCO2 experimental results is close to that outlined by Righi et al. (2020), with the difference that we use warming estmated from T_s instead of T_2m and net TOM radiative balance $\mathcal{N}$ rather than the change in TOA radiative imbalance $\mathrm{\Delta }{\mathcal{N}}_{TOA}$. Table A1 shows iECS estimates derived in four different ways.

The procedure for deriving ECS‐SOM estimates from SOM configurations is less well established. We would like to use multiyear averages of $\bar{T}$ from well equilibrated control and 2xCO2 or 4xCO2‐SOM experiments to define ECS‐SOM. In practice, the choice of averaging periods is somewhat subjective and can lead to small differences in estimates of ECS‐SOM. For example, in Figure A1a we show time series from three SOM experiments using CESM2.0 (1xCO2 in black, 2xCO2 in green, and 4xCO2 in red). Note that all of these experiments are initialized from the same unequilibrated atmospheric state.

A1 Figure. (a) Time series of T‾ from SOM integrations. (b) Inferred iECS derived from SOM runs.

Gettelman, Hannay, et al. (2019) used averages over Years 40–60 for both the control and 2xCO2 simulations to derive an ECS‐SOM of 5.3 K for CESM2. If a later period is used for the CESM‐2xCO2‐SOM (green curve) this estimate will increase since a small additional warming occurs after Year 60. The ECS‐SOM of 5.5 K for CESM2 in Table 1 is calculated using an average of Years 70–100 for the 2xCO2 experiment and a reference temperature averaged over Years 20–75 of CESM2‐1xCO2‐SOM. We do not advocate either value but simply present both to illustrate the level of uncertainty that may exist in published numbers for ECS‐SOM. ECS‐SOM(4x) is calculated using the same reference temperature and an average temperature over Years 70–100 of CESM2‐4xCO2‐SOM. The difference between these values is divided by 2 to account for the 4xCO2 versus 2xCO2 increase.

Another approach to estimating ECS‐SOM is to apply the Gregory et al. (2004) approach to $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ from the SOM runs. Results of this approach are shown in Figure A1b. Interestingly, the results of this method for CESM2‐2xCO2‐SOM (green) appear to converge on an ECS‐SOM value of around 5.2 K, closer to the Gettelman, Hannay, et al. (2019) value, even though this number is based on what appears to be slightly unequilibrated $\bar{T}$ from the 2xCO2‐SOM experiment. We note, however, that the Gregory et al. (2004) method suffers from the same pitfalls when applied to SOM $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ results as it does when applied to full ESM results; that is, rapid initial adjustment can affect the regression estimate of ${\bar{\lambda }}_{\mathcal{N}}$. As with full ESM results, better estimates of ECS may be obtained if initial rapid adjustment in $\bar{\mathcal{N}}(\mathrm{\Delta }\bar{T})$ is excluded (Dunne et al., 2020).

The calculation details discussed in this appendix have only small impacts on estimates of TCR and ECS, generally less than a few tenths of a degree K. We present them to explain possible discrepancies in published numbers of TCR and ECS for CESM.

B AppendixRegional Feedback Analysis

Regional or local analyses of radiation feedbacks have been discussed by Armour et al. (2013) and Andrews et al. (2015). Andrews et al. (2015) in particular questions whether nonlocal effects in the atmosphere interfere with a local analysis. To probe the validity of a regional approach we have performed regional “checkerboard” decompositions into rectangular regions of 60° × 45°, 30° × 22.5°, 15° × 11.25°, and 7.5° × 5.125°. Figure B1 shows gridded maps of ${\delta }_T{\mathcal{S}}_k$ the regression of regional ${\mathcal{S}}_k$ with global‐mean temperature T, as well as $A_k{\lambda }_{\mathcal{S};k}$ the local feedback weighted by the local warming amplification factor (section 3.1) and the difference of these quantities for CESM1‐4xCO2 Years 100–800 using a decomposition into 15° × 11.25° regions. The area a_k weighted sum of ${\delta }_T{\mathcal{S}}_k$ will give the global‐mean shortwave feedback. If the regional decomposition on this scale is valid, $A_k{\lambda }_{\mathcal{S};k}$ should be close to ${\delta }_T{\mathcal{S}}_k$.

B1 Figure. Regional shortwave feedbacks in CESM1‐4xCO2 in a 15° × 11.25° resolution “checkerboard” domain decomposition: (a, b, c) Results without decadal averaging; (d, e, f) results with decadal averaging. Panels (a) and (d) show regression of shortwave radiation with respect to global T‾, that is, δTSk. Panels (b) and (e) show AkλS;k. Panels (c) and (f) show the difference of the leftmost and middle panels.

Figure B1 show results with and without decadal averaging applied to the time series of ${\mathcal{S}}_k$ and A_k. Without decadal averaging (top row) we see significant differences between $A_k{\lambda }_{\mathcal{S};k}$ and ${\delta }_T{\mathcal{S}}_k$ in the tropics. Decadal averaging (bottom row) largely removes these differences suggesting that with decadal averaging a regional analysis is valid even at scales of 15° × 11.25°.

Figures B2 and B3 summarize the results of our tests. In Figure B2 the first four pairs of bars in each panel show the direct sums of $a_k{A}_k{\lambda }_{\mathcal{S};k}$ (Figures B2a and B2c) and $a_k{A}_k{\lambda }_{\mathcal{L};k}$ (Figures B2b and B2d) over all of the rectangular regions in the four decompositions applied. The fifth pair of bars shows independent global regressions. For the regional analysis to be useful all bars should match those in the fifth position.

B2 Figure. Comparison of direct sums of akAkλS;k and akAkλN;k over all domains in checkerboard decompositions with independent global regressions of S and L with T. (a) Results for shortwave forcing in CESM1‐4xCO2 Years 100–800. Black bars show results with decadal averaging applied to time series. Gray bars show results without decadal averaging. (b) As in (a) except for longwave. (c, d) As in (a) and (b) except for CESM2‐4CxO2. Red bars show results with decadal averaging; lighter red shows results without decadal averaging.

B3 Figure. Normalized global RMS values of ϵX;k=AkλX;k−δTX, where X=S or L, for “checkerboard” decompositions. Values are normalized by RMS values of S and L in each decomposition: (a, b) CESM1‐4xCO2; (c, d) CESM2‐4xCO2. Results with and without decadal averaging are shown. Bars are shaded as in Figure B2.

Results with and without decadal averaging are shown. It is striking that decadal averaging has little impact on CESM2 results (Figures B2c and B2d) while remaining critical for the regional approach to work at any scale in CESM1. In the longwave analysis especially, we see that results in CESM1 also improve as larger regions are used. By contrast, the regional analysis for CESM2 appears to work reasonably well for all scales.

Good agreement between the summed regional decompositions in Figure B2 and the independent global regressions could still occur if there is cancelation of regional errors. A stricter metric is shown in Figure  B3. The figure shows the global RMS values of the quantities ${ϵ}_{\mathcal{S};k}=A_k{\lambda }_{\mathcal{S};k}-{\delta }_T{\mathcal{S}}_k$ (Figures B3a and B3c) and ${ϵ}_{\mathcal{L};k}=A_k{\lambda }_{\mathcal{L};k}-{\delta }_T{\mathcal{L}}_k$ (Figures B3b and B3d). The quantities have been normalized by the RMS values of ${\delta }_T{\mathcal{S}}_k$ and ${\delta }_T{\mathcal{L}}_k$, respectively. A value of 1 means that the global‐mean variance in ϵ is as large as the variance in ${\delta }_T{\mathcal{S}}_k$ and ${\delta }_T{\mathcal{L}}_k$. We do not have a well established criterion for success using this metric. However, we note that the bars in Figure B3a for 15° × 11.25° correspond to the maps shown in Figure B1. This suggests values of 0.2 or less in Figure B3 indicate satisfactory performance by the regional decompostion.

This analysis suggests that a regional analysis is useful for sufficiently large regions and for atmospheric time series that have had shorter interannual time scales removed. The precise meaning of “sufficiently large” is not yet clear, however, we speculate it is related to the size of circulation anomalies induced by tropical deep convection in warming climates. The wider range of scales for which the approach is valid in CESM2 suggests that tropical deep convection in CESM2‐4xCO2 is producing predominantly local changes.

Acknowledgments

The authors are grateful for interesting and constructive comments on this work from Tim Andrews and an anonymous reviewer. These comments led to substantial improvements in the paper. We would also like to acknowledge helpful discussions with Julie Arblaster, Gokhan Danabasoglu, Dave Lawrence, and Jerry Meehl. The CESM project is supported primarily by the National Science Foundation (NSF). This material is based upon work supported by the National Center for Atmospheric Research, which is a major facility sponsored by the NSF under Cooperative Agreement No. 1852977. Computing and data storage resources, including the Cheyenne supercomputer (https://doi.org/10.5065/D6RX99HX), were provided by the Computational and Information Systems Laboratory (CISL) at NCAR.

Data Availability Statement

Model output from CESM simulations used in this study (Table 2) along with output from preindustrial control simulations is available online (at https://doi.org/10.5065/8406‐p773).