We present a theoretical framework that integrates the dynamics of glaciers with and without the topographic confinement. This Part 1 paper concerns the former, which may exhibit surge cycles when subjected to thermal switches associated with the bed condition. With the topographic trough setting the glacier width and curbing the lateral drainage of the meltwater, the problem falls under the purview of the undrained plastic bed (UPB) formalism. Employing the UPB, we shall examine the external controls of the glacial behavior and test them against observations. Through our non-dimensionalization scheme, we construct a 2-D regime diagram, which allows a ready prognosis of the glacial properties over the full range of the external conditions, both climate- and size-related. We first discern the boundaries separating the glacial regimes of steady-creep, cyclic-surging and steady-sliding. We then apply the regime diagram to observed glaciers for quantitative comparisons. These include the Svalbard glaciers of both normal and surge types, Northeast Greenland Ice Stream characterized by steady-sliding, and Hudson Strait Ice Stream exhibiting cyclic surges. The quantitative validation of our model containing no free parameters suggests that the thermal switch may unify the dynamics of these diverse glaciers.

This study involves examination of glaciological mass-balance time series, glacier and climatic descriptors, the application of machine learning methods for glaciological clustering, and computation of mass-balance time series based upon the clustering and statistical analyses relative to gridded air temperature datasets. Our analysis revealed an increasingly coherent mass-balance trend but a latitudinal bias of monitoring programs. The glacier classification scheme delivered three clusters, suggesting these correspond to climate-based first-order regimes, as glacier morphometric characteristics weighed little in our multivariate analysis. We combined all available surface mass-balance data from in situ monitoring programs to study temperature sensitivity for each cluster. These aggregated mass-balance time series delivered spatially different statistical relationships to temperature. Results also showed that surface mass balance tends to have a temporal self-correlation of ~20 years. Using this temporal window to analyze sensitivity since ~ 1950, we found that in all cases temperature sensitivity, while generally negative, tended to fluctuate through time, with the largest absolute magnitudes occurring in the 1980s and becoming less negative in recent years, revealing that glacier sensitivity is non-stationary. These findings point to a scenario of a coherent signal of change no matter the glacier regime. This work provides new insights into glacier–climate relationships that can guide observational and analytical strategies.

To understand the mechanism of simultaneous drainage event related to supraglacial lakes on a debris-covered glacier, we investigated water-level variations of supraglacial lakes on the southern Inylchek Glacier in Kyrgyzstan. To examine these variations, we used daily aerial images for 2017–2019 from an uncrewed aerial vehicle that were converted to 15 cm-digital surface models and ortho-images. Our main results are as follows: (1) When one lake drained, the water levels of other lakes simultaneously increased, indicating that drainage water is shared with several lakes through a main englacial conduit. In one drainage event, a branched off englacial conduit clearly connected to a main englacial conduit. (2) Sometimes several lakes discharged simultaneously, indicating that several lakes had connected to a main englacial conduit that had opened. Such cases can cause larger-scale drainage than that from the opening of a branched off englacial conduit. (3) Simultaneous drainage occurred twice in the same year, each time through a different conduit, indicating that the main englacial conduit can be abandoned and reused. (4) In some lakes, the water level on the hydraulic gradient line increased gradually with nearly the same increase rate just before drainage. Such an increase may be an indicator of a possible simultaneous drainage event.

This paper presents a physics-based macroscale model for the densification of dry snow which provides for a smooth transition between densification by grain-boundary sliding (stage 1) and densification by power-law creep (stage 2). The model uses established values of the stage 1 and 2 densification rates away from the transition zone and two transition parameters with a simple physical basis: the transition density and the half-width of the transition zone. It has been calibrated using density profiles from the SUMup database and physically based expressions for the transition parameters have been derived. The transition model produces better predictions of the depth of the nominal bubble close-off horizon than the Herron and Langway model, both in its classical form and in a recent version with re-optimised densification rates.

