Exploring climate stabilisation at different global warming levels in ACCESS-ESM-1.5

Under the Paris Agreement, signatory nations aim to keep global warming well below 2 °C above pre-industrial levels and preferably below 1.5 °C. This implicitly requires achieving net-zero or net-negative greenhouse gas emissions to ensure long-term global temperature stabilisation or reduction. Despite this requirement, there have been few analyses of stabilised climates, and there is a lack of model experiments to address our need for understanding the implications of the Paris Agreement. Here, we describe a new set of experiments using the Australian Community Climate and Earth System Simulator Earth system model (ACCESS-ESM-1.5) that enables the analysis of climate evolution under net-zero emissions, and we present initial results. Seven 1000-year-long simulations were run with global temperatures stabilising at levels in line with the Paris Agreement and at a range of higher global warming levels (GWLs). We provide an overview of the experimental design and use these simulations to demonstrate the consequences of delayed attainment of global net-zero carbon dioxide emissions. We show that there are substantial differences between transient and stabilising climate states and differences in stabilisation between GWLs. As the climate stabilises under net-zero emissions, we identify significant and robust changes in temperature and precipitation patterns including continued Southern Ocean warming and changes in regional precipitation trends. Changes under net-zero emissions differ greatly between regions, including contrasting trajectories of sea ice extent between the Arctic and Antarctic. We also examine the El Niño–Southern Oscillation (ENSO) and find evidence of reduced amplitude and frequency of ENSO events under climate stabilisation relative to projections under transient warming. An analysis at specific GWLs shows that significant regional changes continue for centuries after emission cessation and that these changes are stronger at higher GWLs. Our findings suggest substantial long-term climate changes are possible even under net-zero emission pathways. These simulations are available for use in the community and will hopefully motivate further experiments and analyses based on other Earth system models.

1 Introduction

The world has warmed by over 1 °C above pre-industrial levels due to anthropogenic greenhouse gas emissions (Haustein et al., 2017; IPCC, 2021a). Anthropogenic emissions of carbon dioxide have had major impacts on the climate already, and these impacts are worsening as the rate of emissions continues at record-high levels (Friedlingstein et al., 2023; Liu et al., 2024). The growing impacts of climate change motivated the Paris Agreement of 2015 in which signatory nations agreed to aim to limit global warming to well below 2 °C above pre-industrial levels and preferably to below 1.5 °C. Limiting global warming in line with the Paris Agreement implicitly requires at least stabilising global temperatures if not cooling back towards pre-industrial conditions (Rogelj et al., 2017). Thus, it is imperative that humanity gains a greater understanding of the long-term consequences of low global warming levels and, in contrast, the implications of failing to meet the Paris Agreement through continued global warming or stabilisation at higher global warming levels.

Model analyses suggest that there is a near-linear relationship between cumulative carbon dioxide emissions and global average temperature change in a transient climate (Allen et al., 2022; IPCC, 2021a; Seneviratne et al., 2016). Furthermore, Earth system model (ESM) simulations and simpler model runs performed as part of the Zero Emissions Commitment Model Intercomparison Project (ZECMIP; Jones et al., 2019) suggest that the cessation of carbon dioxide emissions would result in an almost immediate halt to global warming and near-zero global temperature change for the following century, albeit with uncertainty between models (MacDougall et al., 2020; Palazzo Corner et al., 2023). These studies suggest that an emission level very near zero is required to halt global warming in line with the Paris Agreement.

Since the Paris Agreement was signed, there has been a significant effort to understand what 1.5 and 2 °C global warming levels (GWLs) would entail in terms of climate change impacts, but this body of work has not focused on net-zero emission scenarios. This is primarily due to a lack of availability of suitable model simulations in the immediate aftermath of the Agreement being developed. Studies on the 1.5 and 2 °C GWLs were in part motivated by the urgent need for literature ahead of the Intergovernmental Panel on Climate Change (IPCC) Special Report on 1.5 °C global warming (SR1.5; Masson-Delmotte et al., 2018). Work was undertaken to examine 1.5 and 2 °C GWLs, and rapidly inform SR1.5, using existing transient multi-model ensembles (e.g. King et al., 2017; Schleussner et al., 2016), pattern scaling (King et al., 2018; Seneviratne et al., 2016; Tebaldi and Knutti, 2018), single coupled model ensembles (Sanderson et al., 2017), and multi-model atmosphere-only ensembles (Mitchell et al., 2017). However, these methods are all, to some extent, based on transient climate states. The characteristics and relative merits of widely used methods for examining GWLs are detailed by James et al. (2017), with a brief summary provided here.

