Mass loss of the Greenland ice sheet until the year 3000 under a sustained late-21st-century climate

Abstract

We conduct extended versions of the ISMIP6 future climate experiments for the Greenland ice sheet until the year 3000 with the model SICOPOLIS. Beyond 2100, the climate forcing is kept fixed at late-21st-century conditions. For the unabated warming pathway RCP8.5/SSP5-8.5, the ice sheet suffers a severe mass loss, which amounts to ~ 1.8 m SLE (sea-level equivalent) for the 12-experiment mean, and ~ 3.5 m SLE (~ 50% of the entire mass) for the most sensitive experiment. For the reduced emissions pathway RCP2.6/SSP1-2.6, the mass loss is limited to a two-experiment mean of ~ 0.28 m SLE. Climate-change mitigation during the next decades will therefore be an efficient means for limiting the contribution of the Greenland ice sheet to sea-level rise in the long term.

1. Introduction

It is established scientific consensus that the Earth's climate system is warming, and that human influence has been the dominant cause of the observed warming since the mid-20th century (e.g. IPCC, 2021). A major consequence of global warming is sea-level rise, currently (for the period 2006–2018) occurring at a global mean rate of 3.69 ± 0.48 mm a⁻¹. The main sources are melting/discharge of ice sheets, ice caps and glaciers ($\sim {}45\percnt$), thermal expansion of ocean water ($\sim {}38\percnt$), and changes in land water storage ($\sim {}17\percnt$) (Fox-Kemper and others, 2021). In the long term, the two ice sheets of Antarctica (AIS) and Greenland (GrIS) are the largest potential contributors to global sea-level rise because of their enormous volumes, together amounting to ~ 65 m SLE (sea-level equivalent) (Morlighem and others, 2017, 2020). The ice sheets have therefore been the focus of intensive observational as well as modelling efforts.

The Coupled Model Intercomparison Project Phase 6 (CMIP6) is a major international climate modelling initiative (Eyring and others, 2016). As a part of this project, the Ice Sheet Model Intercomparison Project for CMIP6 (ISMIP6) brought together a consortium of ice-sheet modellers to explore the sea-level-rise contribution from the GrIS and AIS (Nowicki and others, 2016, 2020). ISMIP6 focussed on the CMIP6 period from 2015 until the end of 2100. The main findings for the GrIS, when forced by output from CMIP5 global climate models (GCMs), were contributions of 90 ± 50 and 32 ± 17 mm SLE for the unabated warming pathway RCP8.5 (RCP: Representative Concentration Pathway) and the reduced emissions pathway RCP2.6, respectively (Goelzer and others, 2020). The CMIP6 GCMs tend to feature a warmer atmosphere, which results in higher mass loss due to increased surface melt (Payne and others, 2021). For the AIS and CMIP5 climate forcings, ISMIP6 found a mass loss in the range of −7.8 to 30.0 cm SLE under RCP8.5 (Seroussi and others, 2020). The limited number of results for RCP2.6 fall within this range, and so do the results obtained with CMIP6 climate forcings (Payne and others, 2021). This rather unclear picture for the AIS is a consequence of the counteracting effects of mass loss due to ocean warming and mass gain from increased snowfall.

The full suite of ISMIP6 experiments with both CMIP5 and CMIP6 forcings was carried out with the ice-sheet model SICOPOLIS (SImulation COde for POLythermal Ice Sheets, www.sicopolis.net), as documented in detail by Greve and others (2020a, 2020b). Chambers and others (2022) extend the ISMIP6 simulations for the AIS with SICOPOLIS until the year 3000, assuming a sustained late-21st-century climate beyond 2100. Compared to the uncertain response projected over the ISMIP6 period, a radically different picture emerges, demonstrating that the consequences of the high-emissions scenario RCP8.5/SSP5-8.5 (SSP: Shared Socioeconomic Pathway) are much greater in the long term even if no further climate trend is applied beyond 2100.