The response of ice to applied stress on ice-sheet flow timescales is commonly described by a non-linear incompressible viscous fluid, for which the deviatoric stress has a quadratic relation in the strain rate with two response coefficient functions depending on two principal strain-rate invariants I2 and I3. Commonly, a coaxial (linear) relation between the deviatoric stress and strain rate, with dependence on one strain-rate invariant I2 in a stress formulation, equivalently dependence on one deviatoric stress invariant in a strain-rate formulation, is adopted. Glen's uni-axial stress experiments determined such a coaxial law for a strain-rate formulation. The criterion for both uni-axial and shear data to determine the same relation is determined. Here, we apply Steinemann's uni-axial stress and torsion data to determine the two stress response coefficients in a quadratic relation with dependence on a single invariant I2. There is a non-negligible quadratic term for some ranges of I2; that is, a coaxial relation with dependence on one invariant is not valid. The data does not, however, rule out a coaxial relation with dependence on two invariants.

Sea-ice pore microstructure constrains ice transport properties, affecting fluid flow relevant to oil-in-ice transport and biogeochemical processes. Motivated by a lack of pore microstructural data, in particular for granular ice and across the seasonal cycle, throat size, tortuosity, connectivity, and other microstructural variables were derived from X-ray computed tomography for brine-filled pores in seasonal landfast ice off northern Alaska. Data were obtained for granular and columnar ice during the ice growth, transition, and melt season. While granular ice exhibits a more heterogeneous pore space than columnar ice, pore and throat size distributions are comparable. The greater tortuosity of pores in granular (1.2 < τg < 1.7) compared to columnar ice (1.0 < τc < 1.1) compounded with a less interconnected pore space translates into lower permeability for granular ice during the growth season for a given porosity. The microstructural data explain findings of granular ice hindering vertical oil-in-ice transport during ice growth and transition stage. With granular ice more frequent in the changing Arctic, data from studies such as this are needed to inform improved modeling of porosity-permeability relationships.

In Part 1, we have considered the dynamics of topographically confined glaciers, which may undergo surge cycles when the bed becomes temperate. In this Part 2, we consider the ice discharge over a flatbed, which would self-organize into alternating stream/ridge pairs of wet/frozen beds. The meltwater drainage, no longer curbed by the bed trough, would counter the conductive cooling to render a minimum bed strength at some intermediate width, toward which the stream would evolve over centennial timescale. At this stationary state, the stream width is roughly twice the geometric mean of the stream height and length, which is commensurate with its observed width. Over a flatbed, streams invariably interact, and we deduce that the neighboring ones would exhibit compensating cycles of maximum velocity and stagnation over the centennial timescale. This deduction is consistent with observed time variation of Ross ice streams B and C (ISB/C), which is thus a manifestation of the natural cycle. Moreover, the model uncovers an overlooked mechanism of the ISC stagnation: as ISB widens following its reactivation, it narrows ISC to augment the loss of the meltwater, leading to its stagnation. This stagnation is preceded by ice thickening hence opposite to the thinning-induced surge termination.

Lake Vostok, East Antarctica, represents an extensive water surface at the base of the ice sheet. Snow, ice and atmospheric pressure loads applied anywhere within the lake area produce a hydrostatic response, involving deformations of the ice surface, ice–water interface and particle horizons. A modelling scheme is developed to derive height changes of these surfaces for a given load pattern. It is applied to a series of load scenarios, and predictions based on load fields derived from a regional climate model are compared to observational datasets. Our results show that surface height changes due to snow-buildup anomalies are damped over the lake area, reducing the spatial standard deviation by one-third. The response to air pressure variations, in turn, adds surface height variability. Atmospheric pressure loads may produce height changes of up to $\pm$4 cm at daily resolution, but decay rapidly with integration time. The hydrostatic load response has no significant impact neither on ICESat laser campaign biases determined over the lake area nor on vertical particle movements derived from GNSS observations.

