Ice Flow Model Research Articles

Creep enhancement and sliding in a temperate, hard-bedded alpine glacier

Abstract. Glacier internal deformation is usually described by Glen's flow law using two material parameters: the creep factor (A) and the flow law exponent (n). However, the values of these parameters and their spatial and temporal variability are rather uncertain due to the difficulty in quantifying internal strain and stress fields at natural scales. In this study, we combine 1-year-long continuous measurements of borehole inclinometry and surface velocity with three-dimensional full-Stokes ice flow modeling to infer ice rheologies and sliding velocities for the ablation zone of the Argentière Glacier, a temperate glacier in the French Alps. We demonstrate that the observed deformation rate profile has limited sensitivity to the flow law exponent (n) and instead mainly reflects an increase in the creep factor (A) with depth, with A departing from its surface value by up to a factor of 2.5 below 160 m depth. We interpret this creep factor enhancement as an effect of increasing interstitial water content with depth (from 0 % to 1.3 %), which results in an average value of A=148 MPa−3 a−1. We further observe that internal ice deformation exhibits seasonal variability similar to that concerning surface velocity, indicating that the local basal sliding velocity exhibits no significant seasonal variation. We suggest that these changes in deformation rate are due to variations in the stress field, driven by contrasting changes in subglacial hydrology conditions between the sides and center of the glacier. Our study provides further evidence that borehole inclinometry, combined with full-Stokes flow modeling, allows for the constraining of both ice rheology and basal friction at scales that cannot be inferred from surface velocity measurements alone.

Quantifying the buttressing contribution of landfast sea ice and melange to Crane Glacier, Antarctic Peninsula

Abstract. The January 2022 disintegration of multi-year landfast sea ice in the Larsen B embayment, Antarctic Peninsula, was closely followed by a significant acceleration of ice flow and ice-front retreat of numerous outlet glaciers. Crane Glacier was a notable example of this, with 6 km of its floating ice shelf lost to calving in the first month following the disintegration and a 3.4 % increase in terminus flow speeds over the same time period. In this study we quantify for the first time the buttressing stresses that were transmitted to Crane by the ice melange at the glacier outlet using the ice-flow model Úa. We constrained our model with high-resolution surface elevation profiles of the glacier and ambient melange and reconstructed the observed flow velocities by optimising the rheology rate factor throughout our model domain. This allowed us to quantify the stress regime across both the glacier and ice melange. Results showed that resistive backstresses were imparted to Crane by the ice melange with a mean buttressing ratio of ΘN=0.68 calculated at the glacier terminus (ΘN=1 implies no buttressing). In addition, diagnostic modelling showed an expected 19.2 kPa mean increase in extensional stress at the ice front following the disintegration of the ice melange. This perturbation in stress likely triggered the observed rapid calving over the near-terminus region, leading to the periodic loss of sections of Crane's buttressing ice shelf and thus further acceleration of ice flow in the subsequent months.

Two-way coupling between ice flow and channelized subglacial drainage enhances modeled marine-ice-sheet retreat

Abstract. Ice-sheet models used to predict sea-level rise often neglect subglacial hydrology. However, theory and observations suggest that ice flow and subglacial water flow are bidirectionally coupled: ice geometry affects hydraulic potential, hydraulic potential modulates basal shear stress via the basal water pressure, and ice flow advects the subglacial drainage system. This coupling could impact rates of ice mass change but remains poorly understood. We develop a coupled ice–subglacial-hydrology model to investigate the effects of coupling on the long-term evolution of marine-terminating ice sheets. We combine a one-dimensional channelized subglacial hydrology model with a depth-integrated marine-ice-sheet model, incorporating each component of the coupling listed above, yielding a set of differential equations that we solve using a finite-difference, implicit time-stepping approach. We conduct a series of experiments with this model, using either bidirectional or unidirectional coupling. These experiments generate profiles of channel cross-sectional area, channel flow rate, channel effective pressure, ice thickness, and ice velocity. We discuss how the profiles shape one another, resulting in the effective pressure reaching a local maximum in a region near the grounding line. We also describe the impact of bidirectional coupling on the transient retreat of ice sheets through a comparison of our coupled model with ice-flow models that have imposed static basal conditions. We find that including coupled subglacial hydrology leads to grounding-line retreat that is virtually absent when static basal conditions are assumed. This work highlights the role time-evolving subglacial drainage may have in ice-sheet change and informs efforts to include it in ice-sheet models. This work also supplies a physical basis for a commonly used parameterization which assumes that the subglacial water pressure is set by the bed's depth beneath the sea surface.

