Search

For reconstructing paleoclimate or studying glacial isostatic effects on the Earth’s lithosphere, increasingly more studies focus on modeling the large-scale ice cover in mountainous regions over long time scales. However, balancing model complexity and the spatial extent with computational costs is challenging. Previous studies of large-scale ice cover simulation in mountain areas such as the European Alps, New Zealand, and the Tibetan Plateau, typically used 1-2 km spatial resolution. However, mountains are characterized by high peaks and steep slopes - topographic features that are crucial for glacier mass balance and dynamics, but poorly resolved in coarse resolution topography. The Instructed Glacier Model (IGM) is a novel 3D ice model equipped with a Convolutional Neural Network which is trained from high order ice models to simulate ice flow. This results in a significant acceleration of run times, and thereby opening the possibility of higher spatial resolution runs. We use IGM to model the glaciation of the European Alps with different resolutions (2 km, 1 km and 200 m) over a time period of 160,000 years. We apply a linear cooling rate to today’s climate until 6 °C colder to mimic ice age conditions. Preliminary results indicate systematic, resolution-related differences: At the beginning of cooling the 2 km resolution yields slightly more ice volume. However, this trend reverses after the large valleys are filled with thick ice. When the Alps are fully ice covered, we find up to 15% more ice volume in the higher resolution models., The 28th IUGG General Assembly (IUGG2023) (Berlin 2023)

For more than a century, a major focus of Quaternary research has been the investigation of glacier advances and retreats during ice age cycles. The cradle for these studies are the European Alps. In recent years, glacier modelling permitted to produce long-term transient simulations of the European Alps glacier evolution. However, only sparse empirical data, e.g. geological reconstructions of glacier margins in the foreland of the Alps of the Last Glacial Maximum, are available to validate these simulations.Speleothems from the Alps are a widespread palaeoclimate archive. They provide stable carbon isotope records, which can potentially inform about soil and vegetation conditions above a cave site but also about the lack of soil during times of glacier coverage. In addition, speleothem growth in cold, high-elevation cave sites during glacials are a strong indicator of temperatures in the soil-karst-cave system above the freezing point, which is only likely to occur if the cave is covered by a temperate glacier. We here utilise existing speleothem data (growth histories and stable carbon isotopes) from Alpine caves to infer glacier coverage and thermo-dynamical state during the last glacial cycle and to assess the compatibility with modelled reconstructions. We compare data from the last glacial cycle from multiple cave sites located at different elevations (785 to 2319 m a.s.l.) in the Alps with simulations obtained with the Parallel Ice Sheet Model (PISM). We find a general agreement between speleothem-derived soil presence or absence and modelled glacier coverage. However, the speleothem data provide evidence of temperate glacier coverage which is not shown by all of the PISM simulations. Our work demonstrates the value of speleothem-based reconstructions from the Alps as proxies for assessing performance of palaeo-ice flow models.

Mountain glaciers are a major source of sea-level rise and also represent an important freshwater resource in many mountainous regions. Thus, accurate estimations of their thickness and, therefore, the total ice volume are important both in predicting and mitigating the global and local effects of climate change. However, to date, only 2% of the world’s glaciers outside the ice sheets have any thickness data, due to the logistical difficulties of obtaining such measurements, creating a large and policy-relevant scientific gap.The recent development of a global-scale ice-velocity dataset, however, provides an ideal opportunity to fill this gap and determine ice thickness across the 98% of glaciers for which no thickness data is available. This can be done by inverting an ice-dynamics model to solve for the ice thickness. For accurate thickness results, this needs to be a full-Stokes model, but such a model is far too computationally cumbersome to apply on a global scale, and simpler, quicker methods usually based on the shallow ice approximation (SIA) are too inaccurate, particularly where sliding dominates glacier motion. The only attempt that has been made to leverage the global velocity dataset to retrieve ice thickness has, though, used the SIA, simply because higher-order approaches are not computationally realistic at this scale. Consequently, most of the widely-used global glacier models have made no concerted attempt to invert for global ice thickness, owing to these limitations.As an additional related problem, failing to fully assimilate ice-velocity data into an ice-flow model necessarily introduces a shock when initialising prognostic glacier simulations, resulting in model glaciers and predictions that may diverge substantially from their real-world counterparts.As a solution to these problems, we present results from a deep-learning-driven inversion model that emulates the performance of state-of-the-art full-Stokes models at a thousandth of the computational cost. This model, by solving an optimisation problem, can fully use and assimilate all available input datasets (surface velocity and topography, ice thickness, etc.) as components of its cost function to simultaneously invert for and optimise multiple control parameters (here, we focus on ice thickness). This approach also gives us the possibility of using the same ice-velocity field for inversion and forward modelling, reducing the magnitude of the shock inherent in traditional modelling approaches. With a view to a large-scale application to all the world’s 200,000 glaciers, we present initial thickness-inversion results for the relatively well-documented European Alps to help constrain model parameters and provide a test bed for extension to other glaciated regions, with initial extension to the Caucasus and the Southern Alps.