This study investigates stagnation conditions of the Pensilungpa glacier, western Himalaya. Multiple glacier parameters (length, area, debris extent and thickness, snowline altitude (SLA), velocity, downwasting and ice cliffs) were studied using field measurements (2016–18), high-resolution imagery from GoogleEarth (2013–17) and spaceborne Landsat, ASTER and SRTM data (1993–2017) to comprehend the glacier's current state. Results show a moderate decrease in length (6.62 ± 2.11 m a−1) and area (0.11 ± 0.03% a−1), a marked increase in SLA (~6 m a−1) and debris cover (2.86 ± 0.29% a−1) and a slowdown of ~50% during 1993–2016. Notable thinning of −0.88 ± 0.04 m a−1 was observed between 2000 and 2017 showing a similar trend as field measurements during 2016–17 (−0.88 m) and 2017–18 (−1.54 m). Further, results reveal a stagnation of the lower ablation zone (LAZ). Less mass supply and heterogeneous debris growth (6.67 ± 0.41% a−1) over the previous decade resulted in slowdown, margin insulation and slope-inversion, leading to stagnation. Stagnation of LAZ caused bulging in the dynamic upper ablation zone and favored the development of supraglacial ponds and ice cliffs. Ice cliffs have grown significantly (48% in number; 41% in area during 2013–17) and their back-wasting now dominates the ablation process.

Most glaciers on the Tibetan Plateau have experienced continuous mass losses in response to global warming. However, the seasonal dynamics of glaciers on the southeastern Tibetan Plateau have rarely been reported in terms of glacier surface elevation and velocity. This paper presents a first attempt to explore the seasonal dynamics of the debris-covered Dagongba Glacier within the southeastern Tibetan Plateau. We use the multitemporal unoccupied aerial vehicle images collected over the lower ablation zone on 8 June and 17 October 2018, and 13 May 2019, and then perform an analysis concerning climatic fluctuations. The results reveal that the mean surface elevation decrease of the Dagongba Glacier during the warm season ($2.81\pm 0.44$ m) was remarkably higher than the cold season ($0.72\pm 0.45$ m). Particularly notable glacier surface elevation changes were found around supraglacial lakes and ice cliffs where ice ablation rates were $\sim$3 times higher than the average. In addition, a larger longitudinal decline of glacier surface velocity was observed in the warm season than that in the cold season. In terms of further comparative analysis, the Dagongba Glacier experienced a decrease in surface velocity between 1982–83 and 2018–19, with a decrease in the warm season possibly twice as large as that in the cold season.

Refreezing of meltwater in firn is a major component of Greenland ice-sheet's mass budget, but in situ observations are rare. Here, we compare the firn density and total ice layer thickness in the upper 15 m of 19 new and 27 previously published firn cores drilled at 15 locations in southwest Greenland (1850–2360 m a.s.l.) between 1989 and 2019. At all sites, ice layer thickness covaries with density over time and space. At the two sites with the earliest observations (1989 and 1998), bulk density increased by 15–18%, in the top 15 m over 28 and 21 years, respectively. However, following the extreme melt in 2012, elevation-detrended density using 30 cores from all sites decreased by 15 kg m−3 a−1 in the top 3.75 m between 2013 and 2019. In contrast, the lowest elevation site's density shows no trend. Thus, temporary build-up in firn pore space and meltwater infiltration capacity is possible despite the long-term increase in Greenland ice-sheet melting.

Khovsgol Lake, a seasonal ice-covered lake located in boreal Siberia, plays important economic, transportation and agricultural roles for local residents. On 11 January 2021, an Mw 6.7 earthquake struck the center of this lake, providing a unique opportunity to better understand lake ice movement, ridge distribution, and hydraulic and hydrodynamic processes, as a result of the mainshock. We use a pixel-tracking method on both Sentinel-1 synthetic aperture radar and Sentinel-2 optical remote-sensing data covering this lake before, during and after the mainshock to obtain lake-ice displacement. Combining these measurements with different viewing geometries, we determine the locations of pressure ridges in each period and derive their relative 3D displacement maps with a precision of <2.2 m. Our results indicate a large E-W displacement of up to 10.3 m and a vertical displacement of up to 4.6 m during the coseismic period, which has a distinctly different magnitude and pattern of displacement compared to the pre- and post-event periods. Given the ice motion pattern and magnitude, we speculate that the unique ice displacement of Khovsgol Lake may have been caused by a rare seismic seiche after the mainshock.