Precipitation drives western Patagonian glacier variability and may curb future ice mass loss.

Two-thirds of all glaciers worldwide are projected to disappear by 2100 CE. Large uncertainties however remain in maritime settings, where some glaciers have recently gained mass in response to increased snowfall. One of these regions is southern Patagonia, where increased precipitation since the 1980s seems to have attenuated glacier retreat. Whether this exceptional behavior will continue in a warmer future is unclear. Here, we use a numerical ice-flow model constrained by paleoglaciological data to simulate how climate variability influenced the evolution of three maritime outlet glaciers of the Southern Patagonian Icefield during the last 6000 years. Our experiments suggest that precipitation drove 67% of the centennial-scale fluctuations in the volume of the modeled glaciers. When applied to the temperature projected by 2100 CE, our simulations show that precipitation needs to increase by 10-50% to maintain present-day glacier volumes, depending on the climate scenario (SSP1-2.6 to SSP5-8.5). This implies that if greenhouse-gas emission cuts fail, these glaciers will enter a warmer regime unseen over the last 6000 years, where precipitation cannot offset glacier mass loss. Conversely, if emissions are curtailed, increased precipitation could halt mass loss of some of Patagonia's largest glaciers, and potentially of other maritime glaciers worldwide.

Melt sensitivity of irreversible retreat of Pine Island Glacier

Abstract. In recent decades, glaciers in the Amundsen Sea Embayment in West Antarctica have made the largest contribution to mass loss from the entire Antarctic Ice Sheet. Glacier retreat and acceleration have led to concerns about the stability of the region and the effects of future climate change. Coastal thinning and near-synchronous increases in ice flux across neighbouring glaciers suggest that ocean-driven melting is one of the main drivers of mass imbalance. However, the response of individual glaciers to changes in ocean conditions varies according to their local geometry. One of the largest and fastest-flowing of these glaciers, Pine Island Glacier (PIG), underwent a retreat from a subglacial ridge in the 1940s following a period of unusually warm conditions. Despite subsequent cooler periods, the glacier failed to recover back to the ridge and continued retreating to its present-day position. Here, we use the ice-flow model Úa to investigate the sensitivity of this retreat to changes in basal melting. We show that a short period of increased basal melt was sufficient to force the glacier from its stable position on the ridge and undergo an irreversible retreat to the next topographic high. Once high melting begins upstream of the ridge, only near-zero melt rates can stop the retreat, indicating a possible hysteresis in the system. Our results suggest that unstable and irreversible responses to warm anomalies are possible and can lead to substantial changes in ice flux over relatively short periods of only a few decades.

Inland Summer Speedup at Zachariæ Isstrøm, Northeast Greenland, Driven by Subglacial Hydrology

AbstractThe Northeast Greenland Ice Stream (NEGIS) has experienced substantial dynamic thinning in recent years. Here, we examine the evolving behavior of NEGIS, with focus on summer speedup at Zachariae Isstrøm, one of the NEGIS outlet glaciers, which has exhibited rapid retreat and acceleration, indicative of its vulnerability to changing climate conditions. Through a combination of Sentinel‐1 data, in‐situ GPS observations, and numerical ice flow modeling from 2007, we investigate the mechanisms driving short‐term changes. Our analysis reveals a summer speedup in ice flow both near the terminus and inland, with satellite data detecting changes up to 60 km inland, while GPS data capture changes up to 190 km inland along the glacier center line. We attribute this summer speedup to variations in subglacial hydrology, where surface meltwater runoff influences basal friction over the melt season. Incorporating subglacial hydrology into numerical models makes it possible to replicate observed ice velocity patterns.

