Future fluctuations of Mer de Glace, French Alps, assessed using a parameterized model calibrated with past thickness changes

Abstract

Simulations of glacier evolution are needed to assess future changes in the runoff regime of mountain catchments. A simplified parameterized model is applied here to simulate future thickness changes and glacier retreat of Mer de Glace, French Alps. A normalized thickness change function describing the spatial distribution of surface-elevation changes as a function of elevation has been determined. The model reveals that under present climatic conditions Mer de Glace will continue to shrink dramatically in the coming decades, retreating by 1200 m between now and 2040. The method has certain limitations related to the uncertainties of the normalized function based on thickness change data. An error of 10% in the normalized function leads to uncertainties of 46%, 30% and 18% in Mer de Glace front, surface area and glacier-wide mass-balance changes respectively in 2040. Because the difference of the normalized function largely exceeds 10% from one glacier to another, even within a given glacier size class and elevation range, it would be very risky to extrapolate the normalized function to unmeasured glaciers. Consequently, the method is applicable only on glaciers where past surface elevation changes are well constrained.

1. Introduction

Climate change has already had a great impact on glaciers in the Alps (e.g. Vincent and others, 2004; Lemke and others, 2007), and the expected retreat of glaciers in the future will significantly affect the runoff regime of mountain catchments (Huss and others, 2008a; Casassa and others, 2009; Gabbi and others, 2012; Salzmann and others, 2012). Ice flow models can be used to simulate glacier fluctuations in the future from different mass-balance scenarios (Wallinga and Van de Wal, 1998; Le Meur and others, 2007). Numerous ice-flow model studies with different degrees of complexity have been carried out on glaciers (e.g. Bindschadler, 1982; Hubbard and others, 1998; Span and Kuhn, 2003; Leysinger Vieli and Gudmundsson, 2004; Oerlemans 2007; Schäfer and Le Meur, 2007; Adhikari and Marshall, 2013). The fluctuations of alpine glaciers have been assessed in many glaciological studies using simple flowline models (e.g. Greuell, 1992; Oerlemans and others, 1998; Vincent and others, 2000; Sugiyama and others, 2007; Adhikari and Huybrechts, 2009). However, these models do not include the overall three-dimensional (3-D) mechanical processes and therefore cannot accurately simulate glacier fluctuations. On the other hand, complex 3-D ice flow models solve the nonlinear Stokes equations but require extensive bedrock topography data and considerable computational resources (e.g. Gudmundsson, 1999; Le Meur and others, 2007; Jouvet and others, 2008). Moreover, dynamic changes are also driven by the subglacial drainage system that influences basal motion. It is very difficult to take these complex subglacial hydrological processes into account in numerical ice-flow models (Mair and others, 2003; Bartholomaus and others, 2008; Pimentel and Flowers, 2011). This is probably the main cause of large discrepancies between reconstructed and measured long-term temporal changes in ice flow velocities and discrepancies relative to the response time of snout fluctuations (Le Meur and Vincent, 2003). For these reasons, conceptual approaches can be useful (Huss and others, 2008a; Luthi and others, 2010).

The A h parameterization proposed by Huss and others (2008a) can be used to simulate surface elevation changes using a function relating the surface elevation change A h occurring over a given time interval to the elevation of the glacier surface h. The parameterization takes into account the conservation of mass and is thus well suited for transient hydrological modeling of ice storage changes. In addition, this methodology is suitable for glaciers for which the upper bedrock topography remains unknown, as is the case for Mer de Glace. This method has been used to simulate glacier retreat and assess both the future high-mountain hydrology and the contribution of glacier storage change to runoff from large drainage basins (Huss and others, 2010; Huss, 2011; Giesen and Oerlemans, 2013; Salzmann and others, 2013). The parameterization has been shown to yield results for future glacier geometry change that agree well with 3-D finite-element ice flow modeling for two Alpine glaciers (Huss and others, 2010).