Outlet glaciers and ice streams transport ice from the ice sheets to the ocean, where the glaciers lose mass by iceberg calving. Sermeq Kujalleq (SKK, Jakobshavn Isbræ) is one of the largest and most dynamic ice streams of the Greenland Ice Sheet with ice flow velocities up to 40 m/day. With extensive fieldwork and detailed repeated UAV surveys we aim at understanding the complex processes occuring at the ice stream margins and at the calving front of SKK. Such processes are often neclected in numerical models inducing uncertainties in projections of the ice sheet evolution.Within the framework of the COEBELI project we conducted drone photogrammetry surveys in July 2022 at SKK along other field measurements including in-situ GPS, GPRI, seismometers, and time-lapse imagery. Despite challenging weather conditions and constraints due to flying restrictions, we acquired 17 repeated flight surveys over the calving front of SKK during two weeks. As a result, we produced a large imagery data set, which was processed to infer high-resolution ortho-images and digital elevation models (DEM). Comparing the different products enables us to estimate changes in surface topography and ice dynamics. During the observation period several large calving events occurred allowing us to investigate the interaction between frontal processes and ice flow dynamics. With the very detailed data we can study crevasse opening, acceleration at the front, weaknesses in the ice and their origin, and the reaction of the glacier to large calving events.

Outlet glaciers and ice streams are the main channels through which ice sheets transport their mass towards the ocean. One of Greenland’s largest outlet glaciers Sermeq Kujalleq in Kangia (Jakobshavn Isbrae) has been broadly researched after experiencing a rapid retreat of the terminus and accompanying speedup to up to 40 m/day in the early 2000’s. However, such short-term ice dynamic variations remain poorly understood making numerical models difficult to constrain and predictions on future sea-level rise uncertain.The short-term ice dynamics of Sermeq Kujalleq consists in transient states and can only be captured by in-situ measurements of high spatial and temporal resolution. Glacier seismology has proven to be a valuable tool to study these dynamics, it provides data with a high temporal resolution and can provide information on processes happening below the ice surface. Within the COEBELI project we combine passive glacier seismology with global navigation satellite system (GNSS) receivers, long-range drones, time-lapse cameras and terrestrial radar interferometry to capture processes such as calving and basal sliding at their respective timescales.Here, we present results from a multi-array seismic deployment at Sermeq Kujalleq in Summer 2022. From May until September two arrays were deployed in the upstream part of the fast-flowing ice stream (>22 km from calving front) and one array on slower moving ice North of the main trunk. For a 3-week period in July, four more arrays were deployed on the fast-flowing ice stream closer to the calving front (50 Hz) icequakes. By studying these different signals, we are able to better constrain the processes and forces that control fluctuating ice-flow velocity and calving events.

Sermeq Kujalleq in Kangia (Jakobshavn Isbræ), Greenland is one of the most studied glaciers in the world mainly due to its recent retreat associated with extremely fast ice stream flow and high solid ice discharge. However, large limitations remain in the understanding of its short-term ice dynamics as the study of sub-daily variations, generally undetectable in spaceborne observations, requires high-rate field measurements that are challenging to acquire. Here, we present glacier surface velocities determined in Post-Processed Kinematic (PPK) mode from eight autonomous Global Navigation Satellite System (GNSS) stations deployed in July 2022 along the ice stream at a distance of 4 to 30 kilometers from the calving front. During this field campaign, we identified an 8-hour-long glacier speedup which was recorded at all GNSS stations and reached up to 11% of the pre-event velocity, followed by a 12-hour-long slowdown of similar magnitude. We further found the peak velocity was first measured at a GNSS station 16 kilometers away from the calving front, then recorded consecutively at each of the three other downstream GNSS stations with a positive time lag corresponding to a ~3 km/h wave propagation speed. At the station closest to the calving front, the timing of peak velocity corresponded to the occurrence of large-scale calving events. We further present line-of-sight glacier surface velocities measured along three shear margin transects with a terrestrial radar interferometer deployed simultaneously with the GNSS array. Across all profiles, we observed a widespread and simultaneous response of fast- and slow-moving ice suggesting a strong coupling between the main trunk and the shear margins of the ice stream.

