Search

This study investigates the uncertain future contribution to sea-level rise in response to global warming of Upernavik Isstrøm, a tidewater glacier in Greenland. We analyse multiple sources of uncertainty, including Shared Socioeconomic Pathways (SSPs), climate models (global and regional), ice–ocean interactions, and ice sheet model (ISM) parameters. We use weighting methods based on spatio-temporal velocity and elevation data to reduce ice flow model uncertainty and evaluate their ability to prevent overconfidence. Our developed initialization method demonstrates the capability of Elmer/Ice to accurately replicate the hindcast mass loss of Upernavik Isstrøm. Future mass loss predictions in 2100 range from a contribution to sea-level rise from 1.5 to 7.2 mm, with an already committed sea-level contribution projection from 0.6 to 1.3 mm. At the end of the century, SSP-related uncertainty constitutes the predominant component of total uncertainty, accounting for 40 %, while uncertainty linked to the ISM represents 15 % of the overall uncertainty. We find that calibration does not reduce uncertainty in the future mass loss between today and 2100 (+2 %) but significantly reduces uncertainty in the hindcast mass loss between 1985 and 2015 (−32 % to −61 % depending on the weighting method). Combining calibration of the ice sheet model with SSP weighting yields uncertainty reductions in future mass loss in 2050 (−1.5 %) and in 2100 (−32 %).

Observations of ocean-driven grounding-line retreat in the Amundsen Sea Embayment in Antarctica raise the question of an imminent collapse of the West Antarctic Ice Sheet. Here we analyse the committed evolution of Antarctic grounding lines under the present-day climate. To this aim, we first calibrate a sub-shelf melt parameterization, which is derived from an ocean box model, with observed and modelled melt sensitivities to ocean temperature changes, making it suitable for present-day simulations and future sea level projections. Using the new calibration, we run an ensemble of historical simulations from 1850 to 2015 with a state-of-the-art ice sheet model to create model instances of possible present-day ice sheet configurations. Then, we extend the simulations for another 10 000 years to investigate their evolution under constant present-day climate forcing and bathymetry. We test for reversibility of grounding-line movement in the case that large-scale retreat occurs. In the Amundsen Sea Embayment we find irreversible retreat of the Thwaites Glacier for all our parameter combinations and irreversible retreat of the Pine Island Glacier for some admissible parameter combinations. Importantly, an irreversible collapse in the Amundsen Sea Embayment sector is initiated at the earliest between 300 and 500 years in our simulations and is not inevitable yet – as also shown in our companion paper (Part 1, Hill et al., 2023). In other words, the region has not tipped yet. With the assumption of constant present-day climate, the collapse evolves on millennial timescales, with a maximum rate of 0.9 mm a−1 sea-level-equivalent ice volume loss. The contribution to sea level by 2300 is limited to 8 cm with a maximum rate of 0.4 mm a−1 sea-level-equivalent ice volume loss. Furthermore, when allowing ice shelves to regrow to their present geometry, we find that large-scale grounding-line retreat into marine basins upstream of the Filchner–Ronne Ice Shelf and the western Siple Coast is reversible. Other grounding lines remain close to their current positions in all configurations under present-day climate.

Theoretical and numerical work has shown that under certain circumstances grounding lines of marine-type ice sheets can enter phases of irreversible advance and retreat driven by the marine ice sheet instability (MISI). Instances of such irreversible retreat have been found in several simulations of the Antarctic Ice Sheet. However, it has not been assessed whether the Antarctic grounding lines are already undergoing MISI in their current position. Here, we conduct a systematic numerical stability analysis using three state-of-the-art ice sheet models: Úa, Elmer/Ice, and the Parallel Ice Sheet Model (PISM). For the first two models, we construct steady-state initial configurations whereby the simulated grounding lines remain at the observed present-day positions through time. The third model, PISM, uses a spin-up procedure and historical forcing such that its transient state is close to the observed one. To assess the stability of these simulated states, we apply short-term perturbations to submarine melting. Our results show that the grounding lines around Antarctica migrate slightly away from their initial position while the perturbation is applied, and they revert once the perturbation is removed. This indicates that present-day retreat of Antarctic grounding lines is not yet irreversible or self-sustained. However, our accompanying paper (Part 2, Reese et al., 2023a) shows that if the grounding lines retreated further inland, under present-day climate forcing, it may lead to the eventual irreversible collapse of some marine regions of West Antarctica.

Due to increasing surface melting on the Greenland ice sheet, better constraints on seasonally evolving basal water pressure and sliding speed are required by models. Here we assess the potential of using inverse methods on a dense time series of surface speeds to recover the seasonal evolution of the basal conditions in a well-documented region in southwest Greenland. Using data compiled from multiple satellite missions, we document seasonally evolving surface velocities with a temporal resolution of 2 weeks between 2015 and 2019. We then apply the inverse control method using the ice flow model Elmer/Ice to infer the basal sliding and friction corresponding to each of the 24 surface velocity data sets. Near the margin where the uncertainty in the velocity and bed topography are small, we obtain clear seasonal variations that can be mostly interpreted in terms of an effective-pressure-based hard-bed friction law. We find for valley bottoms or “troughs” in the bed topography that the changes in modelled basal conditions directly respond to local modelled water pressure variations, while the link is more complex for subglacial “ridges” which are often non-locally forced. At the catchment scale, in-phase variations in the water pressure, surface velocities, and surface runoff variations are found. Our results show that time series inversions of observed surface velocities can be used to understand the evolution of basal conditions over different timescales and could therefore serve as an intermediate validation for subglacial hydrology models to achieve better coupling with ice flow models.