The objective of our study is to test this method on a well-measured mountain glacier and to analyze its performance and limitations. For this, we apply the method to Mer de Glace for which numerous measurements of thickness changes of the tongue are available going back to the beginning of the 20th century. In this way, the empirical function relating surface elevation change to elevation is well constrained. We simulate glacier fluctuations over the coming decades using different surface mass-balance scenarios, carefully analyzing the uncertainties. We also analyze the sensitivity of simulations to the uncertainties related to the parameterization and compare it to the sensitivity of the model to climatic scenarios.

2. Study Site and Data

Mer de Glace (458550 N, 68570 E), the largest glacier in the French Alps, with an area of -32 km², is located in the Mont Blanc massif (Fig. 1). The maximum elevation of its upper accumulation area reaches ~4300ma.s.l. From this accumulation region, the ice flows rapidly through a narrow, steep portion (an icefall between 2700 and 2400 m) before feeding the lower 7 km of the glacier down to a front located at ~1500m. Mer de Glace includes several tributaries. Numerous geodetic and glaciological measurements have been performed on this glacier. The first topographic measurements were performed on the ablation area of Mer de Glace by Joseph Vallot at the end of the 19th century (Vallot, 1905; Reynaud, 1973; Lliboutry and Reynaud, 1981; Nussbaumer and others, 2007). From 1910 to 1960, ‘Les Eaux et Forets’ institute, which was in charge of forest and hydrology measurements (Mougin, 1933), measured the surface elevations of three cross sections. They added another cross section named ‘Trelaporte’ in 1922. In the 1960s-70s, the hydroelectric power company Electricite de France started to collect subglacial water and measured thickness fluctuations on the same cross sections. A new cross section, named ‘Tacul’, was added in 1967. Since 1985, regular thickness variation measurements have been carried out by the Laboratoire de Glaciologie et Geophysique de l’Environnement (LGGE), Grenoble, along the same five cross sections (Fig. 1) using theodolite or differential GPS methods with an uncertainty of less than ±0.20ma–¹.

Fig. 1. Map of Mer de Glace. Locations of cross sections for calculated thickness changes.

The bedrock topography was determined below 2300 ma.s.l. using mechanical borehole drillings and seismic soundings (Süsstrunk, 1951; Vallon, 1961, 1967; Gluck, 1967).

In addition, annual surface mass balances were measured at 10–30 stakes in the ablation area between 1979 and 1993. However, these measurements were only carried out systematically at the end of the ablation period from 1990 to 1993. Since 1993, systematic winter and summer mass-balance measurements (May and September respectively) have been performed at 30 sites covering the entire glacier surface (Vincent, 2002). Moreover, total cumulative mass balances have been calculated using an old map with elevation contours and repeated photogrammetric measurements. This old map of Mer de Glace was surveyed between 1900 and 1905 by Henri, Joseph and Charles Vallot at 1 :20000 scale. The accuracy of this map is unknown but, from the comparison of the elevations of points located outside glaciated areas, it can be assumed that coordinates are accurate to within several meters. Photogrammetric measurements were performed using aerial photographs by IGN (French National Geographic Institute) in 1958 and by LGGE in 2003 and 2008. Photogrammetric measurements are accurate to within several meters (vertical and horizontal accuracy) for aerial photographs taken prior to 1970 and down to 1 m for photographs taken later.

In addition, results from satellite-derived digital elevation models (DEMs) have been used in this study (Berthier and others, 2004; Berthier and Vincent, 2012). DEMs are calculated from a pair of 10m resolution SPOT (Satellite Pour l’Observation de la Terre) satellite images acquired in 1994. Moreover, a DEM for 1979, made by IGN using aerial photographs, was used to estimate the elevation changes of the tongue between 1979 and 1994. The elevation accuracy of the map of 1979 was assessed by comparison with field topographic measurements (Berthier, 2005). The difference can reach 10m locally but is generally <5 m. Comparison with elevation profiles surveyed in the field each year showed, after Gaussian filtering and averaging over elevation ranges, that the DEM-derived elevation changes from satellite images acquired in 1994 and aerial photographs of 1979 were accurate to within ±5 m in the ablation area (Berthier and others, 2004; Berthier, 2005).