Convolutional Neural Networks trained from physical high-order ice flow model realizations have proven to be outstanding emulators in terms of computational-to-fidelity performance. However, the dependence on an ensemble of realizations of an instructor model renders this strategy difficult to generalize to the variety of glacier shape and ice flow that can be met in nature. To overcome this issue, we adopt the approach of physics-informed deep learning, which fuses traditional numerical solving by finite differences/elements and deep-learning approaches. Here we train a CNN to minimize the energy associated with high-order ice flow equations either offline over a glacier catalogue, or online directly within the time iterations of a glacier evolution model. As a result, we obtain a generic emulator that features simplified and customizable training, computational efficiency (especially on GPU), and high fidelity to the original model. Embedded into the ``Instructed Glacier Model’’ (IGM) framework, the potential of the emulator is illustrated with two applications including a large-scale high-resolution (3000x4000) forward glacier evolution model run and an inverse modelling exercise for data assimilation. Our emulator is a promising alternative to traditional solvers thanks to its high computational efficiency (especially on GPU), its capability to handle various ice flow and memorize previous solutions, and its relatively simple and open-source implementation.

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.

Future mass loss projections of the Greenland ice sheet require understanding of the processes at a glacier terminus, especially of iceberg calving. We present detailed and high-rate terrestrial radar interferometer observations of Eqip Sermia and Bowdoin Glacier, two outlet glaciers in Greenland with comparable dimensions and investigate iceberg calving, surface elevation, velocity, strain rates and their links to air temperature, tides and topography. The results reveal that the two glaciers exhibit very different flow and calving behaviour on different timescales. Ice flow driven by a steep surface slope with several topographic steps leads to high velocities, areas of extension and intense crevassing, which triggers frequent but small calving events independent of local velocity gradients. In contrast, ice flow under smooth surface slopes leaves the ice relatively intact, such that sporadic large-scale calving events dominate, which initiate in areas with high shearing. Flow acceleration caused by enhanced meltwater input and tidal velocity variations were observed for terminus sections close to floatation. Firmly grounded terminus sections showed no tidal signal and a weak short-term reaction to air temperature. These results demonstrate reaction timescales to external forcings from hours to months, which are, however, strongly dependent on local terminus geometry., Journal of Glaciology, 69 (275), ISSN:0022-1430, ISSN:1727-5652