On glaciers and ice sheets, identifying the relationship between velocity and traction is critical to constrain the bed physics that controls ice flow. Yet in Greenland, these relationships remain unquantified. We determine the spatial relationship between velocity and traction in all eight major drainage catchments of Greenland. The basal traction is estimated using three different methods over large grid cells to minimize interpretation biases associated with unconstrained rheologic parameters used in numerical inversions. We find the relationships are consistent with our current understanding of basal physics in each catchment. We identify catchments that predominantly show Mohr–Coulomb-like behavior typical of deforming beds or significant cavitation, as well as catchments that predominantly show rate-strengthening behavior typical of Weertman-type hard-bed physics. Overall, the traction relationships suggest that the flow field and surface geometry of the grounded regions in Greenland is mainly dictated by Weertman-type hard-bed physics up to velocities of approximately 450 m yr−1, except within the Northeast Greenland Ice Stream and areas near floatation. Depending on the catchment, behavior of the fastest-flowing ice (∼ 1000 m yr−1) directly inland from marine-terminating outlets exhibits Weertman-type rate strengthening, Mohr–Coulomb-like behavior, or is not confidently resolved given our methodology. Given the complex basal boundary across Greenland, the relationships are captured reasonably well by simple traction laws which provide a parameterization that can be used to model ice dynamics at large scales. The results and analysis serve as a first constraint on the physics of basal motion over the grounded regions of Greenland and provide unique insight into future dynamics and vulnerabilities in a warming climate.

Alpine glaciers are shrinking and rapidly loosing mass in a warming climate. Glacier modeling is required to assess the future consequences of these retreats on water resources, the hydropower industry and risk management. However, the performance of such ice flow modeling is generally difficult to evaluate because of the lack of long-term glaciological observations. Here, we assess the performance of the Elmer/Ice full Stokes ice flow model using the long dataset of mass balance, thickness change, ice flow velocity and snout fluctuation measurements obtained between 1979 and 2015 on the Mer de Glace glacier, France. Ice flow modeling results are compared in detail to comprehensive glaciological observations over 4 decades including both a period of glacier expansion preceding a long period of decay. To our knowledge, a comparison to data at this detail is unprecedented. We found that the model accurately reconstructs the velocity, elevation and length variations of this glacier despite some discrepancies that remain unexplained. The calibrated and validated model was then applied to simulate the future evolution of Mer de Glace from 2015 to 2050 using 26 different climate scenarios. Depending on the climate scenarios, the largest glacier in France, with a length of 20 km, could retreat by 2 to 6 km over the next 3 decades.

Marine-based sectors of the Antarctic Ice Sheet are increasingly contributing to sea level rise. The basal conditions exert an important control on the ice dynamics and can be propitious to instabilities in the grounding line position. Because the force balance is non-inertial, most ice flow models are now equipped with time-independent inverse methods to constrain the basal conditions from observed surface velocities. However, transient simulations starting from this initial state usually suffer from inconsistencies and are not able to reproduce observed trends. Here, using a synthetic flow line experiment, we assess the performance of an ensemble Kalman filter for the assimilation of transient observations of surface elevation and velocities in a marine ice sheet model. The model solves the shallow shelf equation for the force balance and the continuity equation for ice thickness evolution. The position of the grounding line is determined by the floatation criterion. The filter analysis estimates both the state of the model, represented by the surface elevation, and the basal conditions, with the simultaneous inversion of the basal friction and topography. The idealised experiment reproduces a marine ice sheet that is in the early stage of an unstable retreat. Using observation frequencies and uncertainties consistent with current observing systems, we find that the filter allows the accurate recovery of both the basal friction and topography after few assimilation cycles with relatively small ensemble sizes. In addition it is found that assimilating the surface observations has a positive impact on constraining the evolution of the grounding line during the assimilation window. Using the initialised state to perform century-scale forecast simulations, we show that grounding line retreat rates are in agreement with the reference; however remaining uncertainties in the basal conditions may lead to significant delays in the initiation of the unstable retreat. These results are encouraging for the application to real glacial systems.