3. Methodology: ΔH Parameterization

The parameterization consists in determining the spatial distribution of thickness changes in response to a change in mass balance (Huss and others, 2008a, 2010). First of all, the thickness changes are calculated from mean elevations across cross sections determined from DEMs acquired in the past. For this purpose, for each available DEM (1905, 1958, 2003 and 2008) we extracted the mean surface elevation at 50 fixed cross sections (Fig. 1) covering all altitude ranges between the head and the terminus of the glacier. Mean elevation differences were then computed at those cross sections for 1905-2008, 1958-2008 and 2003-08. An alternative method would be to calculate the averaged DEM-derived elevation changes within a given elevation range (Berthier and others, 2004; Huss and others, 2008b). However, elevation changes show a high spatial variability, even within the same elevation range, for several reasons. First, some areas were not measured accurately in the past because the terrain was inaccessible, as is the case for field measurements performed for the old maps. Second, the average elevations of rough terrain with crevasses are not measured accurately whatever the method used. For these reasons, we thought it was preferable to determine the elevation changes on selected well-measured cross sections, with exactly the same locations.

In this way, we obtained a function relating the surface elevation change to the elevation of the glacier surface. The surface elevations z are normalized with respect to the elevation range according to (z _max-z)/(z _max-z _min) where z _max and z _min are maximum and minimum surface elevations (Huss and others, 2010). Given that z _min changes with time, it is calculated from the average of the two measurement periods. The thickness changes are normalized with respect to the maximum elevation change Ah _max as Ah/Ahmax. The parameterization is approximated using a polynomial function. Secondly, this function was used in an opposite way to calculate the distribution of thickness changes over the whole surface of the glacier. Forced by the glacier-wide mass balance, the parameterization makes it possible to calculate the thickness changes for each elevation range and, consequently, the snout fluctuations. The integrated elevation changes over the whole area should be equal to the glacier-wide mass balance. We assume that the mass loss is glacier ice (density 900 kg m–³).

The glacier-wide mass balance is assessed for the future using three climatic scenarios and the sensitivity of surface mass balance to air temperatures. In the first scenario, the future distribution of the surface mass balances with elevation is assumed to be equal to the average mass balance observed over the last two decades (1992-2012). In the second scenario, we assume atmospheric warming of 0.02°Ca^–1. In a third scenario, we assume an atmospheric warming of 0.04°Ca–¹. Precipitation is assumed to be constant with time. The purpose here is not to provide accurate mass-balance values for climatic scenarios but to assess the ability of the parameterization to simulate the retreat of the glacier.

4. Results

4.1. Glacier-wide mass balance since the beginning of the 20th century

Cumulative mean specific net balances are shown in Figure 2a. Prior to direct measurements, the annual mass balances were reconstructed using hypsometric and meteorological data (Vincent, 2002). The long-term mass balance between 1905 and 2008 is constrained by glacier volume variations deduced from a map (1905) and aerial photographs (in 1958, 2003 and 2008) whereas the interannual variability is deduced from annual glaciological measurements or annual reconstructed data. As shown by Vincent (2002), it would be illusive to reconstruct cumulative mass balances in the past without topographic or photogrammetric maps. For the purpose of our study, the goal is not to test the ability of the model to reconstruct a mass balance, but rather to obtain realistic results. In this way, reconstructed and observed glaciological mass balances have been adjusted so that cumulative mass balances match volumetric mass balances from geodetic measurements (Vincent, 2002). The cumulative mass balances for the periods 1905–2003 and 1905–2008 have been estimated at –22.9 and –28.8 m w.e. respectively from geodetic data. All annual mass-balance measurements available between hydrological years 1989/90 and 2011/12 were processed using the linear mass-balance model (Lliboutry, 1974; Thibert and Vincent, 2009) to obtain centered mass balances. These centered mass balances represent the annual deviation from the mean. The glaciological mass balances calculated between 1989 and 2012 were adjusted over the 2003–08 period using photogrammetric results whereas the reconstructed mass balances were adjusted over the 1905–99 period.