A 3-year record of weekly snow water equivalent (SWE) accumulation at Summit Camp, central Greenland ice sheet, obtained by direct sampling, is presented. While the overall SWE accumulation of 24.2 cm w.e. a−1 matches long-term ice core estimates, variability increases at shorter timescales. Half of the annual SWE accumulation occurs during a few large events, with the average accumulation rate decreasing 35% between the first and second halves of the record coinciding with exceptional anticyclonic conditions in the spring and summer of 2019. No seasonality in accumulation is detected. Rather, local accumulation rates appear to be significantly impacted by wind redistribution that obscures temporal patterns in snowfall. Surface snow density is consistent, on average, with previously measured values but does not correlate with near surface temperature or wind speed. Two surface mass balance reanalysis models significantly underestimate accumulation rates at Summit Camp. This is concerning because such models are often used to estimate ice-sheet mass loss.

This paper presents a feasibility study of in situ field measurements of unsaturated meltwater percolation flux within the vertical profile of a snowpack, using the self-potential (SP) method. On-site snowmelt column tests calibrated the SP measurements. The SP data measured electrical field strength with an electrode spacing of 20 cm, and coincident water saturation (Sw) measurements using time domain reflectometry allowed calculation of SP-modeled vertical percolation flux (qsp), expressed as Darcy velocity. The results reflected transient diurnal snowmelt dynamics, with peak flux lagging arrival of a saturation wetting front. Peak daily qsp was 60 to >300 mm d−1, whereas daily snowmelt was 20–50 mm w.e. Surface refreezing events appeared to cause upward flow, possibly representing water redistribution toward the freezing boundary. Calculated fluxes were comparable to actual fluxes, although average errors ranged from −15 to +46% compared to average of melt expected from surface energy-balance and ablation stake measurements. By advancing method development to measure unsaturated meltwater percolation flux in snowpacks this study creates opportunities to study fundamental snowmelt processes, may improve mathematical modeling and may supplement glacier mass-balance studies and studies of snowmelt interactions with avalanches, groundwater and surface water.

We investigate the errors caused by neglecting the crystal-orientation fabric when inferring the basal friction coefficient field, and whether such errors can be alleviated by inferring an isotropic enhancement factor field to compensate for missing fabric information. We calculate the steady states that arise from ice flowing over a sticky spot and a bedrock bump using a vertical-slab numerical ice-flow model, consisting of a Weertman sliding law and the anisotropic Johnson flow law, coupled to a spectral fabric model of lattice rotation and dynamic recrystallisation. Given the steady or transient states as input for a canonical adjoint-based inversion, we find that Glen's isotropic flow law cannot necessarily be used to infer the true basal drag or friction coefficient field, which are obscured by the orientation fabric, thus potentially affecting vertically integrated mass fluxes. By inverting for an equivalent isotropic enhancement factor, a more accurate mass flux can be recovered, suggesting that joint inversions for basal friction and the isotropic flow-rate factor may be able to compensate for mechanical anisotropies caused by the fabric. Thus, in addition to other sources of rheological uncertainty, fabric might complicate attempts to relate subglacial conditions to basal properties inferred from an inversion relying on Glen's law.

The once-contiguous Ellesmere Ice Shelf, first reported in writing by European explorers in 1876, and now almost completely disintegrated, has rolling, wave-like surface topography, the origin of which we investigate using a viscous buckling instability analysis. We show that rolls can develop during a winter season (~ 100 d) if sea-ice pressure (depth-integrated horizontal stress applied to the seaward front of the Ellesmere Ice Shelf) is sufficiently large (1 MPa m) and ice thickness sufficiently low (1–10 m). Roll wavelength initially depends only on sea-ice pressure, but evolves over time depending on amplitude growth rate. This implies that a thinner ice shelf, with its faster amplitude growth rate, will have a shorter wavelength compared to a thicker ice shelf when sea-ice pressure is equal. A drawback of the viscous buckling mechanism is that roll amplitude decays once sea-ice pressure is removed. However, non-Newtonian ice rheology, where effective viscosity, and thus roll change rate, depends on total applied stress may constrain roll decay rate to be much slower than growth rate and allow roll persistence from year to year. Whether the viscous-buckling mechanism we explore here ultimately can be confirmed as the origin of the Ellesmere Ice Shelf rolls remains for future research.