The majority of Antarctic ice shelves are bounded by grounded ice rises. These ice rises exhibit local flow fields that partially oppose the flow of the surrounding ice shelves. Formation of ice rises is accompanied by a characteristic upward-arching internal stratigraphy (“Raymond arches”), whose geometry can be analysed to infer information about past ice-sheet changes in areas where other archives such as rock outcrops are missing. Here we present an improved modelling framework to study ice-rise evolution using a satellite-velocity calibrated, isothermal, and isotropic 3-D full-Stokes model including grounding-line dynamics at the required mesh resolution (<500 m). This overcomes limitations of previous studies where ice-rise modelling has been restricted to 2-D and excluded the coupling between the ice shelf and ice rise. We apply the model to the Ekström Ice Shelf, Antarctica, containing two ice rises. Our simulations investigate the effect of surface mass balance and ocean perturbations onto ice-rise divide position and interpret possible resulting unique Raymond arch geometries. Our results show that changes in the surface mass balance result in immediate and sustained divide migration (>2.0 m yr−1) of up to 3.5 km. In contrast, instantaneous ice-shelf disintegration causes a short-lived and delayed (by 60–100 years) response of smaller magnitude ( m yr−1). The model tracks migration of a triple junction and synchronous ice-divide migration in both ice rises with similar magnitude but differing rates. The model is suitable for glacial/interglacial simulations on the catchment scale, providing the next step forward to unravel the ice-dynamic history stored in ice rises all around Antarctica.

Oceanic melting beneath ice shelves is the main driver of the current mass loss of the Antarctic ice sheet and is mostly parameterised in stand-alone ice-sheet modelling. Parameterisations are crude representations of reality, and their response to ocean warming has not been compared to 3-D ocean–ice-sheet coupled models. Here, we assess various melting parameterisations ranging from simple scalings with far-field thermal driving to emulators of box and plume models, using a new coupling framework combining the ocean model NEMO and the ice-sheet model Elmer/Ice. We define six idealised one-century scenarios for the far-field ocean ranging from cold to warm, and representative of potential futures for typical Antarctic ice shelves. The scenarios are used to constrain an idealised geometry of the Pine Island glacier representative of a relatively small cavity. Melt rates and sea-level contributions obtained with the parameterised stand-alone ice-sheet model are compared to the coupled model results. The plume parameterisations give good results for cold scenarios but fail and underestimate sea level contribution by tens of percent for warm(ing) scenarios, which may be improved by adapting its empirical scaling. The box parameterisation with five boxes compares fairly well to the coupled results for almost all scenarios, but further work is needed to grasp the correct number of boxes. For simple scalings, the comparison to the coupled framework shows that a quadratic as opposed to linear dependency on thermal forcing is required. In addition, the quadratic dependency is improved when melting depends on both local and non-local, i.e. averaged over the ice shelf, thermal forcing. The results of both the box and the two quadratic parameterisations fall within or close to the coupled model uncertainty. All parameterisations overestimate melting for thin ice shelves while underestimating melting in deep water near the grounding line. Further work is therefore needed to assess the validity of these melting parameteriations in more realistic set-ups.

Ice sheet numerical modeling is an important tool to estimate the dynamic contribution of the Antarctic ice sheet to sea level rise over the coming centuries. The influence of initial conditions on ice sheet model simulations, however, is still unclear. To better understand this influence, an initial state intercomparison exercise (initMIP) has been developed to compare, evaluate, and improve initialization procedures and estimate their impact on century-scale simulations. initMIP is the first set of experiments of the Ice Sheet Model Intercomparison Project for CMIP6 (ISMIP6), which is the primary Coupled Model Intercomparison Project Phase 6 (CMIP6) activity focusing on the Greenland and Antarctic ice sheets. Following initMIP-Greenland, initMIP-Antarctica has been designed to explore uncertainties associated with model initialization and spin-up and to evaluate the impact of changes in external forcings. Starting from the state of the Antarctic ice sheet at the end of the initialization procedure, three forward experiments are each run for 100 years: a control run, a run with a surface mass balance anomaly, and a run with a basal melting anomaly beneath floating ice. This study presents the results of initMIP-Antarctica from 25 simulations performed by 16 international modeling groups. The submitted results use different initial conditions and initialization methods, as well as ice flow model parameters and reference external forcings. We find a good agreement among model responses to the surface mass balance anomaly but large variations in responses to the basal melting anomaly. These variations can be attributed to differences in the extent of ice shelves and their upstream tributaries, the numerical treatment of grounding line, and the initial ocean conditions applied, suggesting that ongoing efforts to better represent ice shelves in continental-scale models should continue.

Reliable projections of ice sheets' future contributions to sea-level rise require models that are able to accurately simulate grounding-line dynamics, starting from initial states consistent with observations. Here, we simulate the centennial evolution of the Amundsen Sea Embayment in response to a prescribed perturbation in order to assess the sensitivity of mass loss projections to the chosen friction law, depending on the initialisation strategy. To this end, three different model states are constructed by inferring both the initial basal shear stress and viscosity fields with various relative weights. Then, starting from each of these model states, prognostic simulations are carried out using a Weertman, a Schoof and a Budd friction law, with different parameter values. Although the sensitivity of projections to the chosen friction law tends to decrease when more weight is put on viscosity during initialisation, it remains significant for the most physically acceptable of the constructed model states. Independently of the considered model state, the Weertman law systematically predicts the lowest mass losses. In addition, because of its particular dependence on effective pressure, the Budd friction law induces significantly different grounding-line retreat patterns than the other laws and predicts significantly higher mass losses.