Origin of ice in the Medusae Fossae Formation, equatorial Mars

We assess and test the recent radar-sounding findings of huge, kms-thick quantities of ice along the equatorial region of Mars covered by a thick protective cap (the Medusae Fossae Formation (MFF); Watters et al., 2024), using 1) atmospheric general circulation models (GCMs), 2) glacial ice accumulation and flow models, and 3) models for ice ablation-induced accumulation residues. Our results indicate that under Hesperian-era conditions, Mars at ∼40° spin-axis obliquity is predicted to accumulate snow/ice in excess of a km in less than a few million years in the MFF region, producing cold-based glaciers with basal melting over <7% of the deposit. We find that subsequent ablation of this ice deposit surface in the several billion year-long Amazonian during periods of episodic eolian stripping of MFF protective sublimation residues and dust/tephra deposits, provides a mechanism sufficient to form the thick capping layer. Similar, shorter-duration obliquity excursions during this period may have also contributed additional ice-sublation residue layers, consistent with the complex stratigraphy of the MFF. On the basis of the estimated non-ice component of the MFF ice deposit, we suggest that ablation alone could have formed the cap unit in a minimum of ∼350 million years, but was likely to have been episodic, operating in concert with periods of eolian surface lag removal over a much longer time period. The tripartite subdivision of MFF stratigraphy could indicate major periods of Amazonian obliquity excursions that deposited and removed thinner layers of ice and sublimation residue. The very high abundance of non-esker-like fluvial channels in part of the Lower Member of the MFF, combined with the paucity of ice-sheet basal melting in our analysis, suggests that ablation processes were sometimes dominated by top-down ice heating, melting and fluvial runoff. In summary, our three-part modeling approach supports the new findings, and offers new dimensions for the further analysis of the enigmatic MFF.

Weak relationship between remotely detected crevasses and inferred ice rheological parameters on Antarctic ice shelves

Abstract. Over the past decade, a wealth of research has been devoted to the detection of crevasses in glaciers and ice sheets via remote sensing and machine learning techniques. It is often argued that remotely sensed damage maps can function as early warning signals for shifts in ice shelf conditions from intact to damaged states and can serve as an important tool for ice sheet modellers to improve future sea level rise predictions. Here, we provide evidence for the Filchner–Ronne and Pine Island ice shelves that remotely sensed damage maps are only weakly related to the ice rate factor field A derived by an ice flow model when inverting for surface velocities. This technique is a common procedure in ice flow models, as it guarantees that any inferred changes in A relate to changes in ice flow measured through observations. The weak relationship found is improved when investigating heavily damaged shear margins, as observed on the Pine Island Ice Shelf; however, even in this setting, this association remains modest. Our findings suggest that many features identified as damage through remote sensing methods are not of direct relevance to present-day ice shelf flow. While damage can clearly play an important role in ice shelf processes and thus be relevant for ice sheet behaviour and sea level rise projections, our results imply that mapping ice damage directly from satellite observations may not directly help improve the representation of these processes in ice flow models.

Geothermal heat source estimations through ice flow modelling at Mýrdalsjökull, Iceland