The response of the GrIS to longer-term climate change has also been investigated. In addition to an ensemble of projections for the 21st century with their higher-order ice-sheet model, the study by Fürst and others (2015) conducted projections until 2300 for the two low-emission scenarios RCP2.6 and RCP4.5, forced by selected CMIP5 GCMs. Vizcaino and others (2015) carried out simulations until 2300 with a coupled ECHAM5.2/MPI-OM/SICOPOLIS model for the pathways RCP2.6, RCP4.5, and a modified RCP8.5 with a 4 × CO₂ limit. Calov and others (2018) drove an extended version of SICOPOLIS (IGLOO) with RCP4.5 and RCP8.5 surface temperature and surface mass-balance anomalies created by the regional climate model MAR with boundary conditions from simulations with three CMIP5 GCMs. Similar to our approach, prolongation of the climatic forcing beyond 2100 was done by assuming no further warming trend. Aschwanden and others (2019) used projections based on four CMIP5 GCMs until 2300 for RCP2.6, RCP4.5 and RCP8.5, extrapolated until 3000, to force the ice-sheet model PISM. Their climatic forcing was processed by the regional climate model HIRHAM5, delivering the surface-temperature anomaly as the main driver, from which precipitation changes were parameterized, and runoff was computed by a positive-degree-day (PDD) method. Van Breedam and others (2020) projected the response of the GrIS and AIS 10000 years into the future with the Earth system model of intermediate complexity LOVECLIMv1.3 (including the ice-sheet model AGISM), forced by the extended concentration pathways ECP2.6, 4.5, 6.0 and 8.5 until 2300 and zero emissions thereafter. In the present study, we transfer the approach by Chambers and others (2022) to the GrIS. The objective is to assess its long-term response to late-21st-century climatic conditions for the full ensemble of ISMIP6 climate forcings, which consists of 14 scenarios from ten different CMIP5 and CMIP6 GCMs.

2. Methods

The main tool used for this study is the ice-sheet model SICOPOLIS. We apply it to the GrIS with hybrid shallow-ice–shelfy-stream dynamics (Bernales and others, 2017), a Weertman–Budd-type sliding law tuned separately for 20 different regions (Greve and others, 2020b), and ice thermodynamics treated by the one-layer melting-CTS enthalpy scheme (CTS: cold-temperate transition surface; Blatter and Greve, 2015; Greve and Blatter, 2016). The horizontal resolution is 5 km. In the vertical, we use terrain-following coordinates (sigma transformation) with 81 layers in the ice domain and 41 layers in the thermal lithosphere layer below. For details on the set-up, the initialization procedure by a paleoclimatic spin-up, comparisons between the simulated and observed ice thickness and surface velocity for our initialization year 1990, as well as the historical run (hist) that bridges the gap between 1990 and the start date of the projections in January 2015 by employing MIROC5/RCP8.5 surface mass balance (SMB) and surface temperature (ST) forcing, we refer to Greve and others (2020b).

Following the ISMIP6 protocol, climate forcing from 2015 until the end of 2100 has an atmospheric and an oceanic component. The atmospheric forcing consists of a 1960–1989 reference climatology, plus space-time-dependent anomalies for SMB and ST. These were derived from a systematic sampling of CMIP5 GCMs that reflects their spread in future projections (Barthel and others, 2020), while CMIP6 GCMs were added on the basis of availability only (Payne and others, 2021). All GCM forcings were downscaled to the GrIS surface with the regional climate model MAR v3.9.6 (Fettweis and others, 2017; Delhasse and others, 2020). Although MAR uses a static GrIS, it also provides vertical gradients of SMB and ST, thus allowing to include SMB–height and ST–height feedbacks in the ice-sheet simulations (Franco and others, 2012; Nowicki and others, 2020). The oceanic forcing is based on a retreat parameterization for tidewater glaciers, forced by MAR runoff and ocean temperature changes specified for seven ice–ocean sectors around Greenland (Slater and others, 2019, 2020).

For the period from 2101 until the end of 3000, we extend the simulations in a similar way than Chambers and others (2022) do for the AIS. For every year of this extended period, the atmospheric forcing (SMB, ST, vertical gradients) for the 10-year interval 2091–2100 is randomly sampled such that no further trend is applied, but some inter-annual fluctuations remain (similar to Calov and others, 2018). The oceanic forcing (prescribed retreat maps) does not show any notable year-to-year fluctuations, so we simply keep it fixed at 2100 conditions.

