Numerical ice-sheet model experiments sometimes exhibit asymmetries in the solutions, despite the symmetric conditions imposed. By first identifying arithmetic asymmetry in the models as one of the reasons for symmetry-breaking through loss of trailing digits, this paper presents a numerical procedure to preserve the symmetries by restructuring the order of the floating-point evaluation of the equations in the numerical ice-sheet model. Re-examination of the series of experiments in the HEINO topic of the ISMIP demonstrates that small perturbations triggered by arithmetic asymmetries significantly amplify and cause qualitative differences in the simulated ice-sheet evolutions. This study shows that it is imperative to apply a symmetric scheme to maintain overall symmetries in the simulation of ice-sheet evolution, at least under a highly idealized configuration.

Recent observations and modeling studies emphasize the crucial role of fracture mechanics for the stability of ice shelves and thereby the evolution of ice sheets. Here we introduce a macroscopic fracture-density field into a prognostic continuum ice-flow model and compute its evolution incorporating the initiation and growth of fractures as well as their advection with two dimensional ice flow. To a first approximation, fracture growth is assumed to depend on the spreading rate only, while fracture initiation is defined in terms of principal stresses. The inferred fracture-density fields compare well with observed elongate surface structures. Since crevasses and other deep-reaching fracture structures have been shown to influence the overall ice-shelf dynamics, we propose the fracture- density field introduced here be used as a measure for ice softening and decoupling of the ice flow in fracture-weakened zones. This may yield more accurate and realistic velocity patterns in prognostic simulations. Additionally, the memory of past fracture events links the calving front to the upstream dynamics. Thus the fracture-density field proposed here may be employed in fracture-based calving parameterizations. The aim of this study is to introduce the field and investigate which of the observed surface structures can be reproduced by the simplest physically motivated fracture source terms.

Numerical experiments are conducted on a synthetic topography with a three-dimensional thermomechanically coupled ice-sheet model, the Parallel Ice Sheet Model (PISM). Within the model, combined stress balances are connected to evolving thermodynamics and hydrology. The sensitivity of cyclic behaviour to changes in sliding-law parameters and the climate input is studied. Multiple types of oscillations were found, with strong variations in both amplitude and frequency. A physical description is given, in which these variations and transitions from one oscillation type to another are linked to the interplay of stresses, heat transport and hydrological processes. High-frequency oscillations (period 114-169 years), which are shown to have a major impact on ice velocities and a small effect on the ice volume, are related to variations in the water distribution at the base. Low-frequency cycles (period 1000+ years), which have a major impact on both velocities and ice volume, are linked to changes in the thermal regime. Oscillation characteristics are shown to be strongly sensitive to changes in sliding-law parameters and the prescribed surface temperature and mass balance. Incorporating a surface-height dependence of the mass balance is shown to provide an additional feedback, which may induce long- period oscillations.

Glacier mass-balance parameters such as the equilibrium-line altitude (ELA) play an important role when working with large glacier samples. While the number of observational mass-balance series to derive such parameters is limited, more and more modeled data are becoming available.Here we explore the possibilities of analyzing such 'synthetic' mass-balance data with respect to mass-balance parameters. A simplified energy-balance model is driven by bias-corrected regional climate model output to model mass-balance distributions for 94 glaciers in the Swiss Alps over 15 years. The modeling results in realistic interannual variability and mean cumulative mass balance. Subsequently model output is analyzed with respect to 18 topographic and mass-balance parameters and a correlation analysis is performed. Well-known correlations such as for ELA and median elevation are confirmed from the synthetic data. Furthermore, previously unreported parameter relationships are found such as a correlation of the balance rate at the tongue with the accumulation-area ratio (AAR) and of the glacier elevation range with the AAR. Analyzing modeled data complements in situ observations and highlights their importance: the small number of accurate mass-balance observations available for validation is a major challenge for the presented approach.

We extend through 2011 an ice-sheet elevation and surface velocity record across three measurement networks established in south-central Greenland by The Ohio State University in 1980/81. Surface parameters are derived from repeat GPS in situ observations, elevations measured by airborne laser altimetry and by the Ice, Cloud and land Elevation Satellite (ICESat). Elevations at the western network steadily rose early in the record by 0.10 ± 0.02 m a-1, but an eastward-progressing thinning trend began in the mid-1990s followed by a ~1m elevation drop at all stations from 2005 to 2011. Measurements weakly suggest a surface velocity increase at the western cluster from 1980 to 2005. At the central network, elevations rose by 0.08 ± 0.02 m a-1 through 2005 and surface speed increased by 0.5–0.7 m a-1. Surface elevations at the central network remained nearly constant thereafter. Thickening occurred at the southern ice divide by 0.05 ± 0.02 m a-1, while east of the divide the ice sheet thinned with increasing rate from the divide, likely because of decreasing accumulation rate trends and drawdown into rapidly retreating coastal glaciers. Our most recent data show that thinning rates are slowing at several sites just east of the divide and that the elevation at the divide continues to increase.