Abstract. Geothermal heat sources beneath glaciers and ice caps influence local ice-dynamics and mass balance but also control ice surface depression evolution as well as subglacial water reservoir dynamics. Resulting jökulhlaups (i.e., glacier lake outburst floods) impose danger to people and infrastructure, especially in Iceland, where they are closely monitored. Due to hundreds of meters of ice, direct measurements of heat source strength and extent are not possible. We present an indirect measurement method which utilizes ice flow simulations and glacier surface data, such as surface mass balance and surface depression evolution. Heat source locations can be inferred accurately to simulation grid scales; heat source strength and spatial distributions are also well quantified. Our methods are applied to the Mýrdalsjökull ice cap in Iceland, where we are able to refine previous heat source estimates.

Study on ice resistance of Antarctic krill ship with trawl under floating ice sea conditions

IntroductionThis study focused on a Chinese Antarctic krill vessel utilising continuous pumping fishing technology. The resistance characteristics of Antarctic krill ships trawling in floating ice areas is of great significance for the navigation and fishing of krill ships in ice areas.MethodsFirstly, MATLAB programming using discrete elements combined with genetic algorithms was used to construct a normal distribution ice flow model. Secondly, a fluid-structure coupling interface is created through the contact between the fluid and the trawl grid, and the displacement and resistance of the trawl grid are evaluated on the shared interface. Finally, the effects of ice density and ship sailing speed on ice resistance were studied.Results and discussionThe results of the calculations results show that ice resistance is positively related to the concentration and speed of floating ice, moreover, there is a special speed point where ice resistance increases rapidly. As the speed increases, the proportion of trawl resistance to the total resistance continues to increase, while the proportion of ice resistance continues to decrease. This paper provides a reference for the navigation and fishing resistance assessment of Antarctic krill ships in floating ice areas.

Impact of subsurface crevassing on the depth–age relationship of high-Alpine ice cores extracted at Col du Dôme between 1994 and 2012

Abstract. Three seasonally resolved ice core records covering the 20th century were extracted in 1994, 2004, and 2012 at a nearly identical location from the Col du Dôme (4250 m above sea level, m a.s.l.; Mont Blanc, French Alps) drill site. Here, we complete and combine chemical records of major ions and radiometric measurements of 3H and 210Pb obtained from these three cores with a 3D ice flow model of the Col du Dôme glacier to investigate in detail the origin of discontinuities observed in the depth–age relation of the ice cores drilled in 2004 and 2012. Taking advantage of the granitic bedrock at Col du Dôme, which makes the ice core 210Pb records sensitive to the presence of upstream crevasses, and the fact that the depth–age disturbances are observed at depths for which absolute time markers are available, we draw an overall picture of a dynamic crevasse formation. This can explain the non-disturbed depth–age relation of the ice core drilled in 1994 and the perturbations observed in those drilled in 2004 and 2012. Since crevasses are common at high-Alpine glacier sites, our study points to the important need for rigorous investigations of the depth–age scale and glaciological conditions upstream of drill sites before interpreting high-alpine ice core records in terms of atmospheric changes.

Projections of Peak Water Timing From the East Rongbuk Glacier, Mt. Everest, Using a Higher‐Order Ice Flow Model

AbstractIn this study, we apply a three‐dimensional (3D) thermomechanically coupled higher‐order ice flow model to simulate the East Rongbuk Glacier (ERG), Mt. Everest. We first diagnostically investigate its present‐day ice dynamic features in 2009 and then prognostically simulate the glacier during the time period 2010–2100. The ice flow model is initialized based on a Robin‐type inversion method by conducting six sensitivity experiments relating to glacier thermal boundary conditions. We apply two different surface mass balance parameterizations in the model, and both of them can reproduce the observed ice volume loss (around 0.1 km3) during 2010–2020. We find that ERG is likely to experience maximum meltwater runoff at the year 2030 under the SSP‐126 scenario, while under SSP‐370 and ‐585 scenarios, the peak water will both likely occur at around 2060. The ice dynamics may contribute more to ice loss as climate warms in time.

Failure strength of glacier ice inferred from Greenland crevasses

