Tidewater glaciers frequently advance and retreat in ways uncoupled from climate forcing. This complicates the task of forecasting the evolution of individual glaciers and the overall Greenland ice sheet, much of which is drained by tidewater glaciers. Past observational research has identified a set of processes collectively known as the tidewater glacier cycle (TGC) to describe tidewater glacier evolution in four stages: the advancing stage, the extended stage, the retreating stage and the retreated stage. Once glacier retreat is initiated, the TGC is thought to depend largely on the glacier's calving rate, which is controlled by fjord geometry. However, there has been little modeling or systematic observational work on the topic. Measuring calving rates directly is challenging and thus we developed an averaged von Mises stress state at the glacier terminus as a calving rate proxy that can be estimated from surface velocities, ice thickness, a terminus position and subglacial topography. We then analyzed 44 tidewater glaciers in Greenland and assessed the current state in the TGC for them. Of the 44 glaciers, we find that fjord geometry is causing instability in ten cases, vs stability in seven, with 11 in rapid retreat and 16 have been historically stable.

Outlet glaciers in Greenland are undergoing retreat and diffusive thinning in response to external forcings, but the rates and magnitudes of these responses differ from glacier to glacier for unclear reasons. We test how changes in ice overburden pressure and basal lubrication affect diffusive thinning rates and their spatial patterns by conducting numerical experiments over various idealized Greenland-like glacier domains. We find that ~10 km frontal retreat over a decade can produce sustained thinning rates as large as 16 m a−1 due to ice overburden pressure changes, at outlet glaciers with high basal drag (>60 kPa) and lateral resistive stress (>70 kPa). Localized basal lubrication perturbations induce upstream thinning and downstream thickening up to 12 m a−1; the duration of the lubrication forcing generally has a greater effect than its intensity on induced thickness changes. Lastly, episodic grounding line retreats over a rough bed produce a stepped time series of thinning broadly consistent with observations of dynamic elevation change on multiple Greenland glaciers. Our findings highlight the critical role of the total grounding zone – not ice front position – through the resistive stress change in relation to total glacier thinning.

Ruth Glacier is situated in the Central Alaska Range, with the Don Sheldon Amphitheater comprising much of its broad accumulation area, directly adjacent to North America's tallest mountain, Denali. From there it funnels through the ‘Great Gorge,’ flanked by steep valley walls reaching over 1500 m. We combine airborne and ground-based radar measurements of ice thickness with satellite-derived surface velocities to constrain ice flux above and below the gorge, and employ a mass conservation approach to estimate the glacier's thickness within the gorge. We measure ice thickness in the amphitheater to reach 950 m, and estimate centerline thickness in the gorge to range from 610 to 960 m. Our estimates are up to two times greater than those suggested by global models, and allow us to confirm that the Great Gorge rivals Hells Canyon as the deepest gorge in North America. We found that the geometry of the gorge prevents radar measurements of ice thickness there since returns from the subglacial valley walls would precede and potentially occlude nadir bed returns. The same may be true of other unmapped mountain glaciers; however, thickness may be determined using appropriately located flux gates where radar sounding is feasible, combined with mass conservation methods.

The time-evolution of glacier basal motion remains poorly constrained, despite its importance in understanding the response of glaciers to climate warming. Athabasca Glacier provides an ideal site for observing changes in basal motion over long timescales. Studies from the 1960s provide an in situ baseline dataset constraining ice deformation and basal motion. We use two complementary numerical flow models to investigate changes along a well-studied transverse profile and throughout a larger study area. A cross-sectional flow model allows us to calculate transverse englacial velocity fields to simulate modern and historical conditions. We subsequently use a 3-D numerical ice flow model, Icepack, to estimate changes in basal friction by inverting known surface velocities. Our results reproduce observed velocities well using standard values for flow parameters. They show that basal motion declined significantly (30–40%) and this constitutes the majority (50–80%) of the observed decrease in surface velocities. At the same time, basal resistive stress has remained nearly constant and now balances a much larger fraction of the driving stress. The decline in basal motion over multiple decades of climate warming could serve as a stabilizing feedback mechanism, slowing ice transport to lower elevations, and therefore moderating future mass loss rates.