The transfer of surface-generated meltwater to the subglacial drainage system through full ice thickness crevassing may lead to accelerated glacier velocities, with implications for ice motion under future climatic scenarios. Accurate predictions of where surface meltwater accesses the ice/bed interface are therefore needed in fully coupled hydrodynamic ice-sheet models. We present a spatially distributed modelling routine for predicting the location and timing of delivery of surface-derived meltwater to the ice/bed interface through moulins and supraglacial lake drainage. The model is explained as it is applied to the Croker Bay glacial catchment of Devon Ice Cap, Canada. The formation of moulins, drainage of lakes, and the transfer of meltwater through the full ice thickness are modelled for the 2004 and 2006 ablation seasons. Through this case study we assess the model’s sensitivity to degree-day factors, fracture toughness, tensile strength and crevasse width, and confirm that parameters influencing the rate at which water fills a crevasse are the most significant controls on the ability of a crevasse to reach the bed. Increased surface melt production, therefore, has the potential to significantly influence the spatial and temporal transfer of meltwater through surface-to-bed connections in a warmer climate

Radar data (center frequency 150 MHz) collected on the Antarctic plateau near the EPICA deep-drilling site in Dronning Maud Land vary systematically in backscattered power, depending on the azimuth antenna orientation. Backscatter extrema are aligned with the principal directions of surface strain rates and change with depth. In the upper 900m, backscatter is strongest when the antenna polarization is aligned in the direction of maximal compression, while below 900 m the maxima shift by 90° pointing towards the lateral flow dilatation. We investigate the backscatter from elongated air bubbles and a vertically varying crystal-orientation fabric (COF) using different scattering models in combination with ice-core data. We hypothesize that short-scale variations in COF are the primary mechanism for the observed anisotropy, and the 900 m boundary between the two regimes is caused by ice with varying impurity content. Observations of this kind allow the deduction of COF variations with depth and are potentially also suited to map the transition between Holocene and glacial ice.

Snow avalanches often form levees and en-echelon shear planes in the run-out zone. We describe the formation of these depositional structures using a simple model that accounts for the role of granular fluctuations in avalanche motion. A mathematical feature of this model is the existence of a bifurcation saddle point, describing how granular fluctuations control the avalanche velocity in the runout zone. The saddle point discriminates between a flowing and stopping regime and defines the physical boundary between the flow and non-flow regions of the avalanche, i.e. the location of shear planes in the avalanche deposits. The formation of a shear plane depends on the interplay between terrain slope and avalanche mass flux, which varies from avalanche head to tail. Levees can form immediately at the avalanche front or, for steep slopes and low mass fluxes, at the avalanche tail. At ravine and gully shoulders the mass flux is restricted, thus initiating levee formation. We find that the levee lines are parallel to the flow direction when the mass flux is constant; en-echelon shear lines occur when the mass flux is decreasing. We test the model using several case studies where we have accurate laser scans of avalanche deposits. Our results suggest that avalanche flow parameters can be determined from simple levee measurements or, conversely, formation of levees and flow fingers can be predicted once the parameters governing the granular fluctuations are known.

Remote-sensing and GIS techniques in conjunction with field investigations show how glacier mass loss has led to the rapid growth of Linggo Co, a glacier-fed lake on the central Tibetan Plateau, which has expanded by 21.3% in area between 1974 and 2010, with a lake-level rise of ˜11.2m. The lake volume of Linggo Co increased at a rate of 0.02 × 106, 42.67 × 106 and 65.8 × I06m3a-1 during the periods 1974-92, 1992-99 and 1999-2010, respectively. Other nonglacier-fed lakes in the vicinity (i.e. Longwei Co, Amur Co and Darngo Co Ngion) shrank considerably from the early 1970s to 1992 and then expanded from 1992 to 2010. Despite being in the same climate region, Linggo Co and the non-glacier-fed lakes have differed in response to climate change. The glaciers in the catchment of Linggo Co retreated by 2.4% in area between 1974 and 2007, and their mean thickness decreased by 6.19 ± 1.91 m between 1974 and 2000, with an associated glacier meltwater runoff of (7.52 ± 2.32) × 108 m3. The results indicate that glacier mass loss had a significant impact on the growth of Linggo Co over the past 40 years.