Abstract. Ice fractures when subject to stress that exceeds the material failure strength. Previous studies have found that a von Mises failure criterion, which places a bound on the second invariant of the deviatoric stress tensor, is consistent with empirical data. Other studies have suggested that a scaling effect exists, such that larger sample specimens have a substantially lower failure strength, implying that estimating material strength from laboratory-scale experiments may be insufficient for glacier-scale modeling. In this paper, we analyze the stress conditions in crevasse onset regions to better understand the failure criterion and strength relevant for large-scale modeling. The local deviatoric stress is inferred using surface velocities and reanalysis temperatures, and crevasse onset regions are extracted from a remotely sensed crevasse density map. We project the stress state onto the failure plane spanned by Haigh–Westergaard coordinates, showing how failure depends on mode of stress. We find that existing crevasse data are consistent with a Schmidt–Ishlinsky failure criterion that places a bound on the absolute value of the maximal principal deviatoric stress, estimated to be 158±44 kPa. Although the traditional von Mises failure criterion also provides an adequate fit to the data with a von Mises strength of 265±73 kPa, it depends only on stress magnitude and is indifferent to the specific stress state, unlike Schmidt–Ishlinsky failure which has a larger shear failure strength compared to tensile strength. Implications for large-scale ice flow and fracture modeling are discussed.

Ice volume and thickness of all Scandinavian glaciers and ice caps