Glaciers of Baffin Island and nearby islands of Arctic Canada have experienced rapid mass losses over recent decades. However, projections of loss rates into the 21st century have so far been limited by the availability of model calibration and validation data. In this study, we model the surface mass balance of the largest ice cap on Baffin Island, Penny Ice Cap, since 1959, using an enhanced temperature index model calibrated with in situ data from 2006–2014. Subsequently, we project changes to 2099 based on the RCP4.5 climate scenario. Since the mid-1990s, the surface mass balance over Penny Ice Cap has become increasingly negative, particularly after 2005. Using volume–area scaling to account for glacier retreat, peak net mass loss is projected to occur between ~2040 and 2080, and the ice cap is expected to lose 22% (377.4 Gt or 60 m w.e.) of its 2014 ice mass by 2099, contributing 1.0 mm to sea level rise. Our 2015–2099 projections are approximately nine times more sensitive to changes in temperature than precipitation, with an absolute cumulative difference of 566 Gt (90 m w.e.) between +2 and −2°C scenarios, and 63 Gt (10 m w.e.) between +20% and −20% precipitation scenarios.

Rifts are full-thickness fractures that propagate laterally across an ice shelf. They cause ice-shelf weakening and calving of tabular icebergs, and control the initial size of calved icebergs. Here, we present a joint inverse and forward computational modeling framework to capture rifting by combining the vertically integrated momentum balance and anisotropic continuum damage mechanics formulations. We incorporate rift–flank boundary processes to investigate how the rift path is influenced by the pressure on rift–flank walls from seawater, contact between flanks, and ice mélange that may also transmit stress between flanks. To illustrate the viability of the framework, we simulate the final 2 years of rift propagation associated with the calving of tabular iceberg A68 in 2017. We find that the rift path can change with varying ice mélange conditions and the extent of contact between rift flanks. Combinations of parameters associated with slower rift widening rates yield simulated rift paths that best match observations. Our modeling framework lays the foundation for robust simulation of rifting and tabular calving processes, which can enable future studies on ice-sheet–climate interactions, and the effects of ice-shelf buttressing on land ice flow.

We investigate unusual discontinuous glacier motion on Thompson Glacier, Umingmat Nunaat, Arctic Canada, using synthetic aperture radar (SAR) images and ice-flow modeling. A novel intensity-rescaling scheme is developed to reduce errors in high-resolution speckle tracking, resulting in a ~25% improvement in accuracy. Interferometric SAR (InSAR) and speckle tracking using high resolution RADARSAT-2 data indicate velocity discontinuities of up to 1 cm d−1 across deep and longitudinally extensive supraglacial channels on Thompson Glacier. We use a cross-sectional finite-element ice-flow model to determine the conditions under which velocity discontinuities of the observed magnitude and signature are possible. The modeling suggests that discontinuous motion across (long and straight) supraglacial channels can occur without ice fracture and under a wide variety of glacier thermal structures, including in fully temperate glaciers. Despite the wide range of conditions conducive to discontinuous motion, the form we observe requires that the associated channels be deep, longitudinally extensive and located in regions of lateral shearing. We speculate that these combined conditions are rare except on polythermal glaciers, where drainage features such as moulins are comparatively scarce and lower deformation rates allow channels to incise consistently and persist over many years.