The surface restitution method we present reconstructs the evolution of a glacier surface between two time-separated surface topographies using seasonal surface mass balance (SMB) data. A conservative and systematic error analysis is included, based on the availability of surface elevation measurements within the period. The method is applied from 2001 to 2013 at Hurd Glacier (a 4 km2 glacier), where we have sufficient SMB and elevation data. We estimate surface elevation changes in two steps: (1) elevation change due to SMB and (2) elevation change due to glacier dynamics. Four different models of the method are compared depending on whether or not accumulation is memorised at each time step and whether they employ balance profiles or SMB maps. Models are validated by comparing a set of surface measurements retrieved in 2007 with the corresponding restituted elevations. Although surface elevation change between 2001 and 2007 was larger than 10 m, more than 80% of the points restituted by the four models showed errors below ±1 m compared to only 33% when predicted by a linear interpolator. As error estimates between models differ by 0.10 m, we recommend the simplest model, which does not memorise accumulation and interpolates SMB by elevation profiles.

This paper introduces the Instructed Glacier Model (IGM) – a model that simulates ice dynamics, mass balance and its coupling to predict the evolution of glaciers, icefields or ice sheets. The novelty of IGM is that it models the ice flow by a Convolutional Neural Network, which is trained from data generated with hybrid SIA + SSA or Stokes ice flow models. By doing so, the most computationally demanding model component is substituted by a cheap emulator. Once trained with representative data, we demonstrate that IGM permits to model mountain glaciers up to 1000 × faster than Stokes ones on Central Processing Units (CPU) with fidelity levels above 90% in terms of ice flow solutions leading to nearly identical transient thickness evolution. Switching to the GPU often permits additional significant speed-ups, especially when emulating Stokes dynamics or/and modelling at high spatial resolution. IGM is an open-source Python code which deals with two-dimensional (2-D) gridded input and output data. Together with a companion library of trained ice flow emulators, IGM permits user-friendly, highly efficient and mechanically state-of-the-art glacier and icefields simulations.

In polar deserts, depth hoar (hereinafter: DH) growth is not systematic unlike on tundra and this is critical for snowpack properties. Here, we address the spatio-temporal variability of the DH layer in the polar desert at two sites in the Canadian High Arctic: Ward Hunt Island (83° N) and Resolute Bay (75° N). Our data show that, over humid areas, DH represented a larger fraction of the snowpack and was characterized by lower density and coarser crystals than over dry gravelly areas. Increased soil moisture extends the zero-curtain period during freeze-up, leading to stronger temperature gradients in the snowpack and greater kinetic metamorphism. Our results also demonstrate that the large inter-annual variability in DH is primarily driven by wind conditions in the fall since this key variable controls the initial snow density and snow onset date. These strong controls exerted by soil moisture and meteorological conditions on DH growth in polar deserts highlight the possibility of major changes in polar snowpacks physical properties in response to the rapid climate and environmental changes currently affecting these regions.

The englacial and subglacial drainage systems exert key controls on glacier dynamics. However, due to their inaccessibility, they are still only poorly understood and more detailed observations are important, particularly to validate and tune physical models describing their dynamics. By creating artificial glacier moulins – boreholes connected to the subglacial drainage system and supplied with water from surface streams – we present a novel method to monitor the evolution of an englacial hydrological system with high temporal resolution. Here, we use artificial moulins as representations for vertical, pressurised, englacial R-channels. From tracer and pressure measurements, we derive time series of the hydraulic gradient, discharge, flow speed and channel cross-sectional area. Using these, we compute the Darcy–Weisbach friction factor, obtaining values which increase from 0.1 to 13 within five days of channel evolution. Furthermore, we simulate the growth of the cross-sectional area using different temperature gradients. The comparison to our measurements largely supports the common assumption that the temperature follows the pressure melting point. The deviations from this behaviour are analysed using various heat transfer parameterisations to assess their applicability. Finally, we discuss how artificial moulins could be combined with glacier-wide tracer experiments to constrain parameters of subglacial drainage more precisely.