Subglacial topography is an important feature in numerous ice-sheet analyses and can drive the routing of water at the bed. Bed topography is primarily measured with ice-penetrating radar. Significant gaps, however, remain in data coverage that require interpolation. Topographic interpolations are typically made with kriging, as well as with mass conservation, where ice flow dynamics are used to constrain bed geometry. However, these techniques generate bed topography that is unrealistically smooth at small scales, which biases subglacial water flowpath models and makes it difficult to rigorously quantify uncertainty in subglacial drainage patterns. To address this challenge, we adapt a geostatistical simulation method with probabilistic modeling to stochastically simulate bed topography such that the interpolated topography retains the spatial statistics of the ice-penetrating radar data. We use this method to simulate subglacial topography using mass conservation topography as a secondary constraint. We apply a water routing model to each of these realizations. Our results show that many of the flowpaths significantly change with each topographic realization, demonstrating that geostatistical simulation can be useful for assessing confidence in subglacial flowpaths.

Ice shelves restrain flow from the Greenland and Antarctic ice sheets. Climate-ocean warming could force thinning or collapse of floating ice shelves and subsequently accelerate flow, increase ice discharge and raise global mean sea levels. Petermann Glacier (PG), northwest Greenland, recently lost large sections of its ice shelf, but its response to total ice shelf loss in the future remains uncertain. Here, we use the ice flow model Úa to assess the sensitivity of PG to changes in ice shelf extent, and to estimate the resultant loss of grounded ice and contribution to sea level rise. Our results have shown that under several scenarios of ice shelf thinning and retreat, removal of the shelf will not contribute substantially to global mean sea level (<1 mm). We hypothesize that grounded ice loss was limited by the stabilization of the grounding line at a topographic high ~12 km inland of its current grounding line position. Further inland, the likelihood of a narrow fjord that slopes seawards suggests that PG is likely to remain insensitive to terminus changes in the near future.

Gaining knowledge of tidewater glacier (TWG) margin evolution, solid ice flux and their responses to climate over large spatio-temporal scales provides valuable context for the projection of future Greenland ice sheet (GrIS) change. Although studies of sector-wide responses of TWGs exist, studies at an ice-sheet-wide scale have only just become feasible. Here, we present a dataset of 224 annual TWG margins for 1984–2017 (n = 3801), showing that averaged over regional scales, normalised terminus change is linear. Regionally linear retreat trends were identified across most sectors of the GrIS starting in the mid-1990s, although in contrast to previous studies, the northeastern sector is shown to have experienced sustained retreat since the mid-1980s. Through cointegration analyses, individual glaciers are shown to have differing sensitivities to potential climate drivers, though on a sector-wide scale the northwest and southeast are shown to be especially sensitive to annual sea surface temperature and June–July–August air temperature, respectively. Although 92% of the analysed glaciers experience retreat across the GrIS, observed increases in absolute flux for the entire ice sheet can be explained by changes in just 11 of these TWGs.

The recently-advancing Taku Glacier is excavating subglacial sediments at high rates over multi-decadal timescales. However, sediment redistribution over shorter timescales remains unquantified. We use a variety of methods to study subglacial and proglacial sediment redistribution on decadal, seasonal, and daily timescales to gain insight into sub- and proglacial landscape formation. Both excavation and deposition were observed from 2003 to 2015 (2.8 ± 0.9 m a−1 to +2.9 ± 0.9 m a−1). The observed patterns imply that a subglacial conduit has occupied the same site over the past decade. Outwash fans on the subaerial end moraine experience fluvial sediment reworking almost year-round, with net sediment gain in winter and net sediment loss in summer, and an overall mass gain between 2005 and 2015. We estimate that tens of meters of sediment still underlie the glacier terminus, sediments which can be remobilized during future activity. However, imminent retreat from the proglacial moraine will limit its sediment supply, leaving the moraine vulnerable to erosion by bordering rivers. Retreat into an over-deepened basin will leave the glacier vulnerable to increased frontal ablation and accelerating retreat.

We present continuous estimates of snow and firn density, layer depth and accumulation from a multi-channel, multi-offset, ground-penetrating radar traverse. Our method uses the electromagnetic velocity, estimated from waveform travel-times measured at common-midpoints between sources and receivers. Previously, common-midpoint radar experiments on ice sheets have been limited to point observations. We completed radar velocity analysis in the upper ~2 m to estimate the surface and average snow density of the Greenland Ice Sheet. We parameterized the Herron and Langway (1980) firn density and age model using the radar-derived snow density, radar-derived surface mass balance (2015–2017) and reanalysis-derived temperature data. We applied structure-oriented filtering to the radar image along constant age horizons and increased the depth at which horizons could be reliably interpreted. We reconstructed the historical instantaneous surface mass balance, which we averaged into annual and multidecadal products along a 78 km traverse for the period 1984–2017. We found good agreement between our physically constrained parameterization and a firn core collected from the dry snow accumulation zone, and gained insights into the spatial correlation of surface snow density.