Fig. 2. (a) Cumulative mass balance of Mer de Glace. (b) Thickness changes observed on five cross sections of the tongue of Mer de Glace.

The cumulative mean specific balance of Mer de Glace between 1905 and 2008 is negative (loss of about 28.8 mw.e.) but shows strong fluctuations. Note that the glacier net balance increased between 1960 and the beginning of the 1980s as elsewhere in the Alps (Huss, 2008b). After 1987, the mass balance was consistently negative except for 1995 (glaciological year 1994/95). In addition, the cumulative mass balance reveals that the loss of ice has been more pronounced since 2001 (Fig. 2a). lost 207 m from the end of the 19th century until the ice disappeared completely in 2008. Those values are in agreement with similar earlier estimates by Nussbaumer and others (2007).

4.2. Observed thickness changes

Comparison of thickness changes of the five cross sections located at different elevations of the tongue reveals significant differences (Fig. 2b). The tongue of the glacier is located in a narrow gorge, and changes in ice flux lead to large changes in ice thickness and ice flow velocity. In addition, as for other glaciers, the thickness changes of the tongue are driven by the changes in compressing flow (Vincent and others, 2009). Between 1990 and 2012, ‘Tacul’ cross section (2200 m a.s.l.) lost 58m in thickness while the thickness decreased by 77 m at ‘Montenvers’ cross section (1700 m a.s.l.). Over longer periods, the differences in thickness changes with elevation are even more pronounced. For example, ‘Trelaporte’ and ‘Mottets’ cross sections experienced losses of 102 and 201 m respectively between 1922 and 2008, whereas the altitude difference is 500 m between these two cross sections (Fig. 2). The lower cross section (‘Mottets’), located in the vicinity of the snout,

For the 50 selected cross sections, we calculated the thickness changes for the periods 1905-2008, 1958-2008 and 2003-08. In Figure 3a, the thickness changes are plotted against elevation. In addition, the thickness changes for the periods 1994-2008 and 1979-94, obtained in the ablation zone only (Berthier and Vincent, 2012), are plotted, although these data do not cover the entire surface area of the glacier. The thickness changes derived from the DEM comparisons reveal a regular relationship between elevation changes and elevation (Fig. 3a). As shown by other studies (e.g. Schwitter and Raymond, 1993; Arendt and others, 2002; Span and Kuhn, 2003; Vincent and others, 2009; Huss and others, 2010), the thickness changes are smaller in the upper part of the glacier and larger near the terminus.

Fig. 3. (a) Thickness changes observed on 50 selected cross sections distributed over the whole surface area of the glacier for the periods 1905–2008, 1958–2008, 1979–94, 1994–2008 and 2003–08. (b) Normalized thickness variations against normalized elevation.

4.3. Normalized thickness changes

Thickness change data during 1905-2008, 1958-2008 and 2003-08 have been normalized in order to obtain the dimensionless A h parameterization (Fig. 3b). Thickness change data for the periods 1979-94 and 1994-2008 were not used given that these data did not cover the entire surface area. Normalized thickness changes for 2003-08 show significant discrepancies in the lower part of the glacier compared to 1905-2008 and 1958-2008. These discrepancies could be due to the short duration of the period (5 years) which could be insufficient to perform such an analysis. The duration of the period required for the parameterization is discussed below (Section 5). Averaged normalized thickness changes have been calculated from a sixth-order polynomial function fitted to the 1905-2008 and 1958-2008 data in order to obtain a smooth A h function (thick black line in Fig. 3b) and to reduce noise on elevation changes. Another option could be to select the period 1905-58 instead of 1905-2008. However, the uncertainties of the 1905 and 1958 maps are large and the signal of elevation change is small during this period, leading to a signal-to-noise ratio too small to construct a relevant normalized function (Section 5). Note that individual values do not differ by more than 10% from the polynomial fit except for data from the 2003-08 period. The standard deviation between averaged and measured normalized elevation changes is 0.07. Some discrepancies shown for the period 1905-2008 could be due to the inaccuracy of the geodetic measurements of 1905 which were performed using a theodolite and which may not have covered the whole glacier surface, in particular areas with crevasses.