An overview of our extended ISMIP6 experiments is given in Table 1, and the magnitude of the atmospheric forcing is shown in Supplementary Table S1. Twelve experiments are for the 21st-century unabated warming pathway RCP8.5 (CMIP5) / SSP5-8.5 (CMIP6), and two are for the reduced emissions pathway RCP2.6 (CMIP5) / SSP1-2.6 (CMIP6) that is largely in line with the commitments of the Paris Agreement (maintaining the global mean temperature well below a 2^°C increase above pre-industrial levels). In two experiments, the impact of different sensitivities of the retreat parameterization due to oceanic forcing (‘high’ and ‘low’ vs the normal, ‘medium’ sensitivity, thereby exploring the uncertainty of the parameterization; Slater and others, 2019, 2020) is tested. In addition, a projection control run for the period 2015–3000 (ctrl_proj) employs constant climate conditions based on a 1960–1989 climatology and no explicit oceanic forcing.

Table 1. Extended ISMIP6-Greenland Tier-1 and 2 future climate experiments discussed in this study

See Nowicki and others (2020) for references for the GCMs.

3. Results

The simulated mass change of the GrIS, expressed as a sea-level contribution, and ice area are shown in Figure 1. For the control run ctrl_proj, the ice sheet remains nearly stable, showing a slight mass gain of 6.4 mm SLE and area loss of 4.7 × 10³ km² during the 986 years model time, which is of the order of permilles of the present-day values. For all future projections, the ice sheet keeps losing both mass and extent over the entire period. (Values of the resulting sea-level contribution by 2100, 2300 and 3000 relative to ctrl_proj are provided in the Supplementary Table S2.) The largest rate of change occurs typically around the year 2100, beyond which it slows down to some extent; however, without reaching or coming close to a new steady state. This demonstrates that the committed mass loss due to 21st-century climate change extends way beyond the 21st century and impacts the ice sheet on a much longer time scale. Corroborating the findings for the 21st century (Goelzer and others, 2020; Greve and others, 2020b), the GrIS responds much more strongly to the ensemble of RCP8.5/SSP5-8.5 simulations than to the two RCP2.6/SSP1-2.6 simulations. By the year 3000, the mass loss amounts to 1.79 ± 0.80 m SLE (mean ± 1-sigma range) for RCP8.5/SSP5-8.5, while it is limited to 0.28 ± 0.12 m SLE for RCP2.6/SSP1-2.6.

Fig. 1. Extended ISMIP6-Greenland historical run (hist, 1990–2015), projection control run (ctrl_proj) and Tier-1 and 2 future climate experiments: (a) simulated ice mass change (counted positively for loss and expressed as sea-level contribution), (b) ice area. The red and blue boxes to the right show the $\hbox {mean}\pm {}1\hbox {-sigma}$ ranges for RCP8.5/SSP5-8.5 and RCP2.6/SSP1-2.6, respectively; the whiskers show the corresponding full ranges.

The percentage area loss (Fig. 1b) is similar to the mass loss. However, it is striking that the variability in area varies significantly between the different simulations. This results from the different variability of the SMB forcing and affects mainly the thin, near-margin parts of the ice sheet, which do not contribute much to the total ice mass. Therefore, the variability is not paralleled in Figure 1a.

The influence of the ice retreat due to oceanic forcing is explored by Exps. 5, 9, 10 (MIROC5/RCP8.5 with ‘medium’, ‘high’ and ‘low’ sensitivity, respectively). The results are shown by the olive lines and olive-shaded regions in Figure 1. By 3000, the simulated mass loss is $1.62^{ + 0.051}_{ - 0.031}\, {\rm m\, SLE}$. Thus, the uncertainty due to these three calibrations is very small in the long range. Relative to the uncertainty due to the different climate forcings, it is more pronounced for the 21st century (Greve and others, 2020b). This is because the continued retreat of the ice sheet decreases its contact with the ocean, so that the oceanic forcing plays a smaller role in the longer term.