As a result of climate warming, glacial meltwater discharge has been increasing in Greenland. During the summers of 2015 and 2016, there were rapid increases in discharge from Qaanaaq Glacier in northwestern Greenland. These discharges resulted in floods that destroyed the road linking the settlement of Qaanaaq to Qaanaaq Airport. Field measurements were performed and a numerical model of glacier runoff was developed to quantify these discharges. The high discharge associated with the 2015 flood, estimated at 9.1 m3 s−1 (hourly mean), resulted from intensive glacier melting due to warm air temperature and strong winds, while the high discharge associated with the 2016 flood resulted from heavy rainfall (90 mm d−1) that led to a peak discharge estimated at 19.9 m3 s−1. The developed model, when used to investigate future glacier runoff under warming conditions, revealed a nonlinear increase in glacial melt with increasing temperature. Additionally, the model forecasted a threefold increase in total summer discharge, owing to a 4 °C rise in temperature. Thus, this study quantified the impact of a changing climate on glacier runoff, which gives insight into future risks of flood hazards along the coast of Greenland.

Ice streams are warmed by shear strain, both vertical shear near the bed and lateral shear at the margins. Warm ice deforms more easily, establishing a positive feedback loop in an ice stream where fast flow leads to warm ice and then to even faster flow. Here, we use radar attenuation measurements to show that the Siple Coast ice streams are colder than previously thought, which we hypothesize is due to along-flow advection of cold ice from upstream. We interpret the attenuation results within the context of previous ice-temperature measurements from nearby sites where hot-water boreholes were drilled. These in-situ temperatures are notably colder than model predictions, both in the ice streams and in an ice-stream shear margin. We then model ice temperature using a 1.5-dimensional numerical model which includes a parameterization for along-flow advection. Compared to analytical solutions, we find depth-averaged temperatures that are colder by 0.7°C in the Bindschadler Ice Stream, 2.7°C in the Kamb Ice Stream and 6.2–8.2°C in the Dragon Shear Margin of Whillans Ice Stream, closer to the borehole measurements at all locations. Modelled cooling corresponds to shear-margin thermal strengthening by 3–3.5 times compared to the warm-ice case, which must be compensated by some other weakening mechanism such as material damage or ice-crystal fabric anisotropy.

During the past 20 years, multi-channel radar emerged as a key tool for deciphering an ice sheet's internal architecture. To assign ages to radar reflections and connect them over large areas in the ice sheet, the layer genesis has to be understood on a microphysical scale. Synthetic radar trace modelling based on the dielectric profile of ice cores allows for the assignation of observed physical properties’ variations on the decimetre scale to radar reflectors extending from the coring site to a regional or even whole-ice-sheet scale. In this paper we rely on the available dielectric profiling data of the northern Greenland deep ice cores: NGRIP, NEEM and EGRIP. The three records are well suited for assigning an age model to the stratigraphic radar-mapped layers, and linking up the reflector properties to observations in the cores. Our modelling results show that the internal reflections are mainly due to conductivity changes. Furthermore, we deduce fabric characteristics at the EGRIP drill site from two-way-travel-time differences of along and across-flow polarized radarwave reflections of selected horizons (below 980 m). These indicate in deeper parts of the ice column an across-flow concentrated c-axis fabric.

Turbulent fluxes make a substantial and growing contribution to the energy balance of ice surfaces globally, but are poorly constrained owing to challenges in estimating the aerodynamic roughness length (z0). Here, we used structure from motion (SfM) photogrammetry and terrestrial laser scanning (TLS) surveys to make plot-scale 2-D and 3-D microtopographic estimations of z0 and upscale these to map z0 across an ablating mountain glacier. At plot scales, we found spatial variability in z0 estimates of over two orders of magnitude with unpredictable z0 trajectories, even when classified into ice surface types. TLS-derived surface roughness exhibited strong relationships with plot-scale SfM z0 estimates. At the glacier scale, a consistent increase in z0 of ~0.1 mm d−1 was observed. Space-for-time substitution based on time since surface ice was exposed by snow melt confirmed this gradual increase in z0 over 60 d. These measurements permit us to propose a scale-dependent temporal z0 evolution model where unpredictable variability at the plot scale gives way to more predictable changes of z0 at the glacier scale. This model provides a critical step towards deriving spatially and temporally distributed representations of z0 that are currently lacking in the parameterisation of distributed glacier surface energy balance models.