Normalized data from the 1994-2008 period have been added in Figure 3b, assuming no elevation change at the head of the glacier during this period. This dataset, restricted to the lower part of the glacier, is assumed to include the maximum elevation change. Those measurements are in agreement with the polynomial fit. Using the same assumption, we calculated normalized elevation changes for the 1979-94 dataset. These results reveal large discrepancies that are discussed below.

4.4. Glacier retreat based on the parameterized model

First, the parameterized model was tested for the period 1958-2012, for which numerous observations are available. For this purpose, our modeling started with the surface topography from 1958. The model was forced by the reconstructed and observed annual glacier-wide mass balance (Fig. 2a). The yearly thickness changes of each DEM cell were calculated using the elevation change distribution obtained from A h parameterization. The A h parameterization is applied to update the 3-D glacier geometry at the end of each hydrological year. If the thickness change is larger than the total glacier thickness, the elevation of the surface is taken to be the bedrock elevation.

The modeling results are shown in Figure 4 and compared to the observed terminus positions for 2003, 2008 and 2012. These results show good agreement with the measurements. Note that the thickness changes calculated at a yearly scale are not relevant given the delayed response of the glacier geometry to mass-balance fluctuations (Jóhannesson and others, 1989). The appropriate timescale on which the calculated thickness changes and front retreat are relevant is dealt with in Section 5.

Fig. 4. Simulated thickness and snout fluctuations in (a) 2003, (b) 2008 and (c) 2012 (simulations start in 1958). Dashed lines correspond to measurements. The triangle corresponds to the snout position measured in 2012.

Numerical simulations were then carried out to simulate the terminus positions in the future using the three mass-balance scenarios described previously. The first scenario, for which the future distribution of the surface mass balances with elevation is assumed to be equal to the average mass balance observed over the last two decades (1992-2012) corresponds to a glacier-wide mass balance of -1.05mw.ea–¹ for the first year of the simulation in 2013. For the following years, the change in glacier-wide mass balance is mainly influenced by the change in surface area and the surface lowering of the tongue (Huss and others, 2012). The glacier-wide mass balance tends to be less negative due to the decrease in the surface area at low elevation while the surface lowering tends to increase ablation. The first effect is dominant and the simulated glacier-wide mass balance is -0.90mw.ea–¹ in 2040 with the same distribution of surface mass balance with elevation.

For the second scenario (atmospheric warming of 0.02°Ca–¹) and the third scenario (warming of 0.04°Ca–¹), the surface mass balance was calculated for each elevation range using the mass-balance sensitivity to temperature obtained by Vincent (2002). Thanks to numerous measurements, regression coefficients have been calculated between ablation and cumulative positive degree-day (CPDD) as a function of altitude. In order to estimate ablation sensitivity to temperature variations over the whole summer period, each degree-day factor has been multiplied by the mean number of days for which temperature is higher than 0°C at the observation elevation (Vincent, 2002).

The surface mass balance could be calculated more accurately using a degree-day model (Vincent, 2002; Huss and others, 2008a), but this level of sophistication is beyond the scope of this study given that the climatic scenarios used here are arbitrary.

Bedrock topography data are required to update glacier surface area. Unfortunately, bedrock topography is unmeasured in the upper part of the glacier, above 2300 ma.s.l. The thickness distribution could be inferred from the methodology developed in Huss and Farinotti (2012). This approach is very useful to derive ice fluxes. However, it does not allow us to infer ice thicknesses accurately, especially for the edges of the glacier. It would be risky to determine the surface area changes from this method.

Consequently, in this zone, the surface area of the glacier is assumed to be constant and, for example, the growth of rock outcrops has to be neglected although it has been observed in the Alps (Paul and others, 2007).

The simulated terminus fluctuations are plotted in Figure 5 and 6. With the first scenario (constant climatic conditions), the terminus retreats 1200 m between 2012 and 2040. The second and third scenarios show a retreat of 1295 m and 1375 m respectively over the same period.