We present the Instructed Glacier Model (IGM) -- a new framework to model the evolution of glaciers at any scale by coupling ice dynamics, surface mass balance, and mass conservation. The key novelty of IGM is that it models the ice flow by a Convolutional Neural Network (CNN), which is trained from physical high-order ice flow mechanical models. Doing so has major advantages in both forward and inverse modelling. In forward modelling, the most computationally demanding model component (the ice flow) is substituted by a very cheap CNN emulator. Once trained with representative data, IGM permits to model individual mountain glaciers several orders of magnitude faster than high-order ones on CPU with fidelity levels above 90 \% in terms of ice flow solutions leading to nearly identical transient thickness evolution. Switching to Graphics Processing Unit (GPU) permits additional significant speed-ups, especially when modelling large-scale glacier networks and/or high spatial resolutions. In inverse modelling, the substitution by a CNN emulator does not only speed up but facilitates dramatically the data assimilation step, i.e. the search for optimal ice thickness and ice flow parameter spatial distributions that match spatial observations at best (such as ice flow, surface topography or ice thickness profiles) while being consistent with the high-order ice flow mechanics. The reason is that inverting a CNN can take great benefit from the tools used for its training such as automatic differentiation, stochastic gradient methods, and GPU. IGM is an open-source Python code (https://github.com/jouvetg/igm), which deals with two-dimensional gridded input and output data. Together with a companion library of trained ice flow emulators, IGM permits user-friendly, computationally highly-efficient, easy-to-customize, and mechanically state-of-the-art glacier forward and inverse modelling at any scale. We illustrate its potential by replicating a simulation of the great Aletsch Glacier, Switzerland, from 1880 to 2100, based on a Stokes model.

Glacier modelling of the Alpine region during past ice ages has received increasing attention in the recent years. Considering the complexity of the Alpine topography, high spatial resolution climate information is required for running glacier models over the area. However, continuous climatic signal at resolution higher than 10 km and covering several hundred thousand of years cannot be directly derived using dynamical climate models. Alternative strategies must be considered.Here, a climate emulator providing monthly temperature and precipitation over the Alpine region, with a horizontal resolution of 2x2 km and covering the last 400’000 years at intervals of 100 years, is presented. The emulatorcombines a dynamical modelling chain of Earth System Models (ESMs) and a Regional Climate Model (RCM) with different statistical modelling methods. The dynamical modelling chain delivers climate information at highest accuracy based upon physical prognostic equations for specific time slices. The subsequent statistical modelling uses these time slices as physically consistent boundaries and estimates the climate conditions in between them, thus generatinga long-term climate evolution.A total of 19 climate model experiments are conducted for different time-slices of the considered study period, including also sensitivity tests with changes in the ice sheet height of the Northern Hemisphere and land cover type. Starting from a simple linear regression, a series of different statistical approaches is tested for building the emulator. An evaluation of the different versions is then conducted against one of the RCM time-slice experiments, left out in turns from the training set of the statistical model. Results show robust skills of the emulator in the representation of temperature, whose changes are mainly driven by smooth variations in the seasonal pattern of insolation at different time-steps, already using a simple linear regression. For precipitation, non-linearities associated to changes in the large-scale atmospheric circulation seem to dominate, making the use of more complex statistical approaches more appropriate. Additional evaluation tests conducted using glacier modelling driven by the outputs of the developed emulator confirm its potential for reconstructing the ice extent over the Alpine region during the last ice ages.

Sermeq Kujalleq in Kangia (Jakobshavn Isbræ), Greenland has been extensively investigated over the past decades due to its recent retreat associated with extremely fast ice stream flow and high solid ice discharge. However, its short-term dynamics still remain poorly understood as they consist in transient states that can only be captured by high spatial and temporal in situ measurements. In the new COEBELI project, we aim at combining high resolution field data sets from seismic arrays, global navigation satellite system (GNSS) receivers, long-range uncrewed aerial vehicles and terrestrial radar interferometers (TRI) to achieve a comprehensive and detailed study of the short-term ice stream dynamics. Here, we present TRI and GNSS retrievals of surface velocity and elevation acquired during the first, exploratory field campaign of the COEBELI project in summer 2021. Seven kilometers away from the calving front, we specifically identified a slowdown of 1.12 m d-1 within a single day in the main trunk of the ice stream. While the absolute slowdown is larger in the main trunk than in the outer area of the shear margin (1.12 m d-1 versus 0.75 m d-1), it corresponds to a larger fraction of the pre-slowdown velocity in the latter zone (-4.48% versus -7.94%). We further discuss the challenges associated with the acquisition, processing and analysis of high-resolution data sets for the study of such complex and dynamic environments.

We give an overview of the Instructed Glacier Model (IGM) -- a new framework to model the evolution of glaciers at any scale by coupling ice dynamics, surface mass balance, and mass conservation. The key novelty of IGM is that it models the ice flow by a Convolutional Neural Network (CNN), which is trained from physical high-order ice flow mechanical models. Doing so has major advantages in both forward and inverse modelling.In forward modelling, the most computationally demanding model component (the ice flow) is substituted by a very cheap CNN emulator. Once trained with representative data, IGM permits to model individual mountain glaciers several orders of magnitude faster than high-order ones on CPU with fidelity levels above 90 % in terms of ice flow solutions leading to nearly identical transient thickness evolution. Switching to Graphics Processing Unit (GPU) permits additional significant speed-ups, especially when modelling large-scale glacier networks and/or high spatial resolutions.In inverse modelling, the substitution by a CNN emulator does not only speed up but facilitates dramatically the data assimilation step, i.e. the search for optimal ice thickness and ice flow parameter spatial distributions that match spatial observations at best (such as ice flow, surface topography or ice thickness profiles) while being consistent with the high-order ice flow mechanics. The reason is that inverting a CNN can take great benefit from the tools used for its training such as automatic differentiation, stochastic gradient methods, and GPU. IGM is an open-source Python code (https://github.com/jouvetg/igm), which deals with two-dimensional gridded input and output data. Together with a companion library of trained ice flow emulators, IGM permits user-friendly, computationally highly-efficient, easy-to-customize, and mechanically state-of-the-art glacier forward and inverse modelling at any scale. We illustrate its potential by replicating a simulation of the great Aletsch Glacier, Switzerland, from 1880 to 2100, based on a Stokes model. The complete workflow (data assimilation and 220 years long forward modelling) at 100 m of resolution takes about 1-2 min on the GPU of a laptop and can be replicated and adapted easily using an online Colab notebook.

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 nov elty 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 com putationally 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, Journal of Glaciology, 68 (270), ISSN:0022-1430, ISSN:1727-5652

Journal of Glaciology, 65 (253), ISSN:0022-1430, ISSN:1727-5652

Glacier calving plays a key role in the recently observed dynamic massloss of the Greenland ice sheet. Calving waves, generated by thesudden detachment of ice from the glacier terminus, can reach tens ofmeters of height and have devastating effects upon impact onsurrounding shores. In this study, we describe a new method for thedetection of source location and timing of calving waves, and theanalysis of their magnitude and spreading properties using aterrestrial radar interferometer (TRI). This method was applied to11,500 minute-interval TRI acquisitions from Eqip Sermia, Greenland.More than 2,000 calving waves were detected within sevendays. Quantitative assessment with a Wave Power Index (WPI) showedspatially distinctive patterns: the sector of the calving front endingin deep water shows a higher wave activity (+49%) with highercumulative WPI (+34%) than the shallow sector. In combination witha detection of meltwater plume locations, we highlighted a 2.3 timeshigher occurrence of visible meltwater plumes in the deep sector than theshallow one. We found both the cumulated WPI and the number of wavesto increase by more than 80% in the presence of a meltwater plumein the deep sector while only by 30% in the shallow sector. Wetherefore explain the higher calving activity in the deep sector to be strongly related to a combination of higher occurrence of meltwater plumes and more efficient calving enhancement linked to better connections to deep warm waters.

Crystal orientation fabric (COF) analysis provides information about the c-axis orientation of ice grains and the associated anisotropy and microstructural information about deformation and recrystallisation processes within the glacier. This information can be used to introduce modules that fully describe the microstructural anisotropy or at least direction-dependent enhancement factors for glacier modelling. The COF was studied at an ice core that was obtained from the temperate Rhonegletscher, located in the central Swiss Alps. Seven samples, extracted at depths between 2 and 79 m, were analysed with an automatic fabric analyser. The COF analysis revealed conspicuous four-maxima patterns of the c-axis orientations at all depths. Additional data, such as microstructural images, produced during the ice sample preparation process, were considered to interpret these patterns. Furthermore, repeated high-precision global navigation satellite system (GNSS) surveying allowed the local glacier flow direction to be determined. The relative movements of the individual surveying points indicated longitudinal compressive stresses parallel to the glacier flow. Finally, numerical modelling of the ice flow permitted estimation of the local stress distribution. An integrated analysis of all the data sets provided indications and suggestions for the development of the four-maxima patterns. The centroid of the four-maxima patterns of the individual core samples and the coinciding maximum eigenvector approximately align with the compressive stress directions obtained from numerical modelling with an exception for the deepest sample. The clustering of the c axes in four maxima surrounding the predominant compressive stress direction is most likely the result of a fast migration recrystallisation. This interpretation is supported by air bubble analysis of large-area scanning macroscope (LASM) images. Our results indicate that COF studies, which have so far predominantly been performed on cold ice samples from the polar regions, can also provide valuable insights into the stress and strain rate distribution within temperate glaciers., The Cryosphere, 15 (2), ISSN:1994-0416, ISSN:1994-0424

Future mass loss predictions, and thereby sea level rise predictions, are strongly affected by the representation of iceberg calving in numerical ice sheet models. Despite recent advances, gaps in our understanding of calving mechanisms remain and there exists a lack of data to constrain mechanical properties related to ice fracturing. For instance, observed critical strain rates for crevasse initiation span two orders of magnitude.Bowdoin Glacier in Northwest Greenland provides a unique opportunity to conduct in-situ measurements near the calving front due to its accessibility via a crevasse-free walkable moraine. In July 2019, two major calving events were surveyed by 10 GPS stations installed along the front in close vicinity to the calving events. Measurements show glacier uplift prior to the first calving event and horizontal compression prior to the second major calving event.In contrast to previously observed major events, no precursor such as a large surface crack was visible on the field. Our data suggest a change in calving behaviour from surface crevasses due to hydro-fracturing to basal crevasse formation due to buoyancy, which may be favoured by observed thinning (~4 m yr-1 since 2013).

The subglacial hydrology of tidewater glaciers is a key but poorly understood component of the complex ice-ocean system, which affects sea level rise. As it is extremely difficult to access the interior of a glacier, our knowledge relies mostly on the observation of input variables such as air temperature, and output variables such as the ice flow velocities reflecting the englacial water pressure, and the dynamics of plumes reflecting the discharge of meltwater into the ocean. In this study we use a cost-effective Vertical Take-Off and Landing (VTOL) Unmanned Aerial Vehicle (UAV) to monitor the daily movements of Bowdoin Glacier, north-west Greenland, and the dynamics of its main plume. Using Structure-from-Motion photogrammetry and feature-tracking techniques, we obtained 22 high-resolution ortho-images and 19 velocity fields at the calving front for 12 days in July 2016. Our results show a two-day-long speed-up event (up to 170%) – caused by an increase in buoyant subglacial forces – with a strong spatial variability revealing that enhanced acceleration is an indication of shallow bedrock. Further, we used the Particle Image Velocimetry (PIV) method to analyze water flow from successive UAV images taken while flying over the main plume of the glacier. We found that PIV successfully captures the area of radially diverging flow of the plume, and provides information on spatial and time variability as no other remote sensing technique can. Most interestingly, the active part of the plume features pulsating water jets at the time scale of seconds, and is 1 to 5 times smaller than its visual footprint defined by the iceberg-free area. Combined with an ice flow model or a non-steady plume model, our approach has the potential to generate a novel set of input data to gather information about the depth of the bedrock, the discharge of meltwater, or the subglacial melting rate of tidewater glaciers.

Environments along the coast of Greenland are rapidly changing under the influence of a warming climate in the Arctic. To better understand the changes in the coastal environments, we performed researches in the Qaanaaq region in northwestern Greenland as a part of the ArCS (Arctic Challenge for Sustainability) Project. Mass loss of ice caps and marine-terminating outlet glaciers were quantified by field and satellite observations. Measurements and sampling in fjords revealed the important role of glacial meltwater discharge in marine ecosystems. Flooding of a glacial stream in Qaanaaq and landslides in a nearby settlement were investigated to identify the drivers of the incidents. Our study observed rapid changes in the coastal environments, and their critical impact on the society in Qaanaaq. We organized workshops with the residents to absorb local and indigenous knowledge, as well as to share the results and data obtained in the project. Continuous effort towards obtaining long-term observations requiring involvement of local communities is crucial to contribute to a sustainable future in Greenland.

We develop a posteriori ‘mechanical’ error estimators that are able to evaluate the solution discrepancy between two ice flow models. We first reformulate the classical shallow ice flow models by applying simplifications to the weak formulation of the Glen–Stokes model. This approach leads to a unified hierarchical formulation which relates the Glen–Stokes model, the Blatter model, the shallow ice approximation and the shallow shelf approximation. Based on this formulation and on residual techniques commonly used to estimate numerical errors, we derive three a posteriori estimators, each of which compares a pair of models using measures of the velocity field from the simpler (shallower) model. Numerical experiments confirm that these estimators can be used to assess the validity of the shallow ice models that are commonly used in glacier and ice sheet modelling.

In this paper, a multilayer generalisation of the Shallow Shelf Approximation (SSA) is considered. In this recent hybrid ice flow model, the ice thickness is divided into thin layers, which can spread out, contract and slide over each other in such a way that the velocity profile is layer-wise constant. Like the SSA (1-layer model), the multilayer model can be reformulated as a minimisation problem. However, unlike the SSA, the functional to be minimised involves a new penalisation term for the interlayer jumps of the velocity, which represents the vertical shear stresses induced by interlayer sliding. Taking advantage of this reformulation, numerical solvers developed for the SSA can be naturally extended layer-wise or column-wise. Numerical results show that the column-wise extension of a Newton multigrid solver proves to be robust in the sense that its convergence is barely influenced by the number of layers and the type of ice flow. In addition, the multilayer formulation appears to be naturally better conditioned than the one of the first-order approximation to face the anisotropic conditions of the sliding-dominant ice flow of ISMIP-HOM experiments.

A new hybrid model for the dynamics of glaciers, ice sheets and ice shelves is introduced. In this ‘multilayer’ model the domain of ice consists of a pile of thin layers, which can spread out, contract and slide over each other such that the two most relevant types of stresses are accounted for: membrane and vertical shear. Assuming the horizontal velocity field to be vertically piecewise constant in each layer, the model is derived from local depth integrations of the hydrostatic approximation of the Stokes equations. These integrations give rise to interlayer tractions, which can be redefined at zeroth order in the interlayer surface slope by keeping the vertical shear stress components. Furthermore, if the layers are chosen such that they are aligned with the streamlines, then second-order accurate interlayer tractions can replace zeroth-order ones. The final model consists of a tridiagonal system of two-dimensional nonlinear elliptic equations, the size of this system being equal to the number of layers. When running the model for prognostic flowline ISMIP-HOM benchmark experiments, the multilayer solutions show good agreement with the higher-order solutions if no severe depression occurs in the bedrock. As an alternative to three-dimensional models, the multilayer approach offers to glacier and ice sheet modellers a way of upgrading the commonly used shallow shelf approximation model into a mechanically complete but mathematically two-dimensional model.

In this paper we reconstruct the space–time trajectory beneath the surface of Aletschgletscher, Switzerland, of the corpses of three mountaineers that disappeared in March 1926 and reappeared at the glacier surface in June 2012. Our method integrates the time-dependent velocity field of an existing full-Stokes glacier model, starting at the point where the corpses were found at the glacier surface. Our main result is that we were able to localize the immersion location where the brothers presumably died. As a second result, the upstream end point of the computed trajectory emerges very close to the glacier surface in 1926, giving a new and global validation of the glacier model in space and time. Testing the sensitivity of the immersion location obtained with respect to the model and other uncertainties indicates an area of 0.6% of the entire glacier area where the accident could have occurred. Our result suggests that death was not caused by an avalanche or a fall into a crevasse; instead, it is likely that the mountaineers became disoriented in prolonged severe weather conditions and froze to death.

In this paper, we consider a model for the time evolution of three-dimensional marine ice sheets. This model combines the Shallow Ice Approximation (SIA) for the ice deformation, the Shallow Shelf Approximation (SSA) for the basal sliding, and the mass conservation principle. At each time step, we solve a scalar p-Laplace minimization-type problem with obstacle (SIA), a vectorial p-Laplace minimization-type problem (SSA) and a transport equation (mass conservation). The two minimization problems are solved using a truncated nonsmooth Newton multigrid method while the transport equation is solved using a vertex-centred finite volume method. Our approach is combined to an heuristic mesh adaptive refinement procedure to face the large gradients of the solution that are expected between the ice sheet and the ice shelf. As applications, we present some simulations of the Marine Ice Sheet Model Intercomparison Project MISMIP (2D and 3D) and validate our results against an analytic solution (2D) and other participant model results (3D). Further numerical results show that the convergence of our Newton multigrid method is insensitive to local refinements making our overall adaptive strategy fully efficient.

In this paper, we analyse the calving activity of Bowdoin Glacier, north-west Greenland, in 2015 by combining satellite images, UAV (Unmanned Aerial Vehicles) photogrammetry and ice flow modelling. In particular, a high-resolution displacement field is inferred from UAV orthoimages taken immediately before and after the initiation of a large fracture, which induced a major calving event. A detailed analysis of the strain rate field allows us to map accurately the path taken by the opening crack. Modelling results reveal i) that the crack was more than half-thickness deep, filled with water, and getting irreversibly deeper when it was captured by the UAV and ii) that the crack initiated in an area of high horizontal shear caused by a local basal bump immediately behind the current calving front. The asymmetry of the bed at the front explains the systematic calving pattern observed in May and July–August 2015. As a corollary, we infer that the calving front of Bowdoin Glacier is currently stabilized by this bedrock bump and might enter in an unstable mode and retreat rapidly if the glacier keeps thinning in the coming years. Beyond this outcome, our study demonstrates that the combination of UAV photogrammetry and ice flow modelling is a promising tool to track horizontally and vertically the propagation of fractures responsible for large calving events.