Water pressure beneath glaciers influences ice velocity. Subglacial hydrology models are helpful for gaining insight into basal conditions, but models depend on unconstrained parameters, and a current challenge is reproducing elevated water pressures in winter. We eliminate terms related to englacial storage, opening by sliding, and melt due to changes in the pressure-melting-point temperature, to create a minimalist version of the Subglacial Hydrology And Kinetic, Transient Interactions (SHAKTI) model, and apply this model to Helheim Glacier in east Greenland to explore the winter base state of the subglacial drainage system. Our results suggest that meltwater produced at the bed alone supports active winter drainage with large areas of elevated water pressure and preferential drainage pathways, using a continuum approach that allows for transitions between flow regimes. Transmissivity varies spatially over several orders of magnitude from 10−4 to 103 m2s−1, with regions of weak transmissivity representing poorly connected regions of the system. Bed topography controls the location of primary drainage pathways, and high basal melt rates occur along the steep valley walls. Frictional heat from sliding is a dominant source of basal melt; different approaches for calculating basal shear stress produce significantly different basal melt rates and subglacial discharge.

Marine-terminating glaciers cover more than one-fourth of the total glacierized area in the Northern Hemisphere outside the Greenland ice sheet. It is therefore crucial to ensure an adequate representation of these glaciers when projecting large-scale glacier mass changes. We investigate how the introduction of marine frontal processes in the modeling chain influences the results of mass change projections, compared to projections neglecting such processes. We find that including frontal processes reduces the projected glacier mass loss, since incorporating frontal ablation in the model's mass-balance calibration results in a decrease in marine-terminating glaciers’ sensitivity to atmospheric temperatures. We also find that retrograde bed slopes lead to increased frontal ablation as the atmosphere warms, while frontal ablation decreases if bed slopes are prograde. These opposing effects have the potential to partly cancel each other when considering large glacier ensembles. Although we do not account for potential future changes in oceanic climate yet, any effect of these would be moderated by around half of today's marine-terminating glaciers becoming land-terminating in the course of the 21st century. While we find a significant influence of ice flow parameters on our results, boundary conditions remain the largest source of uncertainty in our projections.

Floating ice shelves that fringe the coast of Antarctica resist the flow of grounded ice into the ocean. One of the key factors governing the amount of flow resistance an ice shelf provides is the rigidity (related to viscosity) of the ice that constitutes it. Ice rigidity is highly heterogeneous and must be calibrated from spatially continuous surface observations assimilated into an ice-flow model. Realistic uncertainties in calibrated rigidity values are needed to quantify uncertainties in ice sheet and sea-level forecasts. Here, we present a physics-informed machine learning framework for inferring the full probability distribution of rigidity values for a given ice shelf, conditioned on ice surface velocity and thickness fields derived from remote-sensing data. We employ variational inference to jointly train neural networks and a variational Gaussian Process to reconstruct surface observations, rigidity values and uncertainties. Applying the framework to synthetic and large ice shelves in Antarctica demonstrates that rigidity is well-constrained where ice deformation is measurable within the noise level of the observations. Further reduction in uncertainties can be achieved by complementing variational inference with conventional inversion methods. Our results demonstrate a path forward for continuously updated calibrations of ice flow parameters from remote-sensing observations.

We revise and evaluate frontal ablation fluxes obtained by the Open Global Glacier Model (OGGM) for Greenland's tidewater peripheral glaciers de-coupled from the ice sheet. By making use of new region-wide ice thickness and solid ice discharge data, we re-evaluate model performance and suggest future research directions to improve the ice thickness estimation of glacier models. OGGM is unable to predict individual tidewater glacier dynamics well if it has to rely only on surface mass balance estimates and the assumption of a closed budget to constrain the calving parameterization. Velocity observations are essential to constrain the model and estimate the dynamic mass loss of Greenland's tidewater peripheral glaciers.