As reported by Greve and others (2020b) and Payne and others (2021), for both the 21st-century RCP8.5/SSP5-8.5 and RCP2.6/SSP1-2.6 pathways, the CMIP6 climate models produce a larger response of the ice sheet than the CMIP5 ones. While the significance of this statement is limited in the case of RCP2.6/SSP1-2.6 (only one experiment each), it is more robust for RCP8.5/SSP5-8.5, where the ensemble contains eight and four experiments forced by CMIP5 and CMIP6 models, respectively. By 3000, the mean mass loss for the four CMIP6 SSP5-8.5 experiments is 2.73 m SLE, and the maximum value from Exp. B4 (CESM2/SSP5-8.5) is as large as 3.54 m SLE, almost 50% of the entire present-day ice mass.

We now discuss in more detail the results of Exp. 5 (MIROC5/RCP8.5), which was already focused on in the original ISMIP6-Greenland study (Goelzer and others, 2020). Mainly due to its large SMB forcing, it produces, along with Exp. A1 (IPSL-CM5A-MR/RCP8.5), the strongest mass loss among the CMIP5 forcings, while the mass loss is about average for our combined CMIP5/CMIP6 ensemble. Figure 2 shows the components of the global mass balance (integrated over the ice sheet, all counted as positive for mass gain): surface mass balance (SMB), basal mass balance (BMB), calving and ice volume change (dV/dt). On a mean-annual basis, the residual, Res = |SMB + BMB + Calving − dV/dt|, is always <$10^{6}\, {\rm m^{3}\, a^{-1}}$. This is five to six orders of magnitude smaller than the typical range of values in the figure, so that the model conserves mass very well (see also Calov and others, 2018).

Fig. 2. Main components of the global mass balance for Exp. 5 (MIROC5/RCP8.5): surface mass balance (SMB, purple), basal mass balance (BMB, blue), calving (yellow) and ice volume change (dV/dt, green). Note the shifted, right axis for the latter. The black and green dashed lines indicate the zero levels for the left and right axis, respectively.

As already stated above, the ice sheet keeps losing volume (∝ mass) over the entire period, with maximal rates of change occurring shortly before the year 2100. The SMB is initially positive, but changes its sign in the second half of the 21st century and stays negative beyond that. Calving into the surrounding ocean peaks during 2080–2085, when it contributes approximately the same amount to ice volume loss than negative SMB. After that, calving decreases continuously due to ice-sheet retreat from the coast (loss of ocean contact) and becomes almost negligible towards the end of the 3rd millennium. The decrease is likely accelerated by a limitation of the oceanic forcing approach: the fixed retreat mask beyond 2100 does not follow the retreating ice margin further, making it ineffective as the ice sheet recedes beyond its reach. The negative SMB shows a slightly decreasing trend after 2100 because, as the ice sheet shrinks, less area is available for further melting. Both effects together cause the magnitude of dV/dt (loss rate of ice volume) to decrease gradually, which results in the slight flattening of the curves in Figure 1 (the same mechanisms are effective for the other experiments). BMB is relatively small over the entire model time. The inter-annual variability of the volume change is due to that of the SMB, which reflects the variability of the atmospheric forcing.

Figure 3 shows snapshots of the ice thickness and surface velocity for Exp. 5 for the initial year 2015 and the final year 3000. Comparing the thickness distributions demonstrates nicely that the ice sheet retreats from the coast almost all around its perimeter, and contact to the ocean is very limited by the end of the simulation, which entails the low calving rates mentioned above. By contrast, the ice sheet does not suffer much change in its interior parts north of ~ 68^°N. The large-scale pattern of the ice flow and the organization of the ice sheet into major drainage basins remain largely intact. However, on the regional scale, more pronounced changes occur. The fast-flowing outlet glaciers in south-western Greenland disappear entirely due to the extreme retreat in this area. The north-western outlet glaciers, including Petermann Glacier, also slow down substantially. The central-western Jakobshavn Ice Stream loses its clear delimitation to the surrounding glaciers, but remains an area of fast-flowing ice. The major features in East Greenland, e.g. the North-East Greenland Ice Stream, Kangerdlugssuaq and Helheim glaciers, are less affected and remain well identifiable.

Fig. 3. Ice thickness (panels a, b) and surface velocity (panels c, d) for the initial time (2015; panels a, c) and final time (3001; panels b, d) of Exp. 5 (MIROC5/RCP8.5).

4. Discussion and conclusion