We formulate steady, shallow ice sheet flow as an obstacle problem, the unknown being the ice upper surface and the obstacle being the underlying bedrock topography. This generates a free-boundary defining the ice sheet extent. The obstacle problem is written as a variational inequality subject to the positive-ice-thickness constraint. The corresponding PDE is a highly nonlinear elliptic equation which generalizes the $p$-Laplacian equation. Our formulation also permits variable ice softness, basal sliding, and elevation-dependent surface mass balance. Existence and uniqueness are shown in restricted cases which we may reformulate as a convex minimization problem. In the general case we show existence by applying a fixed point argument. Using continuity results from that argument, we construct a numerical solution by solving a sequence of obstacle $p$-Laplacian-like problems by finite element approximation. As a real application, we compute the steady-state shape of the Greenland ice sheet in a steady pre...

The aim of this paper is to study a nonlinear stationary Stokes problem with mixed boundary conditions that describes the ice velocity and pressure fields of grounded glaciers under Glen's flow law. Using convex analysis arguments, we prove the existence and the uniqueness of a weak solution. A finite element method is applied with approximation spaces that satisfy the inf-sup condition, and a priori error estimates are established by using a quasinorm technique. Several algorithms (including Newton's method) are proposed to solve the nonlinearity of the Stokes problem and are proved to be convergent. Our results are supported by numerical convergence studies.