The Antarctic Peninsula's widespread glacier retreat and ice shelf collapse have been attributed to atmospheric and oceanic warming. Following the initial post-collapse period of retreat, several former tributary glaciers of the Larsen A and B ice shelves have been slowly re-advancing for more than a decade. Here, we use a flowline model of Crane Glacier to gauge the sensitivity of former tributary glaciers to future climate change following this period of long-term dynamic adjustment. The glacier's long-term geometry and speed changes are similar to those of other former Larsen A and B tributaries, suggesting that Crane Glacier is a reasonable representation of regional dynamics. For the unperturbed climate simulations, discharge remains nearly unchanged in 2018–2100, indicating that dynamic readjustment to shelf collapse took ~15 years. Despite large uncertainties in Crane Glacier's past and future climate forcing, a wide range of future climate scenarios leads to a relatively modest range in grounding line discharge (0.8–1.5 Gt a−1) by 2100. Based on the model results for Crane, we infer that although former ice shelf tributaries may readvance following collapse, similar to the tidewater glacier cycle, their dynamic response to future climate perturbations should be less than their response to ice shelf collapse.

Numerical glacier and ice-sheet models compute evolving ice geometry and velocity fields using various stress-balance approximations and boundary conditions. At high spatial resolution, with horizontal mesh/grid resolutions of a few kilometers or smaller, these models usually require time steps shorter than climate-coupling time scales because they update ice thickness after each velocity solution. High-resolution performance is degraded by the stability restrictions of such explicit time-stepping. This short note, which considers the shallow ice approximation and Stokes models as stress-balance end members, clarifies the scaling of numerical model performance by quantifying simulation cost per model year in terms of mesh resolution and the number of degrees of freedom. The performance of current-generation explicit time-stepping models is assessed, and then compared to the prospective performance of implicit schemes. The main results highlight the key roles played by the algorithmic scaling of stress-balance solvers and coupled, implicit-step solvers.

In this study, we investigate the differences between two transient, three-dimensional, thermomechanically coupled ice-sheet models, namely, a first-order approximation model (FOM) and a ‘full’ Stokes ice-sheet model (FSM) under the same numerical framework. For all numerical experiments, we take the FSM outputs as the reference values and calculate the mean relative errors in the velocity and temperature fields for the FOM over 100 years. Four different boundary conditions (ice slope, geothermal heat flux, basal topography and basal sliding) are tested, and by changing these parameters, we verify the thermomechanical behavior of the FOM and discover that the velocity and temperature biases of the FOM generally increase with increases in the ice slope, geothermal heat flux, undulation amplitude of the ice base, and with the existence of basal sliding. In addition, the model difference between the FOM and FSM may accumulate over time, and the spatial distribution patterns of the relative velocity and temperature errors are in good agreement.

We perform Bayesian inference of the parameters of a time-dependent model of ice flow and calving at Helheim Glacier, East Greenland. We find that, while a time-independent calving parameterization can recover the mean observed terminus position, such a model is unable to recover sub-annual variability, even when forced with seasonally varying climate. To address this, we develop a simple stochastic model relating surface runoff rates and calving threshold. Again inferring model parameters from observations, we find that this parameterization is able to reproduce observations with respect to both mean position and characteristic temporal variability. This result demonstrates the importance of considering potential sub-annual controls on calving rates in numerical models, which may include variable undercutting rates or surface runoff-dependent surface crevasse propagation.

Increased surface melt in the percolation zone of the Greenland ice sheet causes significant changes in the firn structure, directly affecting the amount and timing of meltwater runoff. Here we force an energy-balance model with automatic weather stations data at two sites in the percolation zone of southwest Greenland ($2040$ and 2360 m a.s.l.) between spring $2017$ and fall $2019$. Extensive model validation and sensitivity analysis reveal that the skin layer formulation used to compute the surface temperature by closing the energy balance leads to a consistent overestimation of melt by more than a factor of two or three depending on the site. In contrast, model results match the observations well when the model is forced by observed surface temperatures; however, unexplained residuals in the energy balance occur. The sensible and ground heat flux differ markedly in the two simulations accounting largely for the difference in modeled melt amounts. This indicates that the energy available for melt is highly sensitive to small changes in surface temperature. Thus, regional climate models that also use the skin layer formulation may have a bias in surface temperature and melt energy in the percolation zone of the ice sheet.