The future climate simulations carried out in this study for the GrIS over the 3rd millennium confirm and continue the trends that were reported by ISMIP6-Greenland for the 21st century (Goelzer and others, 2020; Greve and others, 2020b; Payne and others, 2021). The response of the ice sheet is mainly governed by a negative SMB due to increased surface melting near the ice margin. Marine-terminating glacier retreat, caused by increased oceanic thermal forcing and increased meltwater runoff, constitutes a further negative contribution to the total mass balance, but becomes less important in the longer term. Under the unabated warming pathway RCP8.5/SSP5-8.5, this leads to a severe mass loss during the 3rd millennium, while the loss is much smaller under the reduced emissions pathway RCP2.6/SSP1-2.6. Results obtained with forcings from the newer CMIP6 climate models consistently produce larger mass losses than those obtained with the older CMIP5 models, for SSP5-8.5 in the range of a ~ 25–50% loss of the present-day ice mass (and area) by 3000. For comparison, Aschwanden and others (2019) reported a mass loss of ~ 75–100% by 3000 for their ensemble of RCP8.5 simulations, for which a warming trend is assumed to continue until 2500. Efficient climate change mitigation during the next decades is therefore crucial for limiting the contribution of the GrIS to long-term sea-level rise.

As for interpreting the stronger response of the GrIS to CMIP6-derived forcings compared to CMIP5-derived ones, the different strategies of sampling the respective GCM ensembles for ISMIP6-Greenland must be considered. As we mentioned in Section 2, the six CMIP5 GCMs were chosen systematically, whereas the four CMIP6 GCMs were the only ones available for downscaling at the time. Subsequent analysis of the entire CMIP6 model ensemble revealed that the results cluster around two groups of climate sensitivities (global mean temperature response to doubled atmospheric CO₂), and the four models that were available for ISMIP6-Greenland all fall in the high-sensitivity group (Meehl and 7 others, 2020; Payne and others, 2021). Future work should therefore aim at conducting additional simulations with a more representative sampling of the CMIP6 GCMs.

Our study is limited to investigating the impact of a sustained late-21st-century climate (without imposing a further trend beyond 2100) on the GrIS. However, climate change is projected to continue beyond 2100 (e.g. Bakker and others, 2016; Lyon and others, 2022), with potentially even more devastating effects on the GrIS than reported here, plus a significant decay of the AIS in the long term (Van Breedam and others, 2020). Further, the unidirectional coupling approach (climate model → ice-sheet model) employed by ISMIP6, and thus here, lacks a detailed accounting of feedbacks of the changing ice sheet on the climate. As we explained in Section 2, the climate forcing for Greenland includes vertical gradients of the surface mass balance and surface temperature. Therefore, the changing ice-sheet geometry acts back on these climatic forcing fields. However, the linearized approach was derived for small perturbations of the present-day state only, and it cannot be validated for large changes of the ice sheet. This shortcoming becomes more severe in our simulations over almost a millennium compared to the 86-year scope of ISMIP6, and it is not possible to judge a priori whether it rather leads to an over- or underestimation of the simulated mass loss. Using the PDD method for modelling runoff, like in the study by Aschwanden and others (2019), allows in principle to handle arbitrarily large changes of the ice sheet. PDD is also easy to implement in an ice-sheet model, which adds to its appeal. However, it is a grossly simplifying parameterization of the complex energy balance at the surface of an ice sheet and difficult to calibrate for climatic conditions other than present-day ones (e.g. Bauer and Ganopolski, 2017), so that we do not consider it superior to the ISMIP6-type approach of our study.

Future work in the direction of long-term simulations of ice-sheet response to climate change should aim at employing GCM projections beyond 2100 and improving the representation of feedback processes. The ultimate solution would be to carry out such simulations in a fully coupled way, with the ice-sheet model integrated in the GCM. This approach has been pursued (e.g. Vizcaino and others, 2015; Sellevold and others, 2019; Gregory and others, 2020; Muntjewerf andand others, 2020; Muntjewerf and others, 2020); however, fully coupled simulations are demanding and computationally expensive, which makes it difficult to run large ensembles, involving many different climate and ice-sheet models, over long time scales and at adequate resolution. Intermediary, more manageable solutions may consist of involving snapshots of climate-model results combined with more refined parameterizations for the climatic forcing, similar to the approach by Abe-Ouchi and others (2013) for the paleo-glaciation of the Northern Hemisphere.