In this chapter, we analyze and approximate a nonlinear stationary Stokes problem that describes the motion of glacier ice. The existence and uniqueness of solutions are proved and an a priori error estimate for the finite element approximation is found. In a second time, we combine the Stokes problem with a transport equation for the volume fraction of ice, which describes the time evolution of a glacier. The accumulation due to snow precipitation and melting are accounted for in the source term of the transport equation. A decoupling algorithm allows the diffusion and the advection problems to be solved using a two-grids method. As an illustration, we simulate the evolution of Aletsch glacier, Switzerland, over the 21st century by using realistic climatic conditions.

The existence of solutions of the non-sliding shallow-ice-sheet equation on a flat horizontal bed with a mass balance linearly depending on altitude is proven for fixed margins. Free-margin solutions for the same mass balance do not exist. Fixed-margin solutions show unbounded shear stress and nonzero mass flux at the margin. Steady-state solutions with realistic margins, vanishing ice flux and vanishing shear stress are found numerically for ice sheets with Weertman-type sliding.

We propose a novel Eulerian algorithm to compute the changes of a glacier geometry for given mass balances. The surface of a glacier is obtained by solving a transport equation for the volume of fluid (VOF). The surface mass balance is taken into account by adding an interfacial term in the transport equation. An unstructured mesh with standard stabilized finite elements is used to solve the non-linear Stokes problem. The VOF function is computed on a structured grid with a high resolution. The algorithm is stable for Courant numbers larger than unity and conserves mass to a high accuracy. To demonstrate the potential of the algorithm, we apply it to reconstructed Late-glacial states of a small valley glacier, Vadret Muragl, in the Swiss Alps.