Fig. 5. Simulated thickness and snout fluctuations in (a) 2020, (b) 2030 and (c) 2040 (simulations start in 2008) with constant climatic conditions and a temperature increase of 0.02°C a^–1 (yellow line) and 0.04°C a^–1 (red line).

Fig. 6. Observed length fluctuations since 1874 (black dots) and simulated changes over the coming decades using the parameterized model with three different climatic scenarios: constant climatic conditions (blue dashed line), a warming of 0.02°C a^–1 (yellow dashed line) and a warming of 0.04°C a^–1 (red dashed line).

5. Discussion

Our assumption that the distribution of normalized elevation changes is constant with time can be subject to question. Indeed, several limitations of the method should be considered. First, due to the delayed response of glacier geometry, the redistribution of surface mass-balance anomalies by ice flow is not instantaneous (Lliboutry and Reynaud, 1981). Consequently, the method cannot be applied over short periods. Normalized thickness changes obtained for the period 2003-08 (Fig. 3b) seem to confirm this conclusion given that they show large departure from the polynomial fit. The difference with the average normalized function exceeds 20%. This tends to show that a 5 year period is not sufficient to derive the A h parameterization.

A second limitation, as mentioned by Huss and others (2010), is that the response of glacier geometry could be different for periods with positive mass balances and negative mass balances. To study this question, additional experiments were performed using data obtained for the period 1979–94 (Berthier and Vincent, 2012) during which the glacier was close to steady-state conditions. During this period, the average glacier-wide mass balance was –0.2m w.e. a^–1 and the observed thickness changes did not exceed 15 m. Although these data are relative to the ablation zone only and do not cover the entire surface area of the glacier, we assumed that no elevation change occurred at the top of the glacier during this period and calculated a normalized function from these data. The results, not shown in this paper, reveal a large departure of >50% from the average normalized function calibrated over the 1905–2008 and 1958–2008 periods. Note that this difficulty in reconstructing a normalized function with an uncertainty lower than 10% during 1979–94 could be related to the small values of observed thickness changes over this period. Similarly, for short periods of data, small thickness changes make it impossible to determine a relevant normalized function.

A third limitation is that the normalized elevation change function could be affected over time by the retreat of the glacier and the large changes of glacier geometry which influence the ice flow and consequently the thickness changes. However, the results obtained from the 1958–2008 and 1905–2008 periods are similar within an uncertainty of <10%.

Assuming that the parameterized relationship is stable over a long period, an important question is to quantify the influence of parameterization uncertainty on simulated thickness changes and on future glacier retreat. We carried out some numerical experiments to study this question. As mentioned previously, the standard deviation of the difference with the polynomial fit obtained for the periods 1958-2008 and 1905-2008 is 0.07. Measurements from the 1994-2008 period show that the differences between normalized data and the polynomial function can reach 10%. Consequently we can assume it is difficult to obtain a parameterized function with an uncertainty less than 10%. The previous polynomial function was shifted by this amount to make new numerical simulations and to assess the impacts on glacier retreat. For this purpose, 0.1 was added to or removed from the normalized function except if the obtained value exceeds 1.0 or is negative.

The results (Fig. 5) reveal that the snout fluctuations have an uncertainty of ±430 m in 2030 and ±550 m in 2040. Previous numerical experiments with different climatic scenarios showed that a 0.04°Ca–¹ temperature increase led to an additional retreat of 175 m in 2040. From these results, we can conclude that the bias on the snout changes caused by the uncertainty of 10% related to the normalized function can reach about three times the influence related to a 0.04°Ca–¹ difference in climatic scenarios.

In summary, the uncertainty of 10% related to the normalized function causes an uncertainty of 46% and 30% on the snout fluctuations and the surface-area changes respectively. It also leads to an uncertainty of 18% on the glacier-wide mass-balance changes, and consequently on the glacier contribution to river runoff, in 2040. The annual runoff, measured at the gauge station ‘Le Bois du Bouchet’ located -2.6 km downstream of the glacier snout, reveals a mean annual runoff of 5.21 m³ s–¹ for a total drainage basin of 78 km². Given that the glacier contribution is close to 1.1 m³ s^–1, the uncertainty of 10% related to the normalized function leads to an uncertainty of 3.8% of this runoff.