Abstract We present a new map of bed topography and ice thickness together with a corresponding ice volume estimate representative of the years ~2010 for all Scandinavian ice caps and glaciers. Starting from surface observations, we invert for ice thickness by iteratively running an innovative ice dynamics model on a distributed grid and updating bed topography until modelled and observed glacier dynamics as represented by their rate of surface elevation change (dh/dt) fields align. The ice flow model used is the instructed glacier model (Jouvet and Cordonnier, 2023, Journal of Glaciology 1–15), a generic physics-informed deep-learning emulator that models higher-order ice flow with high-computational efficiency. We calibrate the modelled thicknesses against &gt;11 000 ice thickness observations, resulting in a final ice volume estimate of 302.7 km3 for Norway, 18.4 km3 for Sweden and 321.1 km3 for the whole of Scandinavia with an error estimate of ~ $\pm 11\%$ . The validation statistics computed indicate good agreement between modelled and observed thicknesses (RMSE = 55 m, Pearson's r = 0.87, bias = 0.8 m), outperforming all other ice thickness maps available for the region. The modelled bed shapes thus provide unprecedented detail in the subglacial topography, especially for ice caps where we produce the first maps that show ice-dynamically realistic flow features.

Impact of boundary conditions on the modeled thermal regime of the Antarctic ice sheet

Abstract. A realistic initialization of ice flow models is critical for predicting future changes in ice sheet mass balance and their associated contribution to sea level rise. The initial thermal state of an ice sheet is particularly important, as it controls ice viscosity and basal conditions, thereby influencing the overall ice velocity. Englacial and subglacial conditions, however, remain poorly understood due to insufficient direct measurements, which complicate the initialization and validation of thermal models. Here, we investigate the impact of using different geothermal heat flux (GHF) datasets and vertical velocity profiles on the thermal state of the Antarctic ice sheet and compare our modeled temperatures to in situ measurements from 15 boreholes. We find that the temperature profile is more sensitive to vertical velocity than to GHF. The basal temperature of grounded ice and the amount of basal melting are influenced by both selection of GHF and vertical velocity. More importantly, we find that the standard approach, which consists of combining basal sliding speed and incompressibility to derive vertical velocities, provides reasonably good results in fast-flow regions (ice velocity &gt;50 m yr−1) but performs poorly in slow-flow regions (ice velocity &lt;50 m yr−1). Furthermore, the modeled temperature profiles in ice streams, where bed geometry is generated using a mass conservation approach, show better agreement with observed borehole temperatures compared to kriging-based bed geometry.

Palaeoglacier and palaeoclimate reconstruction since the Neoglacial of the central Gangdise Mountains, southern Tibetan Plateau

Palaeoglacier modelling is a crucial approach for reconstructing palaeoclimate. Although the timing of glaciations during the Little Ice Age (LIA) and Neoglacial in the central Gangdise Mountains has been partially determined, there are still undated fresh-looking moraines. This study utilized coupled mass balance and ice flow models to simulate the palaeoclimate during the five LIA glacial events (LIA-1 to LIA-5) and two Neoglacial glacial events (Neo-1 and Neo-2). The results revealed that the glacial volumes in the Gangdise Mountains during the LIA and Neoglacial periods ranged from ∼1.10 to 1.57 km3. Comparing the reconstructed climate with other proxies, it was found that Neo-2 and Neo-1 exhibited temperatures 0.54 to 0.82 °C and 0.7 to 1.0 °C lower than the present, respectively, with precipitation levels at 110–140% of the present. The temperatures during the LIA periods were 0.36 to 1 °C lower than the present, with precipitation levels at 80–120% of the present. The interannual variability of climate at decadal or centennial scales may have played a significant role in the formation of the LIA moraines.

Validating ensemble historical simulations of Upernavik Isstrøm (1985–2019) using observations of surface velocity and elevation

Abstract The future of tidewater glaciers in response to climate warming is one of the largest sources of uncertainty in the contribution of the Greenland ice sheet to global sea-level rise. In this study, we investigate the ability of an ice-sheet model to reproduce the past evolution of the velocity and surface elevation of a tidewater glacier, Upernavik Isstrøm, by prescribing front positions. To achieve this, we run two ensembles of simulations with a Weertman and a regularised-Coulomb friction law. We show that the ice-flow model has to include a reduction in friction in the first 15 km upstream of the ice front in fast-flowing regions to capture the trends observed during the 1985–2019 period. Without this process, the ensemble model overestimates the ice flow before the retreat of the front in 2005 and does not fully reproduce its acceleration during the retreat. This results in an overestimation of the total mass loss between 1985 and 2019 of 50% (300 vs 200 Gt). Using a variance-based sensitivity analysis, we show that uncertainties in the friction law and the ice-flow law have a greater impact on the model results than surface mass balance and initial surface elevation.

Graphics-processing-unit-accelerated ice flow solver for unstructured meshes using the Shallow-Shelf Approximation (FastIceFlo v1.0.1)

Abstract. Ice-sheet flow models capable of accurately projecting their future mass balance constitute tools to improve flood risk assessment and assist sea-level rise mitigation associated with enhanced ice discharge. Some processes that need to be captured, such as grounding-line migration, require high spatial resolution (under the kilometer scale). Conventional ice flow models mainly execute on central processing units (CPUs), which feature limited parallel processing capabilities and peak memory bandwidth. This may hinder model scalability and result in long run times, requiring significant computational resources. As an alternative, graphics processing units (GPUs) are ideally suited for high spatial resolution, as the calculations can be performed concurrently by thousands of threads, processing most of the computational domain simultaneously. In this study, we combine a GPU-based approach with the pseudo-transient (PT) method, an accelerated iterative and matrix-free solution strategy, and investigate its performance for finite elements and unstructured meshes with application to two-dimensional (2-D) models of real glaciers at a regional scale. For both the Jakobshavn and Pine Island glacier models, the number of nonlinear PT iterations required to converge a given number of vertices (N) scales in the order of 𝒪(N1.2) or better. We further compare the performance of the PT CUDA C implementation with a standard finite-element CPU-based implementation using the price-to-performance metric. The price of a single Tesla V100 GPU is 1.5 times that of two Intel Xeon Gold 6140 CPUs. We expect a minimum speedup of at least 1.5 times to justify the Tesla V100 GPU price to performance. Our developments result in a GPU-based implementation that achieves this goal with a speedup beyond 1.5 times. This study represents a first step toward leveraging GPU processing power, enabling more accurate polar ice discharge predictions. The insights gained will benefit efforts to diminish spatial resolution constraints at higher computing performance. The higher computing performance will allow for ensembles of ice-sheet flow simulations to be run at the continental scale and higher resolution, a previously challenging task. The advances will further enable the quantification of model sensitivity to changes in upcoming climate forcings. These findings will significantly benefit process-oriented sea-level-projection studies over the coming decades.

A Geodetic-Data-Calibrated Ice Flow Model to Simulate Historical and Future Response of Glaciers in Southeastern Tibetan Plateau

Glaciers play a vital role in the Asian mountain water towers and have significant downstream impacts on domestic, agricultural, and industrial water usage. The rate of glacier mass loss in the Southeastern Tibetan Plateau (SETP) is among the highest in Asia and has intensified in recent decades. However, a comprehensive quantification that considers both spatial and temporal aspects of glacier mass loss across the entire SETP is still insufficient. This study aimed to address this gap by utilizing geodetic datasets specific to each glacier by calibrating the Open Global Glacier Model (OGGM) driven by HAR v2 and reconstructing the glacier mass balance of 7756 glaciers in the SETP from 1980 to 2019 while examining their spatial variability. The findings reveal that the average mass balance during this period was −0.50 ± 0.28 m w.e. a−1, with an accelerated loss observed in the 2000s (average: 0.62 ± 0.24 m w.e. a−1). Notably, central glaciers in the SETP exhibited relatively smaller mass loss, indicating a gradient effect of increased loss from the central region toward the eastern and western sides. By the end of this century, the area, length, and volume of glaciers in the entire SETP region are projected to decrease by 83.57 ± 4.91%, 90.25 ± 4.23%, and 88.04 ± 4.52%, respectively. Moreover, the SETP glacier melt runoff is estimated to decrease by 62.63 ± 6.16% toward the end of the century, with the “peak water” point of glacier melt runoff predicted to occur in 2023 under the SSP585 scenario. Sensitivity experiments demonstrated that the SETP glaciers are more than three times more sensitive to temperature changes than to precipitation variations, and the observed decrease in monsoon precipitation indicates the weakening magnitude of the Indian summer monsoon in recent years. The spatially refined and high-temporal-resolution characteristics of glacier mass loss presented in this study contribute to a better understanding of specific glacier changes in the SETP. Additionally, the prediction results provide valuable references for future water resources management and policy formulation in the SETP region.

Evaluating different geothermal heat-flow maps as basal boundary conditions during spin-up of the Greenland ice sheet

Abstract. There is currently poor scientific agreement on whether the ice–bed interface is frozen or thawed beneath approximately one third of the Greenland ice sheet. This disagreement in basal thermal state results, at least partly, from differences in the subglacial geothermal heat-flow basal boundary condition used in different ice-flow models. Here, we employ seven widely used Greenland geothermal heat-flow maps in 10 000-year spin-ups of the Community Ice Sheet Model (CISM). We perform two spin-ups: one nudged toward thickness observations and the other unconstrained. Across the seven heat-flow maps, and regardless of unconstrained or nudged spin-up, the spread in basal ice temperatures exceeds 10 ∘C over large areas of the ice–bed interface. For a given heat-flow map, the thawed-bed ice-sheet area is consistently larger under unconstrained spin-ups than nudged spin-ups. Under the unconstrained spin-up, thawed-bed area ranges from 33.5 % to 60.0 % across the seven heat-flow maps. Perhaps counterintuitively, the highest iceberg calving fluxes are associated with the lowest heat flows (and vice versa) for both unconstrained and nudged spin-ups. These results highlight the direct, and non-trivial, influence of the heat-flow boundary condition on the simulated equilibrium thermal state of the ice sheet. We suggest that future ice-flow model intercomparisons should employ a range of basal heat-flow maps, and limit direct intercomparisons with simulations using a common heat-flow map.