Time sampling methods, based on extracting GWLs from existing simulations, generally use a combination of the Scenario Model Intercomparison Project (ScenarioMIP; O'Neill et al., 2016) projections with continued increases in greenhouse gas concentrations and global warming. Time sampling has become commonly used for examining GWLs (King et al., 2017; Nangombe et al., 2018; Schleussner et al., 2016) and was the principal method employed in the IPCC Sixth Assessment Report Working Group 1 (IPCC, 2021a). This method samples from continually warming transient simulations.

Pattern scaling also typically makes use of ScenarioMIP experiments (King et al., 2018; Tebaldi and Arblaster, 2014; Tebaldi and Knutti, 2018), although earlier studies used carbon-dioxide-only forced simulations (Mitchell, 2003). The method involves the identification of well-understood climate changes expected under high greenhouse gas forcings and assumes a linear scaling with temperature changes to estimate patterns under lower greenhouse gas forcing or GWLs. This method typically samples from continually warming simulations, but it has previously been applied to longer slow-warming simulations (Mitchell, 2003).

Bespoke simulations using the Community Earth System Model (Sanderson et al., 2017) achieve slow global temperature change through near-zero carbon dioxide emission pathways in the last few decades of the 21st century. While these represent a climate nearer to stabilisation than the other methods described here, they are only run to 2100 and thus have not had enough time for slow processes to continue and affect the global climate.

Global atmosphere-only model simulations were developed as part of the “Half a degree additional warming, prognosis and projected impacts” (HAPPI) project (Mitchell et al., 2016, 2017). Simulations were run under warmed sea surface temperatures and increased greenhouse gas concentrations derived from coupled model differences between the recent climate and low-end ScenarioMIP simulations in the late 21st century. This experimental setup is derived from transient climate data but based on slower changes than the time sampling method typically uses.

Figure 1 (a, b) Schematics of the proposed design of a multi-model ensemble of simulations at targeted global warming levels, adapted from King et al. (2021a) and implemented with ACCESS-ESM-1.5 here. (c) Atmospheric carbon dioxide (CO₂) concentrations and (d) global mean surface temperature (GMST) anomalies, relative to the 1850–1900 baseline, in the ACCESS-ESM1-5 historical and SSP5-8.5 simulations (black) and the net-zero greenhouse gas emission simulations (grading from yellow to red for later emission cessation). The net-zero emission simulations begin in 2030, 2035, 2040, 2045, 2050, 2055, and 2060. The global warming levels for the last 30 years of each net-zero emission simulation are shown in panel (d).

There is a gap between the analyses previously performed to understand the implications of the Paris Agreement (based on transient climate states) and the intent of the Agreement itself (for a stabilised or cooling climate under net-zero or net-negative emissions). It is known that there are significant differences between transient and stabilised climates with respect to regional temperatures (Joshi et al., 2008; King et al., 2020; Manabe et al., 1991), atmospheric circulation and precipitation (Ceppi et al., 2018; Grose and King, 2023; Sniderman et al., 2019), and ocean characteristics (Armour et al., 2016; Long et al., 2014), including sea level rise (Nauels et al., 2019). While the examination of low GWLs in transient climate states is useful for understanding near-term climate changes, it is imperative that the understanding of stabilised and overshoot climate behaviour, in line with the Paris Agreement goals, is improved. Projections based on transient methods may not be used as a proxy for stabilised climates aligned with the Paris Agreement.

There are few long, stabilised climate simulations that would be consistent with the intent of the Paris Agreement. While the tier 1 experiment in the ZECMIP protocol results in near-stable global average temperatures, ESM simulations are initiated from the same cumulative emission levels resulting in a range of GWLs, and the simulations are typically around 100 years in length (Jones et al., 2019). Ensembles of millennial-length simulations based on fixed concentration levels, such as LongRunMIP (Rugenstein et al., 2019) or single-model experiments (Fabiano et al., 2024), warm continuously as fixed carbon dioxide concentrations are equivalent to a low continuous rate of carbon dioxide emissions. Again, these ensembles result in different amounts of warming between models dependent on individual models' equilibrium climate sensitivity.

There are also limited studies that explore stabilisation at different global warming levels. This is important, as the temperature response after net-zero has been found to be dependent on total cumulative emissions (Allen et al., 2022; Jenkins et al., 2022); thus simulations must be run at multiple branching points to fully understand the post-net-zero climate evolution. ZECMIP aims to achieve this; however, there are few branching points, and only a small set of models have been run with more than one branching point.

Figure 2 Figure showing the extracted periods of 1.5, 2, and 3 °C global warming levels in the transient (black) and net-zero emission (yellow to red) simulations. GWLs are extracted as all years within decades of a global average temperature of 1.5, 2, and 3 °C ±0.2 °C in the transient simulations and net-zero emission simulations where the climate is stabilising. Extracted decades are shown in bold. An 1850–1900 baseline is used as a proxy for a pre-industrial climate. Only one transient case is shown here for illustrative purposes, but all 40 concentration-driven SSP5-8.5 ensemble members are used.

Given the knowledge gap described above, new model simulations have been proposed to further address the needs of policymakers responding to the Paris Agreement (King et al., 2021a). In this framework, model simulations are initiated from different points in the high greenhouse gas emission scenario, Shared Socioeconomic Pathway (SSP)5-8.5, and run forward in emission-driven mode with net-zero carbon dioxide emissions and 1850-level non-CO₂ greenhouse gases and aerosol concentrations (Fig. 1a, b). These simulations would be run for several centuries with the starting point chosen in each ESM so that the simulations stabilise near targeted GWLs, including the 1.5 and 2 °C GWLs referred to in the Paris Agreement. The choice of starting point would be ESM-dependent and based on prior knowledge of model responses to forcings as discussed in King et al. (2021a). Simulations in carbon-emission-driven mode more easily correspond to policy targets and offer other advantages relative to concentration-driven projections (Sanderson et al., 2023).

We must note that stabilised simulations are unlikely to represent plausible future scenarios, since it is unlikely that we would stabilise at net-zero without then immediately going to net-negative emissions (including carbon dioxide removal). However, stabilised simulations are a useful baseline and reference for understanding stabilised climates before then exploring different plausible scenarios of net negative emissions.

Here, we describe a set of seven 1000 year simulations that follows the framework of King et al. (2021a) and aims to address the current knowledge gap around the global and regional climate response to stabilisation. These simulations were run using the Australian Community Climate and Earth System Simulator (ACCESS) ESM version 1.5 (Ziehn et al., 2020) described in Sect. 2.1, and the experiments are described in Sect. 2.2. In Sect. 2.3 we describe our analysis methods. In Sect. 3.1, we discuss the results based on temporal evolution of the simulations, and in Sect. 3.2 we examine projections through the framework of GWLs. There are many applications of these simulations, and here we perform a preliminary examination of climate extremes in Sect. 3.3 and the El Niño–Southern Oscillation (ENSO) in Sect. 3.4. In Sect. 4 we provide a summary and conclusions.

Figure 3 Time series of (a) decadal average land–ocean temperature difference anomalies, (b) decadal average global mean precipitation anomalies, (c) annual average global mean sea surface temperature anomalies, (d) annual average global mean ocean temperature anomalies, (e) decadal average Arctic sea ice extent in March, (f) decadal average Arctic sea ice extent in September, (g) decadal average Antarctic sea ice extent in March, and (h) decadal average Antarctic sea ice extent in September. Values are shown in the ACCESS-ESM1-5 historical and SSP5-8.5 simulations (black) and in the net-zero greenhouse gas emission simulations (grading from yellow to red). Anomalies in panels (a)–(d) are computed from an 1850–1900 baseline.

2 Data and methods

2.1 ACCESS-ESM-1.5

ACCESS-ESM-1.5 is Australia's latest-generation ESM and is a participant model in the sixth phase of the Coupled Model Intercomparison Project (CMIP6; Eyring et al., 2016). It comprises the Unified Model version 7.3 (Martin et al., 2010), the Community Atmosphere Biosphere Land Exchange model (CABLE; Wang et al., 2011), the Modular Ocean Model version 5 (Griffies, 2012) with the Whole Ocean Model of Biogeochemistry And Trophic-dynamics (WOMBAT; Oke et al., 2013), and the CICE sea ice model version 4.1 (Hunke and Lipscomb, 2010), and it uses the OASIS-MCT coupler (Craig et al., 2017). The ACCESS-ESM-1.5 model has an atmospheric grid spacing of 1.875° longitude by 1.25° latitude and a variable ocean grid with coarsest spacing of 1° by 1°. ACCESS-ESM-1.5 exhibits small biases in radiative flux terms, atmospheric and oceanic properties, and carbon cycle characteristics (Ziehn et al., 2020). Further details on the configuration and basic model evaluation may be found in Ziehn et al. (2020).

ACCESS-ESM-1.5 has an equilibrium climate sensitivity of 3.87 °C and a transient climate response of 1.95 °C (Ziehn et al., 2020), both of which are well within the range of CMIP6 ESMs (IPCC, 2021b) and within estimated likely ranges from a recent community assessment (Sherwood et al., 2020). Overall, ACCESS-ESM-1.5 evaluates favourably compared with many other CMIP6 models for many important climate system properties, including for recent mean climate conditions (Xu et al., 2021), some ENSO characteristics (Planton et al., 2021), and extremes indices (Kim et al., 2020). The model also exhibits unusually small drift in its pre-industrial control simulation (Irving et al., 2020), so drift is less likely a factor in explaining the projected climate evolution after emission cessation than in other ESMs. The model performs well in relation to the Australian climate and produces a particularly dry projection for the region, so it has been selected for downscaling and detailed analysis for national projections (Grose et al., 2023; Rashid et al., 2022). As with all CMIP6 ESMs, it has shortcomings and deficiencies; for example, the representation of different ENSO flavours appears to be worse in ACCESS-ESM-1.5 than in some other CMIP6 models with respect to the spatial structure of central and eastern Pacific ENSO events (Hou and Tang, 2022).

Figure 4 Time series of number of ice-free months per decade in (a) the Arctic in September and (b) Antarctica in March. Time series are shown for the first concatenated historical and SSP5-8.5 simulation and for each net-zero simulation. Ice-free months are defined as months with average sea ice extent under 1 million km².

2.2 Experiments

For this particular analysis, the ACCESS-ESM-1.5 model was run in emission-driven mode with net-zero carbon dioxide emissions and 1850 levels of other greenhouse gases and anthropogenic aerosols for 1000 years, following King et al. (2021a). These simulations were initialised from seven different timestamps in the esm-SSP5-8.5 r10i1p1f1 simulation: 2030, 2035, 2040, 2045, 2050, 2055, and 2060. These simulations are hereafter referred to as “NZ” and the year of emission cessation (e.g. NZ2030). The net-zero simulations exhibit substantial reductions in atmospheric carbon dioxide concentrations (Fig. 1c) due to uptake by the land and ocean, similar to findings based on ZECMIP ESMs (MacDougall et al., 2020). For comparison, Fig. A1 shows the cumulative CO₂ emissions associated with these net-zero emission simulations and also the historical and SSP5-8.5 simulations. The seven 1000-year-long simulations exhibit very slow changes in global mean temperature such that they are suitable for use in examining the effects of climate stabilisation and differences with transient warming (Fig. 1d). The global mean surface temperature (GMST) slightly decreases in the first 20–50 years after the rapid shift to net-zero emissions and 1850-level non-CO₂ and aerosol concentrations has taken place. The initial cooling is due to the rapid return to 1850-level methane concentrations. Much of the analysis described in Sect. 2.3 is conducted on later periods to avoid the effect of this rapid change. These simulations under net-zero emissions achieve different GWLs ranging from 1.6 to 3.3 °C. We note that the concept of “net-zero” is not often related to global warming that exceeds the Paris Agreement targets when used in public discourse. We retain this descriptor because it is accurate with respect to the experimental design.

The net-zero emission ACCESS-ESM-1.5 simulations were analysed and compared against an ensemble of 40 historical and corresponding concentration-driven SSP5-8.5 simulations run with the same model. There is little difference between the emission-driven and concentration-driven SSP5-8.5 simulations, but a large ensemble exists for the concentration-driven version of the model. This large ensemble of transient simulations and the seven 1000-year-long net-zero emission simulations allow robust comparison of climate states as a function of both time and GWL.

We note that the experimental design here does not represent a plausible scenario given the rapid rate of greenhouse gas emissions and instantaneous cessation, but it allows analysis that is helpful for informing policymaking. Future experiments that follow emission pathways that are more plausible would be complementary to this framework.

2.3 Analysis

2.3.1 Global and local changes

Firstly, global and large-scale annual average characteristics of the stabilised runs were computed. GMST and global annual average precipitation were calculated for the stabilised and transient ensembles. Also, land and ocean characteristics were examined separately. Global annual average land and sea surface temperatures were computed and the difference taken, as this is a commonly used climate change metric (e.g. Braganza et al., 2003). Global annual average temperatures through the entire depth of the ocean were calculated to better understand the inertia of the ocean in response to rapid emissions followed by emission cessation. These metrics are computed as anomalies from the 1850–1900 period. The 1850–1900 period averaged across the 40 ACCESS-ESM-1.5 historical simulations was used as a proxy for pre-industrial conditions. Note that this is consistent with the IPCC Sixth Assessment Report (IPCC, 2021a), but previous work suggests this is around 0.1 °C warmer than estimates of 18th-century temperatures representative of earlier in the Industrial Revolution (Hawkins et al., 2017; Schurer et al., 2017). The land–ocean temperature difference rather than ratio is plotted because the baseline 1850–1900 period is included in the plot and the ratio is unstable.

Additionally, sea ice extent in the Arctic and Antarctic regions was computed for March and September to examine seasonally extreme sea ice extent evolution in the post-net-zero runs relative to the transient case. Several of the metrics analysed exhibit very high interannual variability, including the land–ocean temperature difference, average precipitation rate, and sea ice extent metrics, so moving decadal averages of these indices are plotted.

Other analyses were conducted to examine local temperature responses to net-zero emissions. These include the following:

• Analysis of local temperature changes under net-zero emissions relative to a rapid warming scenario. Local warming relative to global warming was computed in the years 200–399 and 800–999 after emission cessation in the millennial-length simulations. The local warming relative to global average warming was also computed in the 40 SSP5-8.5 simulations for 2030–2069. These patterns were then compared.

• Analysis of local temperature changes during the net-zero emission simulations. The difference in annual average temperatures between the years 800–999 and 200–399 after emission cessation was computed in each of the 1000-year simulations. The effect of delay in emission cessation on local temperature response to net-zero emissions was also investigated.

• Analysis of the timing of peak local temperatures relative to net-zero emissions. The timing of peak 11-year annual average temperatures at each location was identified for concatenated historical, esm-SSP5-8.5, and 1000-year-long net-zero emission simulations.

2.3.2 Global warming level projections

In addition, analysis based on GWLs was conducted in the transient SSP5-8.5 and the stabilising climate in the 1000-year-long net-zero emission simulations (Fig. 2). The method used to extract the GWLs is a time-slicing approach (King et al., 2017; Schleussner et al., 2016). Transient 1.5, 2, and 3 °C GWLs were defined by extracting all years within decades where the average temperature in the SSP5-8.5 simulations was within a ±0.2 °C range of the target GWL (e.g. 1.3–1.7 °C for the 1.5 °C GWL) relative to a pre-industrial climate baseline. The use of decades rather than longer averaging windows reduces overlap between the 1.5 and 2 °C GWLs in the transient simulations and is consistent with King et al. (2017).

Figure 5 Maps of annual local temperature change relative to the global average for (a) SSP5-8.5 in 2030–2069, (b) all net-zero simulations 200–399 years after emission cessation, and (c) all net-zero simulations 800–999 years after emission cessation. Maps of the difference between the patterns of warming in (d) all net-zero simulations 200–399 years after emission cessation relative to SSP5-8.5 in 2030–2069 and (e) all net-zero simulations 800–999 years after emission cessation relative to SSP5-8.5 in 2030–2069.

The years chosen to represent each GWL differ slightly between SSP5-8.5 simulations. By centring time slices on the target GWL, the ensemble average is very close to the target GWL in all cases and adequately represents that climate state. A similar approach was taken for an “early stable” period (100–450 years post-emission cessation) using NZ2035 to generate a 1.5 °C GWL ensemble, NZ2045 to generate a 2 °C ensemble, and NZ2060 to generate a 3 °C ensemble. Corresponding “late stable” (650–1000 years post-emission cessation) GWLs were extracted from the NZ2030 run to generate a 1.5 °C ensemble, from NZ2040 to generate a 2 °C ensemble, and from NZ2055 to generate a 3 °C ensemble. Note that each of the “early stable” and “late stable” GWL ensembles could have been compiled using data from multiple simulations, but we chose not to do this to avoid inadvertently weighting towards earlier or later within the time periods. This methodology accounts for the slow global warming in these simulations under net-zero emissions, but results may differ relative to using simulations where global temperature trends are nearer zero for longer, and the same simulations could be used to define the early and late stable GWL ensembles. The length of the millennial-scale simulations and having 40 available SSP5-8.5 simulations result in large ensembles being generated (Table 1). All of the GWL ensembles have average temperatures within 0.1 °C of the target of 1.5, 2, or 3 °C.

The patterns of warming at the 1.5, 2, and 3 °C GWLs under transient (SSP5-8.5), in the early stable period, and in the late stable period were compared. This analysis was performed on boreal summer (June–August or JJA) and boreal winter (December–February or DJF) separately, as previous work has identified differences in warming patterns by season (King et al., 2020, 2021b). Patterns of precipitation change at these GWLs were similarly analysed with percentage differences in seasonal precipitation between the GWLs computed for the transient and the early and late stable period. The significance of the differences was estimated using a Kolmogorov–Smirnov test and a threshold p-value of 0.05. The large sample sizes and low autocorrelation between a given season from one year to the next means that the effective degrees-of-freedom values remain high (Wilks, 2011).

Following this, areas of significant changes in precipitation trends between transient and stabilised climates were examined. This was done by identifying locations where the precipitation distribution is significantly different between the late stable period and SSP5-8.5 for a given GWL. The direction (or sign) of this difference was then compared to the sign of difference between the precipitation distributions during the 1850–1900 baseline and the SSP5-8.5 runs at that GWL. Locations with opposing signs were classified as regions with a significant change in precipitation trend.

Table 1 Number of years and associated cumulative carbon emissions for each ensemble of global warming level by category.

Figure 6 (a) Map of the annual average temperature difference between the years 800–999 and 200–399 averaged across all net-zero simulations. (b) Map of the average difference between net-zero simulations where emission cessation is delayed 5 years (i.e. NZ2035–NZ2030, NZ2040–NZ2035, NZ2045–NZ2040, NZ2050–NZ2045, NZ2055–NZ2050, and NZ2060–NZ2055) in the years 800–999 post-emission cessation.

Differences between the reference GWLs in the Paris Agreement framework were also calculated. Seasonal average temperature and precipitation differences between the 1.5 and 2 °C ensembles in SSP5-8.5 and in the early stable and late stable periods were computed. Difference maps were compared between the transient and stabilising climates.

2.3.3 Changes in extremes and variability

In addition to looking at mean-state changes as a function of time and GWL, changes in extremes and variability were also examined to demonstrate the utility of these millennium-long simulations under near-stable global temperatures. Firstly, changes in monthly temperature extremes were examined by identifying the hottest calendar month at each grid cell in the 1850–1900 baseline period. As we used the 40-member ensemble of historical simulations, there are 2040 monthly mean temperature values for the hottest calendar month at each grid cell, so the standard deviation was computed and is well constrained with a large sample of data. A temperature threshold of the mean plus 2 times the standard deviation was used to compute the frequency of temperature extremes in the early stable and late stable 1.5 and 2 °C GWLs. The difference in frequency of extremes between GWLs and between the early and late stable periods was compared. This method for examining changes in extremes follows that of Cassidy et al. (2023) and is similar to other previous studies (e.g. Thompson et al., 2022).

Finally, a basic analysis of ENSO frequency and variability changes was conducted. Previous work has suggested that ENSO variability is projected to increase under continued global warming in transient simulations (e.g. Cai et al., 2014, 2015), but a study on longer-term ENSO behaviour under slowed global warming in fixed greenhouse gas concentration simulations found a decrease in ENSO amplitude (Callahan et al., 2021). Here, the Niño-3.4 region sea surface temperatures (SSTs) were averaged over the June–July period in the transient (historical and SSP5-8.5 concatenated) and stabilised simulations. The Niño-3.4 SST was detrended using a centred 15-year moving window in all transient and stabilised time series. The standard deviation of the detrended series was calculated in 100-year blocks which were compiled into four groupings: 20th century (based on historical simulations), 21st century (based on historical and SSP5-8.5 simulations), NZ1-500 (based on the first five 100-year blocks in the stabilised simulations), and NZ501-1000 (based on the last five 100-year blocks in the stabilised simulations). The standard deviation was used to examine changes in ENSO amplitude. Changes in the frequency of El Niño and La Niña events were examined by counting events above +0.5 °C and below −0.5 °C, respectively, in the same 100-year blocks. Kolmogorov–Smirnov test statistics were used to identify significant differences between distributions using a p-value of 0.05.

Figure 7 (a) Map of average timing of local peak 11-year smoothed average annual temperatures relative to the point of emission cessation. (b–e) Time series of 11-year smoothed temperatures relative to the point of emission cessation in each net-zero simulation from NZ2030 (yellow) to NZ2060 (red) in the grid cells over (b) New York City, (c) Kinshasa, (d) Delhi, and (e) Melbourne. Locations of these cities are shown in panel (a). Note that the grid cell for Melbourne is north of the city, so it is a land-only grid cell.

3 Results

3.1 Time-varying projections

The seven 1000-year-long simulations exhibit very slow changes in global mean temperature such that they are suitable for use in examining the effects of climate stabilisation and differences with transient warming (Fig. 1d). After the initial change in the first few decades of the simulations, due to the large decrease in methane concentrations, GMST slowly increases over the remainder of these simulations at a rate of around 0.03–0.05 °C per century (Fig. 1d). This is about 1/40 of the rate of observed global warming over the last 30 years. The lack of long-term global cooling despite reduced atmospheric carbon dioxide concentrations (Fig. 1c) is primarily due to slow ocean processes (Armour et al., 2016; MacDougall et al., 2022). The rate of global warming in these ACCESS-ESM-1.5 simulations is a slower rate of global temperature change than would be achieved with fixed atmospheric carbon dioxide concentrations (Rugenstein et al., 2019), but it is slightly higher than could be achieved within the adaptive emission reduction approach framework, where emissions are allowed to vary with the goal of achieving near-zero trends in global mean surface temperature (Terhaar et al., 2022). Beyond GMST, other changes are apparent. The global average land–ocean temperature difference increases strongly under SSP5-8.5, but this trend slowly reverses in all of the net-zero emission simulations (Fig. 3a). Note that, even after 1000 years under net-zero emissions, an increase in land–ocean temperature contrast relative to pre-industrial levels remains and appears as if it would persist beyond the end of the simulations. This is consistent with Joshi et al. (2008), who found persistence of increased land warming relative to ocean warming over simpler, slowly evolving climate model simulations.

As the planet as a whole and the global ocean surface on average warm very slowly under net-zero emissions, continued increases in global average precipitation are simulated in all seven model runs (Fig. 3b). The global average precipitation rate is a function of both global warming and time (Andrews and Forster, 2010; Mitchell et al., 2016) and is illustrative of the continuing global-scale changes projected under net-zero emissions and global average temperature stabilisation. Global average SSTs are projected to continue increasing (Fig. 3c), consistent with the continued slow global average warming and reducing land–ocean temperature contrast. The deeper ocean is projected to warm even centuries after emission cessation, causing the global average vertically averaged ocean temperature to increase (Fig. 3d). The average temperature throughout the entire ocean increases to the extent that, even in the simulation, where emissions cease at 2030, the average ocean temperature exceeds that of SSP5-8.5 in 2100 by 2351. The difference in average global ocean temperatures grows over time between the simulations highlighting the long-term effects of delay in achieving net-zero emissions. This is expected due to the timescales involved in ocean overturning circulation (Armour et al., 2016).

Arctic sea ice extent at the peak and trough of the seasonal cycle exhibits rapid reductions in the transient SSP5-8.5 simulations. Under net-zero emissions, Arctic sea ice change is relatively small with little consistent sign of deterioration or recovery, following the rapid decline before emissions cease (Fig. 3e, f). While there is little change in sea ice extent over the full length of the net-zero simulations, all simulations show apparent substantial trends of different signs over multi-decadal timescales in Arctic sea ice extent in March and September (Fig. 3e, f), although none show a return to pre-industrial-level sea ice extent. These multi-decadal trends are indicative of the significant internal variability and memory in sea ice extent that can give rise to the appearance of significant trends over shorter simulations, although there is evidence that ESMs generally overestimate persistence in Arctic sea ice extent anomalies (Giesse et al., 2021).

Antarctic sea ice extent also decreases under transient global warming but, unlike in the Arctic, continues to decline through the net-zero simulations as well (Fig. 3g, h). In September, around the peak of the seasonal cycle, Antarctic sea ice extent decreases faster under the simulations where emission cessation is achieved later in the 21st century. While many CMIP6 models suffer from significant model biases in sea ice properties, ACCESS-ESM-1.5 performs well in capturing the mean and variability in Antarctic sea ice coverage over the historical period (Roach et al., 2020). Given reasonable model performance, the incidence of low ice, or “ice-free” events, was examined in the Arctic and Antarctic regions using the 1 million km² threshold applied in previous work (Kim et al., 2023; Screen and Williamson, 2017). The likelihood of sea-ice-free Septembers in the Arctic is predominantly a function of global warming level, with little change through the net-zero simulations (Fig. 4a). In contrast, the occurrence of sea-ice-free March events in Antarctica is a function of both global warming and time (Fig. 4b). Under net-zero emissions, in the Antarctic, there is an increasing chance of sea-ice-free events over time to the point that there are more sea-ice-free events in the last century of NZ2030 than in the first century of NZ2060, despite GMST being more than 1 °C higher in the latter.

Figure 8 Warming in boreal summer (JJA) average temperatures at (a) 3 °C, (d) 2 °C, and (g) 1.5 °C global warming levels in SSP5-8.5 simulations relative to the 1850–1900 baseline. The difference in JJA average temperatures between global warming levels extracted between 100 and 450 years after emission cessation and SSP5-8.5 at (b) 3 °C, (e) 2 °C, and (h) 1.5 °C global warming levels. The difference in JJA average temperatures between global warming levels extracted between 650 and 1000 years after emission cessation and SSP5-8.5 at (c) 3 °C, (f) 2 °C, and (i) 1.5 °C global warming levels. (b, c, e, f, h, i) Stippling shows where distributions are not significantly different at the 5 % level using a Kolmogorov–Smirnov test.

Time series analysis of global indices, including land–ocean temperature difference and sea ice extent, suggests that there are complex changes occurring through the net-zero simulations despite minimal change in GMST. To explore this further, we examine local-scale temperature changes. Under transient global warming, a clear pattern of faster warming over land and strong Arctic amplification is identified, whereas warming in the Southern Ocean is very slow (Fig. 5a). As the global temperature stabilises under net-zero emissions, this pattern evolves with land cooling relative to the global average and the Southern Ocean warming (Fig. 5b–e). This is consistent with Fig. 3a, but the near-global nature of land cooling relative to GMST is remarkable. The northern Indian subcontinent is an exception, with warming relative to the global mean in the net-zero simulations compared to the transient case. This is likely related to the reduced atmospheric anthropogenic aerosol forcing in this region in the net-zero emission runs relative to SSP5-8.5 but is also perhaps a response to Southern Ocean warming, as a reduced interhemispheric temperature gradient may reduce water vapour transport to the Indian subcontinent (Oh et al., 2022). The continued Southern Ocean warming also helps to explain the continued decrease in Antarctic sea ice extent under net-zero emissions (Figs. 3g–h and 4). More detailed analysis is needed to understand the role of Southern Ocean warming, atmospheric processes, and ice sheet processes in the large sea ice loss around Antarctica.

While land is projected to cool and ocean is projected to warm relative to GMST, there is also slow global average warming in these net-zero simulations. This results in projected absolute warming over most of the planet, both land and ocean (Figs. 6a, A2). The warming is weak over most land areas but exceeds 2 °C over much of the sea ice region in the Southern Ocean between the early and late periods after emission cessation. We have also already seen that the effect of a delay in achieving net-zero emissions continues through the millennium-length simulations and even grows in some variables. It is remarkable that in the last 2 centuries of the simulations there are still higher temperatures locally across the planet when there is only a 5-year delay in reaching net-zero emissions (Figs. 6b, A3). The difference is particularly large in the planet's sea ice regions, emphasising the role of ice–albedo feedbacks in amplifying small temperature anomalies. Longer delays in emission cessation result in larger amounts of local warming even centuries afterwards (Fig. A4). These results suggest there are major long-term costs associated with even a short delay in achieving net-zero emissions. We note these 5-year delays are based on SSP5-8.5, where the cumulative emission difference is quite large (Table 1) compared with a lower emission scenario.

Different areas of the world experience markedly varied climate trajectories beyond net-zero emissions. To investigate this in more detail, the timing of peak temperatures relative to the point of emission cessation was plotted (Fig. 7a). For a few areas, primarily in northern and central Africa, we find that temperatures peak before the point of net-zero emissions. Many Northern Hemisphere land regions experience peak warming within a couple of centuries of emission cessation, but in the Southern Hemisphere land regions there is a substantial delay to local temperatures peaking. Unsurprisingly, given the previous results, the Southern Ocean peak warming is consistently in the last couple of centuries of the simulations and is likely limited by the runs ending 1000 years after emission cessation. This is consistent with the slow changes found using other long simulations in previous analyses (Armour et al., 2016; Grose and King, 2023). For individual city locations (Fig. 7b–e) we see relatively small temperature changes projected after net-zero relative to the preceding period, but, for the example of Melbourne (Fig. 7e), there is still as much as 1 °C of post-net-zero warming in the NZ2060 simulation up to 1 millennium after emission cessation.

Figure 9 As in Fig. 8 but for boreal winter (DJF).

3.2 Global warming level projections

The length and slow-evolving nature of these net-zero emission simulations means they are suitable for investigating the implications of different GWLs and their dependence on the rate of global warming. Prior to this study, such analysis was challenging and relied on statistical techniques with assumptions (King et al., 2020, 2021b). Firstly, the pattern of seasonal average temperature changes between the rapid-warming SSP5-8.5 ensemble and the early and late parts of the net-zero simulations (see Fig. 2) was examined.

There are differences in the transient warming pattern between boreal summer (JJA; Fig. 8) and boreal winter (DJF; Fig. 9), with Arctic amplification absent in JJA (Lu and Cai, 2009) and a tendency for greater land surface warming in the summertime hemisphere. Generally, in SSP5-8.5, there is scalability in local temperature changes with global changes such that the warming pattern at 3 °C global warming is approximately double the warming pattern at the 1.5 °C GWL (Seneviratne et al., 2016). There are differences in the warming pattern between the early post-net-zero period and SSP5-8.5, and these grow for the later post-net-zero period (Figs. A5, A6). The differences are also seasonally dependent (King et al., 2020) and larger for higher GWLs. In DJF, the difference between sampling the 3 °C GWL from the later stable period relative to the SSP5-8.5 simulations exceeds more than 2 °C across much of the Southern Ocean and the Arctic (Fig. 9c). In JJA, there are many mid-latitude land regions for which sampling a 3 °C GWL from a transient or stabilised climate results in differences of more than 1.5 °C in local seasonal average temperatures. When GWLs are sampled from stabilised runs as opposed to transient runs, temperatures in the Southern Ocean are substantially higher, with temperatures greater than 2 °C higher across much of the Southern Ocean regardless of season. This is consistent with the ongoing warming and reduction in sea ice extent identified previously.

Previous work has identified precipitation changes under climate stabilisation (Ceppi et al., 2018; Grose and King, 2023; Sniderman et al., 2019). Using these ACCESS-ESM-1.5 millennium-length simulations, we also identify significant differences in precipitation change patterns between transient and stabilising GWL projections (Figs. 10, 11). In some cases this is indicative of transient trends intensifying under stabilisation (such as the projection of increased austral winter precipitation over Antarctica; Fig. 10), while in other cases there is projected reversal, particularly over sub-tropical ocean regions projected to dry under continuing greenhouse gas emissions but which are significantly wetter under prolonged net-zero emissions at the same GWLs (Figs. 10, 11).

Figure 10 Change in boreal summer (JJA) average precipitation at (a) 3 °C, (d) 2 °C, and (g) 1.5 °C global warming levels in SSP5-8.5 simulations relative to the 1850–1900 baseline. The difference in JJA average precipitation between global warming levels extracted between 100 and 450 years after emission cessation and SSP5-8.5 at (b) 3 °C, (e) 2 °C, and (h) 1.5 °C global warming levels. The difference in JJA average precipitation between global warming levels extracted between 650 and 1000 years after emission cessation and SSP5-8.5 at (c) 3 °C, (f) 2 °C, and (i) 1.5 °C global warming levels. (b, c, e, f, h, i) Stippling shows where distributions are not significantly different at the 5 % level using a Kolmogorov–Smirnov test.

Figure 11 As in Fig. 10 but for boreal winter (DJF).

Areas of significant precipitation trend differences are shown in Fig. 12 and suggest very large areas of the world may exhibit some return towards pre-industrial levels of seasonal average precipitation. This includes many areas of the world which have been extensively studied due to concern about drying under continued climate change, including the Mediterranean region (Giorgi and Lionello, 2008; Gudmundsson and Seneviratne, 2016) and much of Australia, where stabilisation at a high GWL is projected to result in a significant June–August precipitation increase (Fig. 12a). In contrast, a commonly identified projection of rainfall increase over the Sahara in boreal summer (Pausata et al., 2020) is contrasted with significant drying as the climate stabilises (Figs. 10, 12). There is a general pattern of larger trends in precipitation at higher GWLs under transient climate change (Figs. 10, 11) and also a greater area of significant differences in trends at higher GWLs (Fig. 12). In some locations, seasonal precipitation distributions similar to the 1850–1900 baseline period are identified several centuries after emission cessation despite these locations having significantly different precipitation distributions under SSP5-8.5 (Fig. A7). Recent work by Dittus et al. (2024) found a similar difference in precipitation changes in the Mediterranean region under transient or stabilising climate states using a different ESM.

Figure 12 Locations of differences in precipitation trends under transient and stabilised climate changes for (a, b) 3 °C, (c, d) 2 °C, and (e, f) 1.5 °C global warming levels in JJA and DJF, respectively. Locations in blue are drier under SSP5-8.5 but are significantly wetter (p<0.05) by 650–1000 years after emission cessation at the same global warming level. Locations in brown are wetter under SSP5-8.5 but are significantly drier (p<0.05) by 650–1000 years after emission cessation at the same global warming level. Percentages of the global surface with significant precipitation changes are shown for significant increase under stabilisation relative to SSP5-8.5 (blue) and significant decrease under stabilisation relative to SSP5-8.5 (brown).

Much of the interest in climate projections based on global warming level was initially focused on identifying and understanding differences in climate between the GWLs referred to in the Paris Agreement: 1.5 and 2 °C above pre-industrial levels. Thus, here we examine whether differences between these GWLs vary depending on using transient or stabilised simulations. The ACCESS-ESM-1.5 model simulates a pattern of transient warming and precipitation changes between 1.5 and 2 °C GWLs (Fig. 13) that is consistent with previous multi-model ensemble-based findings (King et al., 2017; Masson-Delmotte et al., 2018).

Figure 13 Differences in (a) JJA and (d) DJF temperatures and (g) JJA and (j) DJF precipitation between 1.5 and 2 °C global warming levels in SSP5-8.5. Differences between 1.5 and 2 °C temperature and precipitation changes for GWLs drawn from 100–450 years after emissions cease and from SSP5-8.5 are shown in panels (b), (e), (h), and (k). Differences between 1.5 and 2 °C temperature and precipitation changes for GWLs drawn from 650–1000 years after emissions cease and from SSP5-8.5 are shown in panels (c), (f), (i), and (l). Hatching illustrates where corresponding absolute differences between stabilised and transient temperature or precipitation are smaller than absolute differences between 1.5 and 2 °C global warming levels in SSP5-8.5.

As the climate stabilises, the pattern of differences between 1.5 and 2 °C GWLs evolves. Over most land areas, the difference in seasonal average temperatures between GWLs under transient warming is substantially greater than between corresponding stabilised GWLs. Over parts of the Southern Ocean, the local differences in JJA and DJF average temperatures between stabilised 1.5 and 2 °C GWLs are more than double the differences between transient 1.5 and 2 °C GWLs. These results are similar to those identified by King et al. (2020). For seasonal average precipitation, differences between 1.5 and 2 °C GWLs are more subtle, and, for many parts of the world, the effect of sampling from stabilised or transient GWLs is greater than the effect of a 0.5 °C difference in global temperatures. These results highlight the importance of framing projections based on GWLs very clearly to reduce the likelihood of misinterpretation.

3.3 Climate extremes

We have demonstrated the utility of these millennium-long ACCESS-ESM1.5 simulations for understanding mean climate changes and climate state dependence on the rate of global warming. Here, we also briefly examine two further applications of these simulations for analysis of changes in climate extremes and climate variability.

Figure 14 Frequency of exceeding 2σ above the 1850–1900 average temperature in the climatologically hottest month of the year at each location at (a, b) 2 °C GWL drawn from 100–450 years and 650–1000 years after emissions cease, respectively. Panels (d) and (e) are the same but for the 1.5 °C GWL. (c, f) Differences in frequency of extremely hot months between GWLs drawn from 100–450 years and 650–1000 years after emissions cease for 2 and 1.5 °C, respectively. (g, h) Differences in frequency of extremely hot months between 2 and 1.5 °C GWLs drawn from 100–450 years and 650–1000 years after emissions cease, respectively.

For the extremes analysis, the hottest calendar month of the year at each location in the ACCESS-ESM-1.5 historical 1850–1900 period was identified, and the variability in that month's temperatures was used to define a threshold of the mean plus 2 standard deviations. The frequency of extremes above this threshold is unsurprisingly higher in the 2 °C GWL samples than in the 1.5 °C GWL samples, regardless of whether these are drawn from earlier or later in the net-zero simulations (Fig. 14). A higher frequency of extremes is projected in the tropics and low-latitude ocean areas where signal-to-noise ratios are greater (Harrington et al., 2016; Hawkins et al., 2020; Hawkins and Sutton, 2012; Mahlstein et al., 2011). The low-latitude regions with high signal to noise are also where the greatest increase in frequency of temperature extremes is projected between the 1.5 and 2 °C GWLs (Fig. 14g, h). For a given GWL, there is a marked reduction in the frequency of hot months over many low- to mid-latitude regions with stabilisation but an increase in the frequency of hot months over the Southern Ocean (Fig. 14c, f). Larger differences may be identifiable comparing transient and stabilised 1.5 and 2 °C GWLs rather than different times within the net-zero simulations, but there are challenges in robustly estimating the frequency and intensity of extremes from fast-warming simulations (King et al., 2020).

3.4 El Niño–Southern Oscillation

The El Niño–Southern Oscillation (ENSO) is the most important mode of climate variability on interannual timescales, and changes in ENSO could have far-reaching implications given its teleconnections to regional climates and impacts (Lieber et al., 2022; Yeh et al., 2018). Here, we examined how the ENSO amplitude and frequency of El Niño and La Niña events compare between transient and stabilised climate states. We previously noted that ACCESS-ESM-1.5 performs well with respect to general ENSO characteristics (Planton et al., 2021) but less well in some aspects of ENSO diversity (Hou and Tang, 2022), so we analyse changes in ENSO amplitude and frequency which are well simulated. ENSO amplitude, as measured by the Niño-3.4 SST standard deviation (Callahan et al., 2021), is projected to increase significantly under continued global warming (Fig. 15a). In contrast, after emission cessation (in both the earlier and later 500-year blocks), the amplitude of ENSO is not significantly different from that simulated for the 20th century and is significantly lower than in transient 21st-century simulations. This is in broad alignment with the previous literature. Cai et al. (2022) found an increase in ENSO amplitude under high- and low-emission transient simulations in CMIP6, while Callahan et al. (2021) identified decreased ENSO amplitude under climate stabilisation in a multi-model ensemble of simulations forced by fixed-CO₂ concentrations. There is also evidence of temporal variation in ENSO amplitude which differs between models (Maher et al., 2023), so further analysis and study of model dependence would be useful.

Figure 15 (a) Box plots of standard deviation in Niño-3.4 average temperatures in the 20th and 21st centuries in fast-warming simulations and in the first and second 500 years after emission cessation. (b) Box plots of El Niño (red) and La Niña (blue) frequency per century in these simulations.

The frequency of El Niño events follows a similar trajectory, with increased occurrence under a scenario of continued high greenhouse gas emission and reduced frequency after the achievement of net-zero emissions (Fig. 15b). On the other hand, no significant change in the frequency of La Niña events is projected in the 21st-century fast-warming simulations. La Niña event frequency is projected to decrease slightly, albeit non-significantly, under net-zero emissions.

For both ENSO amplitude and frequency, we do not find significant differences between the net-zero emission scenarios or between the first 500 years and second 500 years of the simulations. This would suggest the ENSO amplitude and frequency are responding more to the rate of global warming than the amount. Overall, these projections point to contrasting ENSO characteristics between transient and stabilised warmer worlds, with likely effects on regional climates beyond the tropical Pacific Ocean. Further work is needed to explore the mechanisms for ENSO changes under transient versus net-zero conditions.

4 Summary and conclusions

A new set of 1000-year-long ACCESS-ESM-1.5 simulations under net-zero carbon dioxide emissions has been run to help inform critical decision-making around the long-term implications of our current policy goals. We believe that these simulations are unique at present with respect to length and low warming-level stabilisation, allowing the analysis of the climate states associated with the Paris Agreement GWLs. This framework complements existing modelling efforts, including ZECMIP (Jones et al., 2019) and LongRunMIP (Rugenstein et al., 2019), but these new simulations are integrated for longer, the parent simulation is a scenario-based climate, and there are a range of branching points so that the simulated climates span a range of warming levels on both short and long timescales.

Within this framework, net-zero carbon dioxide emission is insufficient to prevent further global warming in ACCESS-ESM-1.5, as continued slow global warming beyond the point of emission cessation is identified. This is despite continued carbon uptake, primarily by the ocean, that causes a reduction in atmospheric CO₂ concentrations. Initial analysis of these ACCESS-ESM-1.5 simulations suggests that the regional climate response to net-zero emissions would be diverse and evolve over subsequent centuries. This work, along with previous studies, highlights the limitations of GMST-based climate targets, as different places and populations experience different climate changes under both continued global warming (Harrington et al., 2018) and in a post-net-zero climate state.

Various properties of the climate system evolve differently under sustained net-zero emissions. While there may be weak changes in Northern Hemisphere land temperature means and extremes, ocean warming is projected to continue. The Southern Hemisphere high latitudes show particularly slow responses to net-zero emissions, with long-term warming for many centuries projected beyond emission cessation accompanied by continuing Antarctic sea ice decline. Extratropical locations in the Southern Hemisphere, such as Melbourne, are projected to warm for centuries after emission cessation in stark contrast to other land areas. While this study does not explore sea level rise, it is already known that this will continue for many centuries under net-zero emission pathways (Nauels et al., 2019). It is clear that the climate of a post-net-zero world will pose different regional risks across the world compared with the near-term transient warming climate and that humanity must prepare so that these risks may be mitigated. This study also suggests that any delay in achieving net-zero emissions may have long-lasting consequences and make the achievement of low GWLs more challenging without substantial net-negative emissions. This adds to the large body of evidence (IPCC, 2021a) showing benefits to earlier emission reductions and the achievement of net-zero carbon dioxide emissions.

The results shown here are based on the ACCESS-ESM-1.5 model only. In ZECMIP, ACCESS-ESM-1.5 has a slightly higher zero emissions commitment (ZEC) than most ESMs (ranked third out of nine) and atmospheric CO₂ reduces at a slower rate than other ESMs (MacDougall et al., 2020), so comparable experiments with other models would be expected to yield slightly different results. While the ACCESS-ESM-1.5 model performs reasonably well in general (Ziehn et al., 2020), this is a single-model analysis and the results should be interpreted in this context. The capability of ESMs to simulate slow-evolving changes in the Earth system or potential tipping points has been debated (Armstrong McKay et al., 2022), and this is a relevant concern given recent rapid changes observed in Antarctic sea ice (Purich and Doddridge, 2023). The results presented in this study use one of our best available modelling tools to understand future climates under net-zero emissions, but improved understanding of slow climate processes and the potential for sudden-onset changes is needed.

The hope with this model framework is that other groups might consider running similar simulations (King et al., 2021a). As discussed previously, other experiments are being run with different protocols but may provide complementary data for analysis. The framework employed here is highly idealised, so future work running simulations under more plausible scenarios and including prolonged net-negative emissions would be beneficial. Our experiments have instant emission cessation following a high rate of emissions, which is unrealistic, although they provide useful results. It would be helpful if the next phase of CMIP included more plausible scenarios that enabled better understanding of the consequences of net-zero emission pathways (including the effect of delay in achieving net-zero emissions), so a multi-model ensemble may be generated and policymakers would be better informed.

This model ensemble may be used to answer many critical questions about future climate changes projected under net-zero emissions. The long simulations with near-stable global temperatures and the pre-existing large ensemble of historical and SSP5-8.5 simulations support many potential analyses to improve the understanding of climate change, variability, and extremes under slow- and fast-changing climate states. This study attempts to demonstrate the utility of these simulations for exploring changes in extremes and in variability, but further work is needed to comprehensively understand climate changes beyond emission cessation. We encourage the analysis of these simulations to help improve the understanding of changes in other forms of climate variability, teleconnections, extremes, ocean and cryosphere properties, and the carbon cycle. In particular, investigation of the evolution of land and ocean carbon fluxes in the different simulations, and comprehensive examination of the processes of sea ice change identified here, would be useful for understanding the future long-term state of the Earth system.

Appendix A

Figure A1 Cumulative CO₂ emissions since 1750 in the historical and SSP5-8.5 simulations (black) and the 1000-year-long net-zero emission simulations beginning in 2030, 2035, 2040, 2045, 2050, 2055, and 2060 (yellow to red).

Figure A2 Maps of annual average temperature difference between the years 800–999 and 200–399 in (a) NZ2060, (b) NZ2055, (c) NZ2050, (d) NZ2045, (e) NZ2040, (f) NZ2035, and (g) NZ2030.

Figure A3 Maps of annual average temperature differences between net-zero simulations where emission cessation is delayed 5 years in the years 800–999 post-emission cessation for (a) NZ2060–NZ2055, (b) NZ2055–NZ2050, (c) NZ2050–NZ2045, (d) NZ2045–NZ2040, (e) NZ2040–NZ2035, and (f) NZ2035–NZ2030.

Figure A4 Maps of annual average temperature differences between net-zero simulations where emission cessation is delayed in the years 800–999 post-emission cessation for (a) NZ2060–NZ2030, (b) NZ2055–NZ2030, (c) NZ2050–NZ2030, (d) NZ2045–NZ2030, (e) NZ2040–NZ2030, and (f) NZ2035–NZ2030.

Figure A5 Warming in boreal summer (JJA), with average temperatures at (a–c) 3 °C, (d–f) 2 °C, and (g–i) 1.5 °C global warming levels in (a, d, g) SSP5-8.5 simulations, (b, e, h) net-zero simulation warming levels extracted between 100 and 450 years after emission cessation, and (c, f, i) net-zero simulation warming levels extracted between 650 and 1000 years after emission cessation, relative to the 1850–1900 baseline.

Figure A6 As in Fig. A5 but for boreal winter (DJF).

Figure A7 Locations of differences in precipitation trends under transient and stabilised climate changes where the difference is large enough for the stabilised climate to be similar to the 1850–1900 baseline period (p>0.05), for (a, b) 3 °C, (c, d) 2 °C, and (e, f) 1.5 °C global warming levels in JJA and DJF, respectively. Locations in blue are drier under SSP5-8.5 but are significantly wetter (p<0.05) by 650–1000 years after emission cessation at the same global warming level and have experienced a return to 1850–1900 seasonal precipitation variability. Locations in brown are wetter under SSP5-8.5 but are significantly drier (p<0.05) by 650–1000 years after emission cessation at the same global warming level and have experienced a return to 1850–1900 seasonal precipitation variability.