Another interesting question concerns the thickness change accuracy required to establish a relevant normalized function. Normalized thickness change is obtained from Δh/Δh _max . Considering that errors are random, the error a_norm on the estimate of the normalized thickness change can be expressed as

(1)

Given that Δh<Δh_max and σ _Δh is similar to a_Ahmax, this gives

(2)

Therefore, to obtain an accuracy better than 10% for the normalized function, the uncertainty of the measured thickness change should be less than 12.2 and 7.8 m for maximum elevation changes of 173 and 110 m respectively, corresponding to the maximum elevation changes of the 1905-2008 and 1958-2008 periods.

This explains why the method is not applicable to the 1979-94 period. Given that the maximum thickness change during this period is 15 m, the uncertainty of measured thickness changes should be less than 1.1 m, which is much smaller than the uncertainty of actual measurements (±5 m). A similar conclusion can be drawn from the 2003-08 measurements. Indeed, although the accuracy of elevations coming from photogrammetric measurements obtained in 2003 and 2008 is <1 m, the measurement period is probably too short to obtain a function that is representative for long periods. As seen in Figure 3a, the averaged differences between surface-elevation change at cross sections and the averaged polynomial function can exceed 5 m.

These calculations show some limitations of the method. The method should therefore be used with caution. It cannot be extrapolated to unmeasured glaciers for which the uncertainty of the parameterized function largely exceeds 10%, even for the same glacier size class within the same mountain range (Huss and others, 2010).

6. Conclusions

A simplified parameterized model offers an attractive alternative to ice flow models which require extensive data relative to bedrock topography. In this paper, the A h parameterization was used to simulate the future thickness changes and glacier retreat of Mer de Glace. This method is suitable for the study site because the bedrock topography remains unknown in the accumulation zone. Numerous elevation-change data coming from measurements performed over the past century allowed us to establish a reliable normalized function of thickness change. The method is mass-conserving and suitable for long-term simulations of glacier fluctuations.

Under present climatic conditions, we found that the terminus would retreat by 1200 m between 2012 and 2040.

We have pointed out the limitations of the method. The results obtained are relevant only if the normalized function is based on numerous elevation-change data obtained over the whole surface area of the glacier. The method requires data over long periods and large elevation changes. Moreover, the uncertainty of elevation-change measurements should be small compared to the observed elevation changes. We found that the method is not relevant for short periods or periods with positive mass balances. For example, we conclude that the method is not applicable to the 1979-94 period, either because the glacier-wide mass-balance was close to zero (-0.2 mw.e.a–¹) or because the thickness changes were small (<15m) compared to the uncertainties. Moreover, the method is not really applicable to positive mass balances because it does not allow the transfer of positive elevation changes close to the snout to a snout advance. Indeed, in this case, an assumption concerning the slope of the snout is required.

From our numerous data obtained for Mer de Glace, we concluded that the best estimation of the uncertainty of the normalized function is ~10%. We calculated that this error leads to an uncertainty of 550 m on the Mer de Glace front fluctuation after 30 years of simulation. Given that the spread of the normalized function largely exceeds 10% from one glacier to another, even within a given glacier size class and a given elevation range (Huss and others, 2010), it is not reasonable to extrapolate the normalized function to unmeasured glaciers in order to estimate snout fluctuations or surface-area changes.

In summary, the uncertainty of 10% related to the normalized function causes a large uncertainty in the snout fluctuations and the surface-area changes. It leads to an uncertainty of 18% on the glacier-wide mass-balance changes, and consequently in the glacier contribution to river runoff, in 2040.

We conclude the method is applicable only when numerous elevation-change data are available for long periods and for negative mass balances in order to keep the normalized function uncertainty below 10%.