1. Introduction

Climate change has already had a great impact on glaciers in the Alps (e.g. Reference Vincent, Kappenberger, Valla, Bauder, Funk and Le MeurVincent and others, 2004; Reference LemkeLemke and others, 2007), and the expected retreat of glaciers in the future will significantly affect the runoff regime of mountain catchments (Reference Huss, Farinotti, Bauder and FunkHuss and others, 2008a; Reference Casassa, López, Pouyard and EscobarCasassa and others, 2009; Reference Gabbi, Farinotti, Bauder and MaurerGabbi and others, 2012; Reference Salzmann, Machguth and LinsbauerSalzmann and others, 2012). Ice flow models can be used to simulate glacier fluctuations in the future from different mass-balance scenarios (Reference Wallinga and Van de WalWallinga and Van de Wal, 1998; Reference Le Meur, Gerbaux, Schäfer and VincentLe Meur and others, 2007). Numerous ice-flow model studies with different degrees of complexity have been carried out on glaciers (e.g. Reference BindschadlerBindschadler, 1982; Reference Hubbard, Blatter, Nienow, Mair and HubbardHubbard and others, 1998; Reference Span and KuhnSpan and Kuhn, 2003; Reference Leysinger Vieli and GudmundssonLeysinger Vieli and Gudmundsson, 2004; Reference OerlemansOerlemans 2007; Reference Schäfer and E.Schäfer and Le Meur, 2007; Reference Adhikari and MarshallAdhikari and Marshall, 2013). The fluctuations of alpine glaciers have been assessed in many glaciological studies using simple flowline models (e.g. Reference GreuellGreuell, 1992; Reference OerlemansOerlemans and others, 1998; Reference Vincent, Vallon, Reynaud and Le MeurVincent and others, 2000; Reference Sugiyama, Bauder, Zahno and FunkSugiyama and others, 2007; Reference Adhikari and HuybrechtsAdhikari 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. Reference GudmundssonGudmundsson, 1999; Reference Le Meur, Gerbaux, Schäfer and VincentLe Meur and others, 2007; Reference Jouvet, Picasso, Rappaz and BlatterJouvet 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 (Reference Mair, Willis, Fischer, Hubbard, Nienow and HubbardMair and others, 2003; Reference Bartholomaus, Anderson and AndersonBartholomaus and others, 2008; Reference Pimentel and FlowersPimentel 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 (Reference Huss, Farinotti, Bauder and FunkHuss and others, 2008a; Luthi and others, 2010).

The A h parameterization proposed by Reference Huss, Farinotti, Bauder and FunkHuss 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 (Reference Huss, Jouvet, Farinotti and BauderHuss and others, 2010; Reference HussHuss, 2011; Reference Giesen and OerlemansGiesen and Oerlemans, 2013; Reference SalzmannSalzmann 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 (Reference Huss, Jouvet, Farinotti and BauderHuss 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 (Reference VallotVallot, 1905; Reference ReynaudReynaud, 1973; Reference Lliboutry and ReynaudLliboutry and Reynaud, 1981; Reference Nussbaumer, Zumbühl and SteinerNussbaumer and others, 2007). From 1910 to 1960, ‘Les Eaux et Forets’ institute, which was in charge of forest and hydrology measurements (Reference MouginMougin, 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 (Reference SüsstrunkSüsstrunk, 1951; Reference VallonVallon, 1961, Reference Vallon1967; Reference GluckGluck, 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 (Reference Berthier, Arnaud, Baratoux, Vincent and RémyBerthier and others, 2004; Reference Berthier and VincentBerthier 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 (Reference Berthier, Arnaud, Baratoux, Vincent and RémyBerthier and others, 2004; Reference BerthierBerthier, 2005).

3. Methodology: ΔH Parameterization

The parameterization consists in determining the spatial distribution of thickness changes in response to a change in mass balance (Reference Huss, Farinotti, Bauder and FunkHuss and others, 2008a, Reference Huss, Jouvet, Farinotti and Bauder2010). 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 (Reference Berthier, Arnaud, Baratoux, Vincent and RémyBerthier and others, 2004; Reference Huss, Bauder, Funk and HockHuss 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 (Reference Huss, Jouvet, Farinotti and BauderHuss 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 (Reference VincentVincent, 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 Reference VincentVincent (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 (Reference VincentVincent, 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; Reference Thibert and VincentThibert 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 Reference Nussbaumer, Zumbühl and SteinerNussbaumer 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 (Reference Vincent, Soruco, Six and Le MeurVincent 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 (Reference Berthier and VincentBerthier 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. Reference Schwitter and RaymondSchwitter and Raymond, 1993; Reference Arendt, Echelmeyer, Harrison, Lingle and ValentineArendt and others, 2002; Reference Span and KuhnSpan and Kuhn, 2003; Reference Vincent, Soruco, Six and Le MeurVincent and others, 2009; Reference Huss, Jouvet, Farinotti and BauderHuss 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.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 (Reference Jóhannesson, Raymond and WaddingtonJó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 (Reference Huss, Hock, Bauder and FunkHuss 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; Reference Huss, Farinotti, Bauder and FunkHuss 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 Reference Huss, Hock, Bauder and FunkHuss 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 (Reference Paul, Kääb and HaeberliPaul 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 (Reference Lliboutry and ReynaudLliboutry 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 Reference Huss, Jouvet, Farinotti and BauderHuss 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 (Reference Berthier and VincentBerthier 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 (Reference Huss, Jouvet, Farinotti and BauderHuss 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 (Reference Huss, Jouvet, Farinotti and BauderHuss 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%.