At the Last Glacial Maximum (LGM), the Rhine glacier in the Swiss Alps covered an area of about 16 000 km2. As part of an integrative study about the safety of repositories for radioactive waste under ice age conditions in Switzerland, we modeled the Rhine glacier using a thermodynamically coupled three-dimensional, transient Stokes flow and heat transport model down to a horizontal resolution of about 500 m. The accumulation and ablation gradients that roughly reproduced the geomorphic reconstructions of glacial extent and ice thickness suggested extremely cold (TJuly ∼ 0 °C at the glacier terminus) and dry ( ∼ 10 % to 20 % of today's precipitation) climatic conditions. Forcing the numerical simulations with warmer and wetter conditions that better matched LGM climate proxy records yielded a glacier on average 500 to 700 m thicker than geomorphic reconstructions. Mass balance gradients also controlled ice velocities, fluxes, and sliding speeds. These gradients, however, had only a small effect on basal conditions. All simulations indicated that basal ice reached the pressure melting point over much of the Rhine and Linth piedmont lobes, and also in the glacial valleys that fed these lobes. Only the outer margin of the lobes, bedrock highs beneath the lobes, and Alpine valleys at high elevations in the accumulation zone remained cold based. The Rhine glacier was thus polythermal. Sliding speed estimated with a linear sliding rule ranged from 20 to 100 m a−1 in the lobes and 50 to 250 m a−1 in Alpine valleys. Velocity ratios (sliding to surface speeds) were > 80 % in lobes and ∼ 60 % in valleys. Basal shear stress was very low in the lobes (0.03–0.1 MPa) and much higher in Alpine valleys ( > 0.2 MPa). In these valleys, viscous strain heating was a dominant source of heat, particularly when shear rates in the ice increased due to flow constrictions, confluences, or flow past large bedrock obstacles, contributing locally up to several watts per square meter but on average 0.03 to 0.2 W m−2. Basal friction acted as a heat source at the bed of about 0.02 W m−2, 4 to 6 times less than the geothermal heat flow which is locally high (up to 0.12 W m−2). In the lobes, despite low surface slopes and low basal shear stresses, sliding dictated main fluxes of ice, which closely followed bedrock topography: ice was channeled in between bedrock highs along troughs, some of which coincided with glacially eroded overdeepenings. These sliding conditions may have favored glacial erosion by abrasion and quarrying. Our results confirmed general earlier findings but provided more insights into the detailed flow and basal conditions of the Rhine glacier at the LGM. Our model results suggested that the trimline could have been buried by a significant thickness of cold ice. These findings have significant implications for interpreting trimlines in the Alps and for our understanding of ice–climate interactions.

The search for ice as old as 1.5 Myr requires the identification of places that maximize our chances to retrieve old, well-resolved, undisturbed and datable ice. One of these locations is very likely southwest of the Dome C summit, where elevated bedrock makes the ice thin enough to limit basal melting. A 3-D ice flow simulation is used to calculate five selection criteria, which together delineate the areas with the most appropriate glaciological properties. These selected areas (a few square kilometers) lie on the flanks of a bedrock high, where a balance is found between risks of basal melting, stratigraphic disturbances and sufficient age resolution. Within these areas, several sites of potential 1.5 Myr old ice are proposed, situated on local bedrock summits or ridges. The trajectories of the ice particles towards these locations are short, and the ice flows over a smoothly undulating bedrock. These sites will help to choose where new high-resolution ground radar surveys should be conducted in upcoming field seasons.

Earlier large-scale Greenland ice sheet sea-level projections (e.g. those run during the ice2sea and SeaRISE initiatives) have shown that ice sheet initial conditions have a large effect on the projections and give rise to important uncertainties. The goal of this initMIP-Greenland intercomparison exercise is to compare, evaluate, and improve the initialisation techniques used in the ice sheet modelling community and to estimate the associated uncertainties in modelled mass changes. initMIP-Greenland is the first in a series of ice sheet model intercomparison activities within ISMIP6 (the Ice Sheet Model Intercomparison Project for CMIP6), which is the primary activity within the Coupled Model Intercomparison Project Phase 6 (CMIP6) focusing on the ice sheets. Two experiments for the large-scale Greenland ice sheet have been designed to allow intercomparison between participating models of (1) the initial present-day state of the ice sheet and (2) the response in two idealised forward experiments. The forward experiments serve to evaluate the initialisation in terms of model drift (forward run without additional forcing) and in response to a large perturbation (prescribed surface mass balance anomaly); they should not be interpreted as sea-level projections. We present and discuss results that highlight the diversity of data sets, boundary conditions, and initialisation techniques used in the community to generate initial states of the Greenland ice sheet. We find good agreement across the ensemble for the dynamic response to surface mass balance changes in areas where the simulated ice sheets overlap but differences arising from the initial size of the ice sheet. The model drift in the control experiment is reduced for models that participated in earlier intercomparison exercises.

The basal topography is largely unknown beneath most glaciers and ice caps, and many attempts have been made to estimate a thickness field from other more accessible information at the surface. Here, we present a two-step reconstruction approach for ice thickness that solves mass conservation over single or several connected drainage basins. The approach is applied to a variety of test geometries with abundant thickness measurements including marine- and land-terminating glaciers as well as a 2400 km2 ice cap on Svalbard. The input requirements are kept to a minimum for the first step. In this step, a geometrically controlled, non-local flux solution is converted into thickness values relying on the shallow ice approximation (SIA). In a second step, the thickness field is updated along fast-flowing glacier trunks on the basis of velocity observations. Both steps account for available thickness measurements. Each thickness field is presented together with an error-estimate map based on a formal propagation of input uncertainties. These error estimates point out that the thickness field is least constrained near ice divides or in other stagnant areas. Withholding a share of the thickness measurements, error estimates tend to overestimate mismatch values in a median sense. We also have to accept an aggregate uncertainty of at least 25 % in the reconstructed thickness field for glaciers with very sparse or no observations. For Vestfonna ice cap (VIC), a previous ice volume estimate based on the same measurement record as used here has to be corrected upward by 22 %. We also find that a 13 % area fraction of the ice cap is in fact grounded below sea level. The former 5 % estimate from a direct measurement interpolation exceeds an aggregate maximum range of 6–23 % as inferred from the error estimates here.

Knowledge of the ice thickness distribution of glaciers and ice caps is an important prerequisite for many glaciological and hydrological investigations. A wealth of approaches has recently been presented for inferring ice thickness from characteristics of the surface. With the Ice Thickness Models Intercomparison eXperiment (ITMIX) we performed the first coordinated assessment quantifying individual model performance. A set of 17 different models showed that individual ice thickness estimates can differ considerably – locally by a spread comparable to the observed thickness. Averaging the results of multiple models, however, significantly improved the results: on average over the 21 considered test cases, comparison against direct ice thickness measurements revealed deviations on the order of 10 ± 24 % of the mean ice thickness (1σ estimate). Models relying on multiple data sets – such as surface ice velocity fields, surface mass balance, or rates of ice thickness change – showed high sensitivity to input data quality. Together with the requirement of being able to handle large regions in an automated fashion, the capacity of better accounting for uncertainties in the input data will be a key for an improved next generation of ice thickness estimation approaches.

Ice flow models are now routinely used to forecast the ice sheets' contribution to 21st century sea-level rise. For such short term simulations, the model response is greatly affected by the initial conditions. Data assimilation algorithms have been developed to invert for the friction of the ice on its bedrock using observed surface velocities. A drawback of these methods is that remaining uncertainties, especially in the bedrock elevation, lead to non-physical ice flux divergence anomalies resulting in undesirable transient effects. In this study, we compare two different assimilation algorithms based on adjoints and nudging to constrain both bedrock friction and elevation. Using synthetic twin experiments with realistic observation errors, we show that the two algorithms lead to similar performances in reconstructing both variables and allow the flux divergence anomalies to be significantly reduced.

Three-dimensional ice flow modelling requires a large number of computing resources and observation data, such that 2-D simulations are often preferable. However, when there is significant lateral divergence, this must be accounted for (2.5-D models), and a flow tube is considered (volume between two horizontal flowlines). In the absence of velocity observations, this flow tube can be derived assuming that the flowlines follow the steepest slope of the surface, under a few flow assumptions. This method typically consists of scanning a digital elevation model (DEM) with a moving window and computing the curvature at the centre of this window. The ability of the 2.5-D models to account properly for a 3-D state of strain and stress has not clearly been established, nor their sensitivity to the size of the scanning window and to the geometry of the ice surface, for example in the cases of sharp ridges. Here, we study the applicability of a 2.5-D ice flow model around a dome, typical of the East Antarctic plateau conditions. A twin experiment is carried out, comparing 3-D and 2.5-D computed velocities, on three dome geometries, for several scanning windows and thermal conditions. The chosen scanning window used to evaluate the ice surface curvature should be comparable to the typical radius of this curvature. For isothermal ice, the error made by the 2.5-D model is in the range 0–10 % for weakly diverging flows, but is 2 or 3 times higher for highly diverging flows and could lead to a non-physical ice surface at the dome. For non-isothermal ice, assuming a linear temperature profile, the presence of a sharp ridge makes the 2.5-D velocity field unrealistic. In such cases, the basal ice is warmer and more easily laterally strained than the upper one, the walls of the flow tube are not vertical, and the assumptions of the 2.5-D model are no longer valid.

Radar inference of the bulk properties of glacier beds, most notably identifying basal melting, is, in general, derived from the basal reflection coefficient. On the scale of an ice sheet, unambiguous determination of basal reflection is primarily limited by uncertainty in the englacial attenuation of the radio wave, which is an Arrhenius function of temperature. Existing bed-returned power algorithms for deriving attenuation assume that the attenuation rate is regionally constant, which is not feasible at an ice-sheet-wide scale. Here we introduce a new semi-empirical framework for deriving englacial attenuation, and, to demonstrate its efficacy, we apply it to the Greenland Ice Sheet. A central feature is the use of a prior Arrhenius temperature model to estimate the spatial variation in englacial attenuation as a first guess input for the radar algorithm. We demonstrate regions of solution convergence for two input temperature fields and for independently analysed field campaigns. The coverage achieved is a trade-off with uncertainty and we propose that the algorithm can be "tuned" for discrimination of basal melt (attenuation loss uncertainty ∼ 5 dB). This is supported by our physically realistic ( ∼ 20 dB) range for the basal reflection coefficient. Finally, we show that the attenuation solution can be used to predict the temperature bias of thermomechanical ice sheet models and is in agreement with known model temperature biases at the Dye 3 ice core.

The dynamical contribution of marine ice sheets to sea level rise is largely controlled by grounding line (GL) dynamics. Two marine ice sheet model intercomparison exercises, namely MISMIP and MISMIP3d, have been proposed to the community to test and compare the ability of models to capture the GL dynamics. Both exercises are known to present a discontinuity of the friction at the GL, which is believed to increase the model sensitivity to mesh resolution. Here, using Elmer/Ice, the only Stokes model which completed both intercomparisons, the sensitivity to the mesh resolution is studied from an extended MISMIP experiment in which the friction continuously decreases over a transition distance and equals zero at the GL. Using this MISMIP-like setup, it is shown that the sensitivity to the mesh resolution is not improved for a vanishing friction at the GL. For the original MISMIP experiment, i.e. for a discontinuous friction at the GL, we further show that the results are moreover very sensitive to the way the friction is interpolated in the close vicinity of the GL. In the light of these new insights, and thanks to increased computing resources, new results for the MISMIP3d experiments obtained for higher resolutions than previously published are made available for future comparisons as the Supplement.

In ice flow modelling, the use of control methods to assimilate the dynamic and geometric state of an ice body has become common practice. These methods have primarily focussed on inverting for one of the two least known properties in glaciology, namely the basal friction coefficient or the ice viscosity parameter. Here, we present an approach to infer both properties simultaneously for the whole of the Antarctic ice sheet. After the assimilation, the root-mean-square deviation between modelled and observed surface velocities attains 8.7 m a−1 for the entire domain, with a slightly higher value of 14.0 m a−1 for the ice shelves. An exception in terms of the velocity mismatch is the Thwaites Glacier Ice Shelf, where the RMS value is almost 70 m a−1. The reason is that the underlying Bedmap2 geometry ignores the presence of an ice rise, which exerts major control on the dynamics of the eastern part of the ice shelf. On these grounds, we suggest an approach to account for pinning points not included in Bedmap2 by locally allowing an optimisation of basal friction during the inversion. In this way, the velocity mismatch on the ice shelf of Thwaites Glacier is more than halved. A characteristic velocity mismatch pattern emerges for unaccounted pinning points close to the marine shelf front. This pattern is exploited to manually identify seven uncharted features around Antarctica that exert significant resistance to the shelf flow. Potential pinning points are detected on Fimbul, West, Shackleton, Nickerson and Venable ice shelves. As pinning points can provide substantial resistance to shelf flow, with considerable consequences if they became ungrounded in the future, the model community is in need of detailed bathymetry there. Our data assimilation points to some of these dynamically important features not present in Bedmap2 and implicitly quantifies their relevance.

Understanding the response of fast flowing ice streams or outlet glaciers to changing climate is crucial in order to make reliable projections of sea level change over the coming decades. Motion of fast outlet glaciers occurs largely through basal motion governed by physical processes at the glacier bed, which are not yet fully understood. Various subglacial mechanisms have been suggested for fast flow but common to most of the suggested processes is the requirement of presence of liquid water, and thus temperate conditions. We use a combination of modelling, field, and remote observations in order to study links between different heat sources, the thermal regime and basal sliding in fast flowing areas on Vestfonna ice cap. A special emphasis lies on Franklinbreen, a fast flowing outlet glacier which has been observed to accelerate recently. We use the ice flow model Elmer/Ice including a Weertman type sliding law and a Robin inverse method to infer basal friction parameters from observed surface velocities. Firn heating, i.e. latent heat release through percolation of melt water, is included in our model; its parameterisation is calibrated with the temperature record of a deep borehole. We found that strain heating is negligible, whereas friction heating is identified as one possible trigger for the onset of fast flow. Firn heating is a significant heat source in the central thick and slow flowing area of the ice cap and the essential driver behind the ongoing fast flow in all outlets. Our findings suggest a possible scenario of the onset and maintenance of fast flow on the Vestfonna ice cap based on thermal processes and emphasise the role of latent heat released through refreezing of percolating melt water for fast flow. However, these processes cannot yet be captured in a temporally evolving sliding law. In order to simulate correctly fast flowing outlet glaciers, ice flow models not only need to account fully for all heat sources, but also need to incorporate a sliding law that is not solely based on the basal temperature, but also on hydrology and/or sediment physics.

Fabric (distribution of crystallographic orientations) along the full NEEM ice core, Greenland was measured in the field by an automatic ice texture analyzer every 10 m, from 33 m down to 2461 m depth. The fabric evolves from a slightly anisotropic fabric at the top, toward a strong single maximum at about 2300 m, which is typical of a deformation pattern mostly driven by uniaxial compression and simple shearing. A sharp increase in the fabric strengthening rate is observed at the Holocene to Wisconsin (HW) climatic transition. From a simple model we estimate that this depth is located at a transition from a state dominated by vertical compression to a state dominated by vertical shear. Comparisons are made to two others ice cores drilled along the same ridge; the GRIP ice core, drilled at the summit of the ice sheet, and the NGRIP ice core, drilled 325 km to the NNW of the summit along the ridge, and 365 km upstream from NEEM. This comparison tends to demonstrate that the ice viscosity change with the HW climatic transition must be associated with the shear-dominated state to induce the abrupt fabric strengthening observed at NEEM. This comparison therefore reflects the increasing role of shear deformation on the coring site when moving NW along the ridge from GRIP to NGRIP and NEEM. The difference in fabric profiles between NEEM and NGRIP also evidences a stronger lateral extension associated with a sharper ridge at NGRIP. Further along the core, centimeter scale abrupt texture (fabric and microstructure) variations are observed in the bottom part of the core. Their positions are in good agreement with the observed folding layers in Dahl-Jensen et al. (2013).

The grounding line is a key element of coastal outlet glaciers, acting on their dynamics. Accurately knowing its position is fundamental for both modelling the glacier dynamics and establishing a benchmark for later change detection. Here we map the grounding line of the Astrolabe Glacier in East Antarctica (66°41' S, 140°05' E), using both hydrostatic and tidal methods. The first method is based on new surface and ice thickness data from which the line of buoyant floatation is found. The second method uses kinematic GPS measurements of the tidal response of the ice surface. By detecting the transitions where the ice starts to move vertically in response to the tidal forcing we determine control points for the grounding line position along GPS profiles. Employing a two-dimensional elastic plate model, we compute the rigid short-term behaviour of the ice plate and estimate the correction required to compare the kinematic GPS control points with the previously determined line of floatation. These two approaches show consistency and lead us to propose a grounding line for the Astrolabe Glacier that significantly deviates from the lines obtained so far from satellite imagery.

We apply a new parameterisation of the Greenland ice sheet (GrIS) feedback between surface mass balance (SMB: the sum of surface accumulation and surface ablation) and surface elevation in the MAR regional climate model (Edwards et al., 2014) to projections of future climate change using five ice sheet models (ISMs). The MAR (Modèle Atmosphérique Régional: Fettweis, 2007) climate projections are for 2000–2199, forced by the ECHAM5 and HadCM3 global climate models (GCMs) under the SRES A1B emissions scenario. The additional sea level contribution due to the SMB–elevation feedback averaged over five ISM projections for ECHAM5 and three for HadCM3 is 4.3% (best estimate; 95% credibility interval 1.8–6.9%) at 2100, and 9.6% (best estimate; 95% credibility interval 3.6–16.0%) at 2200. In all results the elevation feedback is significantly positive, amplifying the GrIS sea level contribution relative to the MAR projections in which the ice sheet topography is fixed: the lower bounds of our 95% credibility intervals (CIs) for sea level contributions are larger than the "no feedback" case for all ISMs and GCMs. Our method is novel in sea level projections because we propagate three types of modelling uncertainty – GCM and ISM structural uncertainties, and elevation feedback parameterisation uncertainty – along the causal chain, from SRES scenario to sea level, within a coherent experimental design and statistical framework. The relative contributions to uncertainty depend on the timescale of interest. At 2100, the GCM uncertainty is largest, but by 2200 both the ISM and parameterisation uncertainties are larger. We also perform a perturbed parameter ensemble with one ISM to estimate the shape of the projected sea level probability distribution; our results indicate that the probability density is slightly skewed towards higher sea level contributions.

We present a new parameterisation that relates surface mass balance (SMB: the sum of surface accumulation and surface ablation) to changes in surface elevation of the Greenland ice sheet (GrIS) for the MAR (Modèle Atmosphérique Régional: Fettweis, 2007) regional climate model. The motivation is to dynamically adjust SMB as the GrIS evolves, allowing us to force ice sheet models with SMB simulated by MAR while incorporating the SMB–elevation feedback, without the substantial technical challenges of coupling ice sheet and climate models. This also allows us to assess the effect of elevation feedback uncertainty on the GrIS contribution to sea level, using multiple global climate and ice sheet models, without the need for additional, expensive MAR simulations. We estimate this relationship separately below and above the equilibrium line altitude (ELA, separating negative and positive SMB) and for regions north and south of 77° N, from a set of MAR simulations in which we alter the ice sheet surface elevation. These give four "SMB lapse rates", gradients that relate SMB changes to elevation changes. We assess uncertainties within a Bayesian framework, estimating probability distributions for each gradient from which we present best estimates and credibility intervals (CI) that bound 95% of the probability. Below the ELA our gradient estimates are mostly positive, because SMB usually increases with elevation: 0.56 (95% CI: −0.22 to 1.33) kg m−3 a−1 for the north, and 1.91 (1.03 to 2.61) kg m−3 a−1 for the south. Above the ELA, the gradients are much smaller in magnitude: 0.09 (−0.03 to 0.23) kg m−3 a−1 in the north, and 0.07 (−0.07 to 0.59) kg m−3 a−1 in the south, because SMB can either increase or decrease in response to increased elevation. Our statistically founded approach allows us to make probabilistic assessments for the effect of elevation feedback uncertainty on sea level projections (Edwards et al., 2014).

The Fourth IPCC Assessment Report concluded that ice sheet flow models, in their current state, were unable to provide accurate forecast for the increase of polar ice sheet discharge and the associated contribution to sea level rise. Since then, the glaciological community has undertaken a huge effort to develop and improve a new generation of ice flow models, and as a result a significant number of new ice sheet models have emerged. Among them is the parallel finite-element model Elmer/Ice, based on the open-source multi-physics code Elmer. It was one of the first full-Stokes models used to make projections for the evolution of the whole Greenland ice sheet for the coming two centuries. Originally developed to solve local ice flow problems of high mechanical and physical complexity, Elmer/Ice has today reached the maturity to solve larger-scale problems, earning the status of an ice sheet model. Here, we summarise almost 10 yr of development performed by different groups. Elmer/Ice solves the full-Stokes equations, for isotropic but also anisotropic ice rheology, resolves the grounding line dynamics as a contact problem, and contains various basal friction laws. Derived fields, like the age of the ice, the strain rate or stress, can also be computed. Elmer/Ice includes two recently proposed inverse methods to infer badly known parameters. Elmer is a highly parallelised code thanks to recent developments and the implementation of a block preconditioned solver for the Stokes system. In this paper, all these components are presented in detail, as well as the numerical performance of the Stokes solver and developments planned for the future.

Over the last two decades, the Greenland ice sheet (GrIS) has been losing mass at an increasing rate, enhancing its contribution to sea-level rise (SLR). The recent increases in ice loss appear to be due to changes in both the surface mass balance of the ice sheet and ice discharge (ice flux to the ocean). Rapid ice flow directly affects the discharge, but also alters ice-sheet geometry and so affects climate and surface mass balance. Present-day ice-sheet models only represent rapid ice flow in an approximate fashion and, as a consequence, have never explicitly addressed the role of ice discharge on the total GrIS mass balance, especially at the scale of individual outlet glaciers. Here, we present a new-generation prognostic ice-sheet model which reproduces the current patterns of rapid ice flow. This requires three essential developments: the complete solution of the full system of equations governing ice deformation; a variable resolution unstructured mesh to resolve outlet glaciers and the use of inverse methods to better constrain poorly known parameters using observations. The modelled ice discharge is in good agreement with observations on the continental scale and for individual outlets. From this initial state, we investigate possible bounds for the next century ice-sheet mass loss. We run sensitivity experiments of the GrIS dynamical response to perturbations in climate and basal lubrication, assuming a fixed position of the marine termini. We find that increasing ablation tends to reduce outflow and thus decreases the ice-sheet imbalance. In our experiments, the GrIS initial mass (im)balance is preserved throughout the whole century in the absence of reinforced forcing, allowing us to estimate a lower bound of 75 mm for the GrIS contribution to SLR by 2100. In one experiment, we show that the current increase in the rate of ice loss can be reproduced and maintained throughout the whole century. However, this requires a very unlikely perturbation of basal lubrication. From this result we are able to estimate an upper bound of 140 mm from dynamics only for the GrIS contribution to SLR by 2100.

The dynamics of Vestfonna ice cap (Svalbard) are dominated by fast-flowing outlet glaciers. Its mass balance is poorly known and affected dynamically by these fast-flowing outlet glaciers. Hence, it is a challenging target for ice flow modeling. Precise knowledge of the basal conditions and implementation of a good sliding law are crucial for the modeling of this ice cap. Here we use the full-Stokes finite element code Elmer/Ice to model the 3-D flow over the whole ice cap. We use a Robin inverse method to infer the basal friction from the surface velocities observed in 1995. Our results illustrate the importance of the basal friction parameter in reproducing observed velocity fields. We also show the importance of having variable basal friction as given by the inverse method to reproduce the velocity fields of each outlet glacier – a simple parametrization of basal friction cannot give realistic velocities in a forward model. We study the robustness and sensitivity of this method with respect to different parameters (mesh characteristics, ice temperature, errors in topographic and velocity data). The uncertainty in the observational parameters and input data proved to be sufficiently small as not to adversely affect the fidelity of the model.

Variegated Glacier (Alaska) is known to surge periodically after a sufficient amount of cumulative mass balance is reached, but this observation is difficult to link with changes in the basal conditions. Here, using a 10-yr dataset, consisting of surface topography and surface velocity observations along a flow line for 25 dates, we have reconstructed the evolution of the basal conditions prior to and during the 1982–1983 surge. The model solves the full-Stokes problem along the central flow line using the finite element method. For the 25 dates of the dataset, the basal friction parameter distribution is inferred using the inverse method proposed by Arthern and Gudmundsson (2010). This method is here slightly modified by incorporating a regularisation term in the cost function to avoid short wavelength changes in the friction parameter. Our results indicate that dramatic changes in the basal conditions occurred between 1973 to 1983. Prior to the surge, periodic changes can be observed between winter and summer, with a regular increase of the sliding from 1973 to 1982. During the surge, the basal friction decreased dramatically and an area of very low friction moved from the upper part of the glacier to its terminus. Using a more complex friction law, these changes in basal sliding are then interpreted in terms of basal water pressure. Our results support that dramatic changes took place in the subglacial drainage system of Variegated Glacier, moving from a relatively efficient drainage system prior to the surge to an inefficient one during the surge. By reconstructing the water pressure evolution at the base of the glacier it is possible to propose a scenario for the hydrological history leading to the occurrence of a surge.