Data assimilation in high-order ice flow modeling is a challenging and computationally costly task, yet crucial to find ice thickness and ice flow parameter distributions that are consistent with ice flow mechanics and mass balance while best matching observations. Failing to find these distributions that are required as initial conditions leads to a disequilibrium between mass balance and ice flow, resulting in nonphysical transient effects in the prognostic model. Here we tackle this problem by inverting an emulator of the Stokes ice flow model based on deep learning. By substituting the ice flow equations using a convolutional neural network emulator, we simplify, make more robust and dramatically speed up the solving of the underlying optimization problem thanks to automatic differentiation, stochastic gradient methods and implementation of graphics processing unit (GPU). We demonstrate this process by simultaneously inferring the ice thickness distribution, ice flow parametrization and ice surface of ten of the largest glaciers in Switzerland. As a result, we obtain a high degree of assimilation while guaranteeing an equilibrium between mass-balance and ice flow mechanics. The code runs very efficiently (optimizing one large-size glacier at 100 m takes < 1 min on a laptop) while it is open-source and publicly available.

Tidewater glaciers have been observed to experience instantaneous, stepwise increases in velocity during iceberg-calving events due to a loss of resistive stresses. These changes in stress can potentially impact tidewater glacier stability by promoting additional calving and affecting the viscous delivery of ice to the terminus. Using flow models and perturbation theory, we demonstrate that calving events and subsequent terminus readvance produce quasi-periodic, sawtooth oscillations in stress that originate at the terminus and propagate upstream. The stress perturbations travel at speeds much greater than the glacier velocities and, for laterally resisted glaciers, rapidly decay within a few ice thickness of the terminus. Consequently, because terminus fluctuations due to individual calving events tend to be much higher frequency than climate variations, individual calving events have little direct impact on the viscous delivery of ice to the terminus. This suggests that the primary mechanism by which calving events can trigger instability is by causing fluctuations in stress that weaken the ice and lead to additional calving and sustained terminus retreat. Our results further demonstrate a stronger response to calving events in simulations that include the full stress tensor, highlighting the importance of accounting for higher order stresses when developing calving parameterizations.

Upernavik Isstrøm, a marine glacier undergoing rapid retreat, is simulated by forcing a numerical model with ocean-driven melt. A review of processes driving retreat led us to hypothesize that a glacier undergoing rapid retreat may be less sensitive to perturbations in the balance of forces than a glacier that is undergoing moderate changes or a glacier in steady state. Numerical experiments suggest this is not the case, and that a system in rapid retreat is as sensitive to basal traction perturbations as a system that is near to steady state. This result is important when considering other glacier systems experiencing marine-forced retreat. While the ice–ocean interface is of primary importance, additional perturbations from meltwater-forced decoupling of the glacier from its bed continue to feature in glacier dynamics.

Iceberg calving, the process where icebergs detach from glaciers, remains poorly understood. Moreover, few parameterizations of the calving process can easily be integrated into numerical models to accurately capture observations, resulting in large uncertainties in projected sea level rise. Recent efforts have focused on estimating crevasse depths assuming tensile failure occurs when crevasses fully penetrate the glacier thickness. However, these approaches often ignore the role of advecting crevasses on calculations of crevasse depth. Here, we examine a more general crevasse depth calving model that includes crevasse advection. We apply this model to idealized prograde and retrograde bed geometries as well as a prograde geometry with a sill. Neglecting crevasse advection results in steady glacier advance and ice tongue formation for all ice temperatures, sliding law coefficients and constant slope bed geometries considered. In contrast, crevasse advection suppresses ice tongue formation and increases calving rates, leading to glacier retreat. Furthermore, crevasse advection allows a grounded calving front to stabilize on top of sills. These results suggest that crevasse advection can radically alter calving rates and hint that future parameterizations of fracture and failure need to more carefully consider the lifecycle of crevasses and the role this plays in the calving process.