Ice flow numerical models are essential for predicting the evolution of ice sheets in a warming climate. Recent research emphasizes the need for higher-order and even full-Stokes flow models, instead of the traditional shallow-ice approximation, whose assumptions are not valid in certain critical areas. These higher-order models are, however, computationally intensive and difficult to use at the continental scale. Here we present a new technique, the Tiling method, to couple ice flow models of varying orders of complexity. The goal of the method is to limit the spatial extent of where higherorder models are applied to reduce the computational cost, while maintaining the model precision. We apply this method on synthetic geometries to demonstrate its practical use. We first use a geometry for which all models yield the same results to check the consistency of the method. Then we apply our method to a geometry for which a full-Stokes model is required in the vicinity of the ice front. Our results show that the hybrid models present significant improvements over mono-model approaches and reduce computational times.

The loss of the floating ice tongue on Jakobshavn Isbrӕ, Greenland, in the early 2000s has been concurrent with a pattern of thinning, retreat and acceleration leading to enhanced contribution to global sea level. These changes on decadal timescales have been well documented. Here we identify how the glacier responds to forcings on shorter timescales, such as from variations in surface melt, the drainage of supraglacial lakes and seasonal fluctuations in terminus position. Ice motion and surface melt were monitored intermittently from 2006 to 2008. Dual-frequency GPS were deployed 20–50 km upstream of the terminus along the glacier center line. Gaps in surface melt measurements were filled using a temperature-index model of ablation driven by surface air temperatures recorded during the same time period. Our results corroborate the premise that the primary factors controlling speeds on Jakobshavn Isbrӕ are terminus position and geometry. We also observe that surface speeds demonstrate a complex relationship with meltwater input: on diurnal timescales, velocities closely match changes in water input; however, on seasonal timescales a longer, more intense melt season was observed to effectively reduce the overall ice flow of the glacier for the whole year

An analytical theory is developed for ice flow velocity in a boundary layer couplet at the calving front. The theory has simple quantitative characteristics that relate ice front velocity to thickness, strain rate and shelf width, matching one set of empirically derived relationships (Alley and others, 2008) and implying that these relationships predict ice velocity rather than calving rate. The two boundary layers are where longitudinal and transverse flow fields change from the interior flow to patterns consistent with the calving-front stress condition. Numerical simulations confirm the analytical theory. The quantitative predictions of the theory have low sensitivity to unmeasured parameters and to shelf plan aspect ratio, while its robustness arises from its dependence on the scale invariance of the governing equations. The theory provides insights into calving, the stability of ice-shelf calving fronts, the stability of the grounding line of laterally resisted ice streams, and also suggests that the calving front is an instructive dynamical analogue to the grounding line.

We interpreted flow dynamics of the Kahiltna Pass Basin accumulation zone on Mount McKinley, Alaska, USA, using 40,100 and 900 MHz ground-penetrating radar profiles and GPS surface velocity measurements. We found dipping, englacial surface-conformable strata that experienced vertical thickening as the glacier flowed westward from a steep, higher-velocity (60 m a-1) region into flat terrain associated with a 90° bend in the glacier and lower velocities (15 m a-1) to the south. Stratigraphy near the western side of the basin was surface-conformable to ˜170m depth and thinned as flow diverged southward, down-glacier. We found complex strata beneath the conformable stratigraphy and interpret these features as buried crevasses, avalanche debris and deformed ice caused by up-glacier events. We also suggest that basin dimensions, bed topography and the sharp bend each cause flow extension and compression, significantly contributing to conformable and complex strata thickness variations. Our findings show that surface-conformable stratigraphy continuous with depth and consistent strata thicknesses cannot be assumed in accumulation basins, because local and up- glacier terrain and flow dynamics can cause structural complexities to occur under and within surface- conformable layers.

A surface roughness measurement system for snow is presented. It is based on a background board with scales on the edges and a digital camera. Analysis software is developed for automatic processing of images to produce calibrated profiles. The image analysis and calibration was fully automatic in >99% of the studied cases. In the others, the intensity adjustment or board detection needed manual intervention. Profile detection, control point picking and calibration always worked autonomously. The accuracy of the system depends on the photographing configuration, and is typically of the order of 0.1 mm vertically and 0.04 mm horizontally. The method tolerates relatively well cases of snowfall, traces of wiping the black background dry, uneven shading, reflected sunlight, reflected flash light, litter on the snow surface and a tilted plate. The repeatability of the system is at least 1%.