Supraglacial debris cover regulates the melt rates of many glaciers in mountainous regions around the world, thereby modifying the availability and quality of downstream water resources. However, the influence of supraglacial debris is often poorly represented within glaciological models, due to the absence of a technique to provide high-precision, spatially continuous measurements of debris thickness. Here, we use high-resolution UAV-derived thermal imagery, in conjunction with local meteorological data, visible UAV imagery and vertically profiled debris temperature time series, to model the spatially distributed debris thickness across a portion of Llaca Glacier in the Cordillera Blanca of Peru. Based on our results, we simulate daily sub-debris melt rates over a 3-month period during 2019. We demonstrate that, by effectively calibrating the radiometric thermal imagery and accounting for temporal and spatial variations in meteorological variables during UAV surveys, thermal UAV data can be used to more precisely represent the highly heterogeneous patterns of debris thickness and sub-debris melt on debris-covered glaciers. Additionally, our results indicate a mean sub-debris melt rate nearly three times greater than the mean melt rate simulated from satellite-derived debris thicknesses, emphasising the importance of acquiring further high-precision debris thickness data for the purposes of investigating glacier-scale melt processes, calibrating regional melt models and improving the accuracy of runoff predictions.

Several global datasets of glacier thickness exist, but the number of observations from western Canada are sparse and spatially biased. To supplement these limited observations, we measured ice thickness with ice penetrating radar on five glaciers in the Columbia Mountains, Canada. Our radar surveys, when combined with previous surveys for two glaciers in the Rocky Mountains, total 182 km of transects that represent 34 672 point measurements of ice thickness. Our measurements are, on average, 38% thicker than previous surface inversion model estimates of glacier thickness. Using our measurements within a cross-validation scheme, we model ice thickness with the Open Global Glacier Model (OGGM) driven with recent observations of surface mass balance and glacier elevation. We calibrated OGGM ice thickness by optimizing the ice creep parameter in the model. The optimized OGGM yields an ice volume for Columbia Basin of 122.5 ± 22.4 km3 for the year 2000, which is 23% greater than the range of previous estimates. At current rates of glacier mass loss for this region, glaciers would disappear from the basin in about 65–80 years. Disappearance of these glaciers will negatively affect the basin's surface hydrology, freshwater availability and aquatic ecosystems.

Assessing the avalanche danger level requires snow stratigraphy and instability data. As such data are usually sparse, we investigated whether distributed snow cover modeling can be used to provide information on spatial instability patterns relevant for regional avalanche forecasting. Using Alpine3D, we performed spatially distributed simulations to evaluate snow instability for the winter season 2016–17 in the region of Davos, Switzerland. Meteorological data from automatic weather stations were interpolated to 100 m horizontal resolution and precipitation was scaled with snow depth measurements from airborne laser scanning. Modeled snow instability metrics assessed for two different weak layers suggested that the weak layer closer to the snow surface was more variable. Initially, it was less stable than the weak layer closer to the ground, yet it stabilized faster as the winter progressed. In spring, the simulated snowpack on south-facing slopes stabilized faster than on north-facing slopes, in line with the regional avalanche forecast. In the winter months January to March 2017, simulated instability metrics did not suggest that the snowpack on south-facing slopes was more stable, as reported in the regional avalanche forecast. Although a validation with field data is lacking, these model results still show the potential and challenges of distributed modeling for supporting operational avalanche forecasting.

We calibrate the calving parameterisation implemented in the Open Global Glacier Model via two methods (velocity constraint and surface mass balance (SMB) constraint) and assess the impact of accounting for frontal ablation on the ice volume estimate of Greenland tidewater peripheral glaciers (PGs). We estimate an average regional frontal ablation flux of 7.38±3.45 Gta−1 after calibrating the model with two different satellite velocity products, and of 0.69±0.49 Gta−1 if the model is constrained using frontal ablation fluxes derived from independent modelled SMB averaged over an equilibrium reference period (1961–90). This second method makes the assumption that most PGs during that time have an equilibrium between mass gain via SMB and mass loss via frontal ablation. This assumption serves as a basis to assess the order of magnitude of dynamic mass loss of glaciers when compared to the SMB imbalance. The differences between results from both methods indicate how strong the dynamic imbalance might have been for PGs during that reference period. Including frontal ablation increases the estimated regional ice volume of PGs, from 14.47 to 14.64±0.12 mm sea level equivalent when using the SMB method and to 15.84±0.32 mm sea level equivalent when using the velocity method.