Due to climatic change, many Alpine glaciers have significantly retreated during the last century. In this study we perform the numerical simulation of the temporal and spatial change of Rhonegletscher, Swiss Alps, from 1874 to 2007, and from 2007 to 2100. Given the shape of the glacier, the velocity of ice u is obtained by solving a 3D nonlinear Stokes problem with a nonlinear sliding law along the bedrock-ice interface. The shape of the glacier is updated by computing the volume fraction of ice @f, which satisfies a transport equation. The accumulation due to snow fall and the ablation due to melting is accounted by adding a source term to the transport equation. A decoupling algorithm allows the two above problems to be solved using different numerical techniques. The nonlinear Stokes problem is solved on a fixed, unstructured finite element mesh consisting of tetrahedrons. The transport equation is solved using a fixed, structured grid of smaller cells. The numerical simulation, from 1874 to 2007, is validated against measurements. Afterwards, three different climatic scenarios are considered in order to predict the shape of Rhonegletscher until 2100. A dramatic retreat of Rhonegletscher during the 21st century is anticipated. Our results contribute to a better understanding of the impact of climatic change on mountain glaciers.

For more than a century Alpine glaciers have been retreating dramatically, and they are expected to shrink even more quickly over the coming decades. This study addresses the future evolution of Grosser Aletschgletscher, Switzerland, the largest glacier in the European Alps. A three-dimensional combined surface mass-balance and glacier dynamics model was applied. The ice flow was described with the full Stokes equations. The glacier surface evolution was obtained by solving a transport equation for the volume of fluid. Daily surface melt and accumulation were calculated on the basis of climate data. The combined model was validated against several types of measurements made throughout the 20th century. For future climate change, scenarios based on regional climate models in the ENSEMBLES project were used. According to the median climatic evolution, Aletschgletscher was expected to lose 90% of its ice volume by the end of 2100. Even when the model was driven using current climate conditions (the past two decades) the glacier tongue experienced a considerable retreat of 6 km, indicating its strong disequilibrium with the present climate. By including a model for the evolution of supraglacial debris and its effect in reducing glacier melt, we show that this factor can significantly slow future glacier retreat.