A dark region in the western ablation zone of the Greenland ice sheet is caused by outcropping ice layers that contain more dust than the surrounding brighter ice. These higher amounts of dust were deposited in the accumulation zone of the ice sheet and travelled with the ice to the ablation zone. To deduce the period and the causes of this higher dust deposition, carbonaceous particles in ice samples from the dark region and from brighter reference ice were analysed and used for dating. Samples including ice from directly below the surface contain high amounts of modern organic carbon, probably from microorganisms on the ice surface. Deeper samples reveal low amounts of carbonaceous particles, which are originally deposited in the accumulation zone. The amount of outcropping carbonaceous particles in the dark region seems significantly higher than in the reference ice. One of the samples that contained material initially deposited in the accumulation zone was dated and revealed Late Holocene ages, coinciding with a period of enhanced eolian activity in the nearby tundra. Therefore, variable eolian activity during the Holocene effected dust fluxes towards the ice and hence leads to albedo variations in the present ablation zone of the ice sheet.

When an ice-flow model is constrained by data that exist over only a section of an ice sheet, it is computationally advantageous to limit the model domain to only that section. For example, a limited domain near an ice-core site might cross an ice divide, and have no termini. Accurately calculating ice-sheet evolution in response to spatial and temporal variations in climate and ice flow depends on accurately calculating the transient ice flux crossing the limited-domain boundaries. In the absence of information from outside the limited domain, this is an ill-posed problem. Boundary conditions based only on information from inside the limited domain can produce ice-sheet evolution incompatible with the full ice sheet within which we suppose it to exist. We use impulse-response functions to provide boundary values that are informed by the external ice sheet, without conventionally 'nesting' the limited domain in a full ice-sheet model. Evolution within a limited domain can then be consistent with evolution of boundary conditions is designed for future use in affected the limited domain can be inferred from the full ice sheet. Our treatment of limited-domain an inverse problem in which external changes that data from within the limited domain.

A numerical model for an interacting ice shelf and ocean is presented in which the ice- shelf base exhibits a channelized morphology similar to that observed beneath Petermann Gletscher’s (Greenland) floating ice shelf. Channels are initiated by irregularities in the ice along the grounding line and then enlarged by ocean melting. To a first approximation, spatially variable basal melting seaward of the grounding line acts as a steel-rule die or a stencil, imparting a channelized form to the ice base as it passes by. Ocean circulation in the region of high melt is inertial in the along-channel direction and geostrophically balanced in the transverse direction. Melt rates depend on the wavelength of imposed variations in ice thickness where it enters the shelf, with shorter wavelengths reducing overall melting. Petermann Gletscher’s narrow basal channels may therefore act to preserve the ice shelf against excessive melting. Overall melting in the model increases for a warming of the subsurface water. The same sensitivity holds for very slight cooling, but for cooling of a few tenths of a degree a reorganization of the spatial pattern of melting leads, surprisingly, to catastrophic thinning of the ice shelf 12 km from the grounding line. Subglacial discharge of fresh water along the grounding line increases overall melting. The eventual steady state depends on when discharge is initiated in the transient history of the ice, showing that multiple steady states of the coupled system exist in general.

Inverse problems are used to estimate model parameters from observations. Many inverse problems are ill-posed because they lack stability, meaning it is not possible to find solutions that are stable with respect to small changes in input data. Regularization techniques are necessary to stabilize the problem. For nonlinear inverse problems, iterative inverse methods can be used as a regularization method. These methods start with an initial estimate of the model parameters, update the parameters to match observation in an iterative process that adjusts large-scale spatial features first, and use a stopping criterion to prevent the overfitting of data. This criterion determines the smoothness of the solution and thus the degree of regularization. Here, iterative inverse methods are implemented for the specific problem of reconstructing basal stickiness of an ice sheet by using the shallow-shelf approximation as a forward model and synthetically derived surface velocities as input data. The incomplete Gauss-Newton (IGN) method is introduced and compared to the commonly used steepest descent and nonlinear conjugate gradient methods. Two different stopping criteria, the discrepancy principle and a recent- improvement threshold, are compared. The IGN method is favored because it is rapidly converging, and it incorporates the discrepancy principle, which leads to optimally resolved solutions.