Modelling rock wall permafrost degradation in the Mont Blanc massif from the LIA to the end of the 21 century

High alpine rock wall permafrost is extremely sensitive to climate change. Its degradation has strong impact on landscape evolution and can trigger rock falls constituting an increasing threat to socio-economical activities of highly frequented areas; quantitative understanding of permafrost evolution is crucial for such communities. This study investigates the long-term evolution of permafrost in three vertical cross-sections of rock wall sites between 3160 and 4300 m a.s.l. in the 15 Mont Blanc massif, since the Little ice Age (LIA) steady-state conditions to 2100. Simulations are forced with air temperature time series, including two contrasted air temperature scenarios for the 21 st century representing possible lower and upper boundaries of future climate change according to the most recent models and climate change scenarios. The 2D Finite Element model accounts for heat conduction and latent heat transfers, and the outputs for the current period (20102015) are evaluated against borehole temperature measurements and an electrical resistivity transect: permafrost conditions 20 are remarkably well represented. Through the past two decades, permafrost has disappeared on faces with a southerly aspect up to 3300 m a.s.l., and possibly higher. Warm permafrost (i.e. > -2°C) has extended up to 3300 and 3850 m a.s.l. in N and S-exposed faces, respectively. Along the 21 st century, warm permafrost is likely to extent at least up to 4300 m a.s.l. on Sexposed rock walls, and up to 3850 m a.s.l. at depth on the N-exposed faces. In the most pessimistic case, permafrost will disappear at depth on the S-exposed rock walls up to 4300 m a.s.l., whereas warm permafrost will extent at depth of the N 25 faces up to 3850 m a.s.l., but possibly disappearing at such elevation under the influence of a close S face. The results are site-specific and extrapolation to other sites is limited by the imbrication of the local topographical and transient effects.


Introduction
The IPCC Fifth Assessment Report (AR5) has drawn a global increase in permafrost temperature since the 1980s (IPCC, 2014).By the end of the 21st century, the near-surface permafrost area is projected to retreat by 37 or 81 % according to RCP 2.6 and RCP 8.5 respectively (Representation Concentration Pathways with a projected increase in radiative forcing of 2.6 W m −2 since 1750; Vuuren et al., 2011).Concerns about natural disasters resulting from mountain permafrost degradation have started to rise during the late 1990s (IPCC, 1996).Haeberli et al. (1997) identified various types of high mountain slope instabilities that could be prepared or triggered by interactive processes between bedrock, permafrost and glaciers.Such processed have then been largely observed, especially with the increase in rockfall activity of high-elevation permafrost rock walls during the past two decades (Ravanel and Deline, 2011).
Since the hot summer of 2003 and the remarkable number of rockfalls observed in the European Alps (Schiermeier, 2003;Ravanel et al., 2011), rock wall permafrost has been intensively studied in various mountain areas (Gruber, 2005;Noetzli, 2008;Allen et al., 2009;Hasler, 2011;Hipp, 2012;Magnin, 2015a).The role of permafrost degradation in rock wall stability is recognised more and more (e.g.Krautblatter et al., 2013;Gjermundsen et al., 2015), and mountain permafrost is of high concern for construction practices (Har-F.Magnin et al.: Modelling rock wall permafrost degradation in the Mont Blanc massif ris et al., 2001a;Bommer et al., 2010).The destabilisation of rock wall permafrost endangers high mountains activities, infrastructure (Duvillard et al., 2015), mountain climbers and workers.Valley floors could be affected by high mountain hazards owing to the possible cascading effects (Deline, 2001;Einhorn et al., 2015).The acceleration of rock wall retreat resulting from rapid permafrost degradation (Haeberli and Burn, 2002) has substantial implication for landscape evolution.Major changes are visible at human timescales, such as the sudden disappearance of the famous Bonatti rock pillar and its climbing routes in the Mont Blanc massif in 2005 (Ravanel and Deline, 2008).
Rock wall permafrost is highly sensitive to air temperature change because (i) it is directly coupled with the atmosphere (absence of debris and seasonal snow cover), (ii) the delaying effect of latent heat processes is reduced due to the low ice content (Smith and Riseborough, 1996), and (iii) it is subject to multi-directional warming from the different summit sides (Noetzli et al., 2007).Therefore, it is prone to much faster changes than any other kind of permafrost (Haeberli et al., 2010).
The monitoring of rock wall permafrost started in the late 1990s in Switzerland with the drilling of two boreholes at the Jungfraujoch site (PERMOS, 2004).A latitudinal transect along European mountains was later installed in the framework of the PACE project (Sollid et al., 2000;Harris et al., 2001bHarris et al., , 2009)).A warming trend clearly appeared over the past decade in most of the existing boreholes (Blunden and Arndt, 2014).
The presence of ice in the fractures of steep alpine bedrock has been demonstrated by engineering work (Keusen and Haeberli, 1983;King, 1996;Gruber and Haeberli, 2007).This ice contributes significantly to rock wall stability because it increases the tensile and shear strengths of the fractures (Davies et al., 2001;Krautblatter et al., 2013).The warming of an ice-filled fracture has two effects on its stability: the loss of bonding and the release of water which increases the hydrostatic pressure.An ice-filled fracture becomes critically unstable by between −1.4 and 0 • C (Davies et al., 2001).In this way, the warming of permafrost and the thickening of the active layer by heat conduction could be responsible for rock wall destabilisation (Gruber and Haeberli, 2007).However, heat advection through the circulation of water supplied by the melting of the interstitial ice, snow or glacier ice could warm permafrost at deeper layers than those reached by heat conduction (Hasler et al., 2011a).Hydraulic and hydrostatic pressures in frozen bedrock are modified under freezing and thawing and can be involved in rock wall destabilisation in a large range of processes (for a review of these processes see Matsuoka and Murton, 2008;Krautblatter et al., 2012).
Historical and recent rockfall events in the Mont Blanc massif have been systematically collated (Ravanel et al., 2010a;Ravanel and Deline, 2013).Their trend revealed a clear relationship with hot climate signals at various timescales from seasonal to decadal (Ravanel et al., 2010b;Ravanel and Deline, 2011;Huggel et al., 2012).In some cases, extreme precipitation events are thought to be the main triggering factor, which increase the hydraulic pressure in an impermeable bedrock permafrost system (Fischer et al., 2010).However, the role of extreme precipitation events in triggering rockfall is less obvious and systematic than that of extreme air temperature.Given recent evidence, one can assert that the magnitude and frequency of these hazards are likely to increase over the 21st century of projected global warming (IPCC, 2011).Knowledge on the current and future thermal state of the Mont Blanc massif rock walls is thus required in order to take into account a risk that threatens activities in this dense, highly frequented high mountain area.
Patterns and processes of long-term permafrost changes in steep mountain flanks were studied in idealised cases for the European Alps (Noetzli et al., 2007;Noetzli and Gruber, 2009) and Norway (Myhra et al., 2015).However, future changes in rock wall permafrost driven by the most recently released RCPs are yet to be addressed, and the sitespecific response to 21st century climate change has not been considered.Furthermore, evaluation of time-dependent rock wall permafrost models have remained limited by the lack of empirical data.To address site-specific long-term changes in rock wall permafrost of the Mont Blanc massif we ran a 2-D finite element model accounting for heat conduction and latent heat transfers on NW-SE cross sections of three sites covering an elevation transect, starting from 3160 up to 4300 m a.s.l., which encompasses currently warm and cold permafrost conditions.Transient simulations are run from the end of the LIA (ca.1850) to the end of the 21st century (2100) based on two different RCPs (4.5 and 8.5) accounting for moderately optimistic and pessimistic scenarios.Bidimensional models of the current period (2010)(2011)(2012)(2013)(2014)(2015) are benchmarked against an independent data set in order to evaluate the model performance.Even though changes in precipitation patterns (seasonality, frequency of extreme events and liquid/solid ratio) may play a marginal role in permafrost degradation pathways and rockfall triggering, only air temperature scenarios are considered since they constitute the dominant controlling factor of rock wall permafrost changes.The underlying research questions are as follows: -Is our modelling approach suitable for reproducing current permafrost conditions at the site scale?
-How has permafrost changed within these sites over the past decades?
-What is the possible evolution of rock wall permafrost by the end of the 21st century considering the latest IPCC projections?
This study provides as much insight into the recent changes of rock wall temperature as into its future evolution and is usable for retrospective analyses of rock wall instability as well as for assessing future hazards.

Study site and available data
The Mont Blanc massif is an external Variscan high mountain range culminating at 4809 m a.s.l., located on the western margin of the European Alps (Fig. 1).Its two major lithological units are a polymetamorphic basement along its western margin and a unit of Mont Blanc granite at its core (Bussy and von Raumer, 1994).It covers 550 km 2 over France, Switzerland and Italy, 30 % of which is glaciated (Gardent et al., 2014;Fig. 1).About 65 % of its rock walls above 2300 m a.s.l. are permanently frozen, according to a first estimation of permafrost distribution on the French side and borders (Magnin et al., 2015a;Fig. 1).For the purpose of this study, we selected three sites at various elevations and under various permafrost conditions: Aiguille du Midi, Grands Montets and Grand Pilier d'Angle.All three sites are located in the granitic area of the massif.Their elevations as well as their permafrost conditions are representative of the Mont Blanc massif rock walls.

Aiguille du Midi and bedrock temperature data
Studies on rock wall permafrost started at the end of 2005 in the Mont Blanc massif with the progressive installation of nine rock surface temperature (RST) sensors at the Aiguille du Midi summit (AdM), a set of three granite pillars.The AdM is accessible by cable car throughout the year (with approximately 500 000 visitors per year).As a pilot site in high-elevated permafrost research, the AdM is now equipped with a variety of instruments to measure rock wall temperature (Magnin et al., 2015c), snow cover (Magnin et al., 2017) and mechanics with extensometers (Ravanel et al., 2016).Three 10 m-deep boreholes of 15-node thermistor chains are installed in the AdM bedrock and have recorded temperature with 3 h time steps since December 2009 (NW and SE faces) and April 2010 (NE face).The AdM was selected because of the possibility to quantitatively evaluate the model outputs.It is characterised by the coexistence of cold permafrost (−4.5 • C at 10 m-depth) on its NW face and warm permafrost (−1.5 • C at 10 m depth) on its SE face (Fig. 2).Thermal effects of snow are observed in the three boreholes.The local cooling effect of a fracture has been detected at 2.5 m depth of the NW borehole.Nevertheless, temperature at 10 m-depth seems to be mainly governed by conductive heat transfer processes and lateral heat fluxes from the warm south face to the cold north face (Magnin et al., 2015c).

Grands Montets and ERT data
The Grands Montets (GM) is a summit culminating at 3296 m a.s.l., to the north and about 800 m below the Aiguille Verte (4122 m a.s.l.).In 1962-1963 a cable car was installed on its 60 • -steep north face to transport skiers up to the glaciated area.In May 2011, a RST logger was installed (GEOPrecision PT1000, sensor accuracy ±0.1 • C) at the foot of the highly fractured NW face (3058 m a.s.l.) in a 85 • -steep rock wall portion.It recorded the rock temperature at depths of 3, 10, 30 and 55 cm until January 2013.The 2012 mean annual rock surface temperature (MARST) at a depth of 3 cm was −1.4 • C. In 2012 and 2013, electrical resistivity tomography (ERT) soundings were conducted along the NW face of the GM and four other sites of the massif (Magnin et al., 2015d).The potential of ERT for qualitative evaluation of 2-D permafrost models has been demonstrated by Noetzli et al. (2008), as it covers a much wider and deeper rock wall portion than borehole.This makes ERT a potentially better approach with which to evaluate distributed models of rock wall permafrost because it has the capacity to represent the spatial variability of rock wall permafrost (Magnin et al., 2015d).Conversely, direct temperature measurements allow for quantitative evaluation, but have the disadvantage of being only representative for the measurement point.We selected the GM site because (i) a 160 m-long and 25 m-deep ERT transect is available for model evaluation, (ii) the site bears socio-economical interests with around 200 000 persons using the cable car every year, and (iii) it is located within the warm permafrost fringe of the massif as revealed by the RST data, permafrost map (Figs. 1 and 2) and the ERT transect.Moreover, this site has been regularly affected by rockfalls during the last decade which furthers the interest in studying its thermal dynamics.

Grand Pilier d'Angle
The third site was chosen based on its elevation in order to include an entirely cold permafrost site.Cold permafrost is likely to be present at the Grand Pilier d'Angle (GPA, 4304 m a.s.l.) on all the rock faces according to the permafrost map (Fig. 2).The east face of the GPA was strongly affected by a rock avalanche in November 1920.About 3 million m 3 of rock detached from the face in several stages and travelled onto the Brenva Glacier on a distance > 5 km, reaching the valley floor (Deline et al., 2015a).Because of its altitude, relief (900 m), steepness (subvertical rock walls) and remoteness, the GPA includes climbing routes among the most difficult and exposed of the Mont Blanc massif.For this last site, located on the Italian side of the massif, no data set is available for model evaluation, and the quality of the 2-D models will be assessed based on the evaluation of the two other sites.
3 Modelling approach

Background and strategy
Rock wall permafrost is mainly governed by air temperature and incoming solar radiations (Gruber et al., 2004).Local snow deposits can, however, either warm or cool the rock surface temperature compared to snow-free conditions, depending on the rock wall aspect and on the snow thickness www.the-cryosphere.net/11/1813/2017/The Cryosphere, 11, 1813-1834, 2017 Figure 1.Location of the Mont Blanc massif, its glaciers and mean annual rock surface temperature (MARST).Areas with MARST < 0 • C can be considered permanently frozen (Magnin et al., 2015a).Background topography: ASTER GDEM.(Haberkorn et al., 2017).However, snow control on permafrost evolution over long timescales is poorly known.This effect is neglected in our modelling approach.At depth, the temperature in hard rock mainly depends on the conductive heat transfer from the surface (Williams and Smith, 1989;Wegmann et al., 1998) and 3-D heat fluxes induced by the aspect-dependent rock surface temperature (RST) variability (Noetzli et al., 2007).Generally, modelling procedures of permafrost rock wall first calculate the RST and then solve the heat conduction equation to simulate subsurface temperature.Pioneer studies used distributed energy balance models to calculate the RST (Gruber et al., 2004) and simulated the subsurface temperature fields with the consideration of only (i) conductive heat transfer within idealised high mountain geometry and (ii) latent heat processes to account for water phase changes in the bedrock interstices (Noetzli et al., 2007;Noetzli and Gruber, 2009).
Due to the high computational efforts in energy balance approaches, statistical approaches were later adopted to compute the RST (Allen et al., 2009;Hipp et al., 2014;Myhra et al., 2015).The increasing amount of available RST time series in the European Alps has permitted the formulation of such statistical model for the entire Alpine range (Boeckli et al., 2012).This last model has been applied on a 4 mresolution DEM of the French part and of the Swiss and Italian borders of the Mont Blanc massif with local air temperature input data to map the mean annual rock surface temperature (MARST; Figs. 1 and 2, Magnin et al., 2015a).
In our modelling procedure, we use the MARST map available for the French part of the Mont Blanc massif to generate the initial RST condition.We run the transient simulations in the commercial hydrogeological software DHI-WASY FEFLOW version 7.0 by forcing the RST with climate time series from 1850 to 2100 and solving the heat conduction equation in 2-D with consideration of freeze and thaw processes in the bedrock interstices.

Conceptual approach
Rock wall permafrost is composed of rock, ice and air in non-saturated conditions.Sass (2003) approached alternating saturated/unsaturated conditions under freezing and thawing of the near-surface pore spaces of a rock wall by mean of geophysical soundings.However, the rates of saturation of alpine rock walls remain poorly understood, and therefore the numerical models of rock wall permafrost have so far considered a saturated, homogeneous and isotropic media (Wegmann et al., 1998;Noetzli et al., 2007;Hipp et al., 2014;Myhra et al., 2015).Nonetheless, such approaches have been shown satisfactory to simulate long-term temperature changes in alpine rock masses.Shorter timescale processes are clearly a hydrogeological problem.
In the scope of this study, we used FEFLOW combined with the plug-in piFreeze 1.0 which accounts for freeze and thaw processes.As a very first use, we adopted existing approaches of long-term simulations (saturated, homogeneous and isotropic media) to simulate transient thermal processes along centennial timescales.Further developments will use the potential of FEFLOW to simulate them with various saturation rates and fluid transfers.

Heat transfers: numerical approach
The conservation equation of energy for advective dispersive-diffusive transport of thermal energy depends on fluid flows (in saturated or unsaturated conditions, i.e.Darcy's law incorporated in continuity equation or Richards' equation), but works in pure conduction when flows are zero.It is usually expressed as follows (Diersch, 2002): with ϕ being the porosity (dimensionless), ρ L c L and ρ S c S the volumetric heat capacities (J m −3 K −1 ) of the liquid and solid phases respectively (obtained by combining the density ρ and the specific heat capacity c), the hydrodynamic thermal dispersion tensor (J m −1 s −1 K −1 ), which includes thermal conductivity, T the temperature (K) and q the apparent flow velocity from Darcy or Richards equation (m s −1 ).Equation ( 1) accounts for only one fluid phase and only one solid phase (water and rock respectively), which is the default use of FEFLOW in saturated conditions such as in this study.
The piFreeze plug-in working with Feflow7.0 adds the ice fraction to the solid phase in order to modify only the parameters of solid thermal conductivity and of solid thermal heat capacity in Eq. ( 1).The addition of the ice in the whole modelled medium is expressed throughout the bulk volume as ε a + ε w + ε i + ε r = 1, with ε a , ε w , ε i and ε r the bulk fractions of air (ε a = 0 in our case), water, ice and rock respectively.A relation is established between the bulk volume of ice and the bulk volume of liquid, which is the mass fraction per bulk volume of the unfrozen liquid to the total liquid mass, also called the freezing function F (Clausnitzer and Mirnyy, 2015): where ρ is the density of the corresponding phase (subscript i for ice and w for water).This function F decreases with the fraction of ice.For the freezing point T 0 , the ice forms gradually within a predefined temperature interval of the length Taking into account the expressions for the bulk volume ε and of the freezing function F described above, the thermal parameters of Eq. ( 1) are modified as follows: a.The solid thermal conductivity λ (W m −1 K −1 ) becomes b. Similarly, the solid volumetric heat capacity becomes where L f is the latent heat of the ice formation.

Thermal parameters
The thermal conductivity of the rock was set to 3 W m −1 K −1 , which stands for a conservative value for saturated granitic rock (Cho et al., 2009).However, the thermal conductivity of a saturated media does not only depend on the mineral properties, but also on the liquid or solid state of water, ice being up to six times more conductive than water at 0 • C (Williams and Smith, 1989).Thermal conductivity variations along freeze and thaw cycles are accounted for by piFreeze.The heat capacity of the rock was set to 1.8 MJ m 3 K −1 .
In addition to the usual adjustable parameters of FE-FLOW, piFreeze allows a user-defined freezing temperature, a customisable temperature interval for freeze/thaw combined with a linear freezing function, adaptable thermal properties of ice, a user-specified residual fluid content and a configurable latent heat.Such possibilities are highly promising for adapting the modelling approach to the natural conditions.
To account for the latent heat processes related to the freeze and thaw of the interstitial ice contained in pores and fractures, we took a 5 % porosity value following the procedure from Noetzli et al. (2007), which is the maximum value for dense crystalline rocks (Domenico and Schwartz, 1997) or lowly fissured crystalline rocks (Banton and Bangoy, 1999).Indeed, dense crystalline rock porosity without any fissure is usually below 1 %, while fractured crystalline rock porosity quickly reaches values greater than 5 %.The 5 % value chosen here accounts then for the ice contained in fractures since bedrock discontinuities are not included.
Water contained in artificial pore spaces is subject to supercooling, i.e. a temperature deviation from the equilibrium freezing point at 0 • C, until it reaches a spontaneous freezing point which depends on pore size and material (Alba-Simionesco et al., 2006).Geophysical experiments on various hard rock samples and under controlled laboratory settings have quantified a freezing point depression of (Krautblatter, 2009) due to the pressure and water salinity.To account for this supercooling characteristic of interstitial water, the freezing point T 0 was set to −1 • C in the piFreeze module, while the temperature interval T of the freezing function was set to 1 • C. The latent heat of fusion was set to 334 kJ kg −1 .

Rock surface temperature
We first extracted the topography and the MARST from the 4 m-resolution DEM (provided by RGD 73-74), and mapped the MARST over it to serve as upper boundary conditions along the NW-SE vertical transects (Fig. 2).The MARST map has been evaluated against 43 measurement points of RST from the multi-year time series of the nine RST log-gers installed around the AdM.The modelled MARST values tend to underestimate measured MARST values of sunexposed rock faces and to overestimate those of the shaded faces.Nevertheless, the mean bias (mean difference between measured and modelled MARST) of −0.21 • C (Magnin et al., 2015a) indicates a generally good approximation of the real-world MARST at this site.
The linear regression used to produce the MARST map has been formulated with the mean air temperature of the 1961-1990 reference period (homogenised by Hiebl et al., 2009), and the measured MARST were adjusted to the reference period.The MARST were adjusted by applying the difference in air temperature between the reference period and the years of the MARST measurements (Boeckli et al., 2012).In our modelling procedure, we considered that the MARST extracted from the map is representative of the year 1961 since the MARST has been mapped using the 1961-1990 MAAT.Then, we lowered this MARST by 1 • C to approximate the MARST at the end of the LIA (Auer et al., 2007;Böhm et al., 2010) and set up the initial RST condition at the upper boundary of the model domain (Fig. 3).MARST differences driven by topo-climatic factors clearly appear along the extracted profiles but are site specific.At the GM, the MARST difference between the SE and NW face is only 1 • C whereas it is 5 • C at the AdM and 6 • C at the GPA.These variable temperature differences for similar aspect differences (180 • ) are attributed to two factors: the differences in slope steepness and the local shading.The GM and AdM are isolated summits with no close shading.Conversely, the GPA NW face is located right below the Mont Blanc which shades the GPA west-exposed faces and lowers their MARST.The less steep NW slope and rounded summit of the GM receive solar radiation for a larger number of daylight hours than in subvertical settings with a more perpendicular incidence of the beams.
Starting from the initial RST representative of the LIA conditions, we first initialise 2-D steady-state temperature fields for the year 1850, and then run transient simulations using reconstructed, measured and projected climate time series up to 2100 (see Sect. 3.3).

Model geometries
Below the topographical profiles, a box with a height of 5000 m was added to shut off the model geometry of each site.A constant geothermal heat flux of 85 mW m −2 (Medici and Rybach, 1995;Maréchal et al., 2002) was set up as a lower boundary condition.Above these boxes, a finite element mesh with triangular elements was generated to discretize the subsurface material (Fig. 3).Even though the spatial resolution of the boundary conditions is 4 m, we refined the meshes close to the surface in order to better represent the near-surface temperature gradient.The spatial resolution of the initial RST, based on the 4 m resolution map, was refined accordingly using linear interpolation.This mesh and RST refinement does not provide much information at depth, nor improve the quality of the models, but it facilitates the model evaluation.At greater depth, we kept a mesh size of 4 m, in coherence with the resolution of the input data.This approach resulted in 8548 nodes and 16 141 mesh elements at the AdM, 5844 nodes and 10 952 elements at the GM and 12 087 nodes and 23 344 elements at the GPA (Fig. 3).
On the AdM site, 37 observation points were defined between the surface and 10 m depth of the NW and SE faces at the location of the boreholes (Sect.2.1).The mesh was refined along these observations points (Fig. 3) and simulated bedrock temperature is extracted at user-defined time steps during the transient simulations to serve for model evaluation (Sect.4.2.1).

Initial condition
To define an initial 2-D temperature field in the model geometries, we ran the model with the upper boundary condition (the 1850 RST) and lower boundary condition (the geothermal heat flux) until steady-state conditions were reached, balancing the respective influence of the atmosphere and geothermal heat fluxes.Steady-state conditions were reached after 80 000 ± 10 000 years.After this initialisation procedure that provides an initial condition for 1850, we ran transient simulations with air temperature forcing time series from 1850 to 2100.

Forcing data
The transient simulations are forced with air temperature time series created from various sources of data.The temporal resolution of these forcing data was gradually refined with the increasing quality of the available data and the pe-riods of interest (Fig. 4).For the period 1850-1961, no continuous air temperature measurements are available for the Mont Blanc massif.Therefore, we assumed a linear increase of 1 • C between 1850 and 1961 (Auer et al., 2007;Böhm et al., 2010), and run the simulations at an annual time step.A sensitivity analysis using larger time steps did not change the final results for the periods of interest.
In the framework of the MARST mapping (Magnin et al., 2015a), a climate time series was created at a monthly time step based on measured temperature values in Chamonix for the period 1961-1990.We used this monthly time series and extended it up to 1993 to force the model between 1961 and 1993, which constituted a first step in temporal resolution refinement of the forcing data.
In 1993, Météo France started continuous records of air temperature at hourly time step.From these hourly records, daily air temperature can be reliably calculated.Therefore, we forced the transient models with this daily air temperature time series between 1993 and 2015.
Finally, two contrasted scenarios were retained for the 21st century.Time series consist of daily 2 m air temperature simulated by the IPSL-CM5A-MR Earth system model (Dufresne et al., 2013)  For RCP4.5 anthropogenic greenhouse gases emissions peak around 2040 and then decline (moderately optimistic scenario), while for RCP8.5 emissions continuously increase throughout the century (pessimistic scenario).The IPSL-CM5A-MR model was chosen because its basic state is very close to the recent observational records during the first years of the 21st century (Fig. 4), and its response to the radiative forcing throughout the century is close to the median of the CMIP5 models, ensuring a representative behaviour to estimate long-term evolutions.Forcing time series are obtained by extracting data at the closest grid-point (1336 m a.s.l.) to the Mont Blanc Massif.

Results
The results of the simulations are presented in three steps.First, we describe the permafrost conditions and changes between the steady-state at the end of the LIA to timedependent conditions during the recent period.The recent conditions are illustrated through model snapshots in September 1992 and September 2015.In a second step, the model outputs for the recent period (2010-2015) are compared to an independent data set of real-world conditions for assessing the quality of the simulations along the 20th and early 21st centuries.Finally, after model evaluation, thermal conditions by the end of the 21st century in response to RCP4.5 and 8.5 forcings are presented.Equilibrium conditions for the end of the LIA (1850) are displayed in Fig. 5a for the three sites.In 1850, the GPA and the AdM showed cold permafrost (< −2 • C).At the GM, a cold permafrost body was characterising the NW subsurface and extending below the SE face between 3260 and 3280 m a.s.l.Warm permafrost was already occupying most of the site, including the summit area.
The shape of the isotherms varies from one site to another, depending on the topographical settings.The steepest site (AdM) shows subvertical isotherms down to 3720 m a.s.l., where they incline downwards and towards the NW to become more oblique.In the top part of the less steep site (GM), the −2 • C isotherm is rather oblique while the −1 • C one is subhorizontal in the lower part.In the GPA, isotherms are vertical in the top part and oblique in the middle part.In the lower part of the SE face, in the > −5 • C area, isotherms obliquity declines to become more parallel to the upwards geothermal heat flow.
The modelled temperature gradients directly depend on the temperature difference between NW and SE flanks (Sect.3.2.1):small temperature gradients are visible in the GM cross section in accordance with the initial RST difference, whereas they are higher in the two other sites with more severe topography and more contrasted RST.

Transient temperature fields along the 20th and early 21st centuries
Figure 5b displays time-dependent conditions for the early 1990s, while Fig. 5c demonstrates those in 2015 after the two past decades of strong air temperature increase (Fig. 4).
During the 20th and early 21st centuries, permafrost has degraded in all the three sites.Warm permafrost has extended over the entire GM site and below the AdM SE face.This degradation shows site-specific features in terms of isotherm shapes and temperature field distributions.At the GM, the cold permafrost body has subsisted until the early 1990s below the NW face, but has narrowed down to two small bodies.Meanwhile, the lower limit of the −1 • C Two decades later, this body has narrowed: its width decreased by 20 m; by depth, its lower limit rose up to 4070 m a.s.l. and the highest limit decreased to 4260 m a.s.l.
The air temperature increase experienced since the end of the LIA had variable effects depending on the site geometry and initial RST.The existing permafrost data in the Mont Blanc massif allow for an evaluation of the model outputs.

Model evaluation
Modelled subsurface temperatures in rock walls are rarely benchmarked against measured borehole values due to the scarcity of subsurface temperature measurements, and because approved heat conduction models are assumed to be accurate enough to reproduce temperature in relatively simple thermal systems such as rock wall permafrost.This assumption is viable only if the upper boundary (the RST) is accurately simulated.Therefore, validation of modelled RST is required before implementing it in heat conduction schemes, especially when it is generated from complex energy balance simulations with many sources of uncertainties related to the high number of input data and parameters (Gruber et al., 2004;Noetzli et al., 2007).The quality of the initial RST was already mentioned in Sect.3.2, based on the study from Magnin et al. (2015a).Here, we propose to evaluate simulated temperature in depth by means of measured borehole temperatures and an electrical resistivity tomogram.For better visibility, only four selected modelled temperatures of the year 2010 for each borehole, encompassing the four seasons, are displayed in Fig. 6.

Daily features
The model reproduces the bedrock temperature below 6 m depth very well, especially in the NW borehole (Fig. 6a,  b).Measured bedrock temperatures at the shallowest depths, which are affected by snow cover, fractures and strong solar irradiation are not taken into account in the modelling approach.This explains the greater discrepancy between modelled and measured ground temperature in winter (presence of snow) whereas summer temperatures show a better fit, especially in the NW face where the effect of solar radiation is weak.
The temperature profiles recorded in the NW borehole are significantly affected by an open fracture at 2.5 m depth which locally cools the bedrock due to air ventilation, especially during winter, whereas above the fracture, the insulating effect of a snow patch accumulating on a ledge at the borehole entrance warms the temperature profile down to the fracture depth (Magnin et al., 2015c).This is visible on the profile from 1 January 2010 (Fig. 6a, b): the upper part shows a small temperature gradient due to snow insulation, while a stronger temperature gradient is visible below the fracture due to its shortcut effect between the air and the subsurface.During the years 2010-2015, the maximum difference between measured and modelled daily temperature at 10 m depth in the NW borehole is 0.5 • C and the mean difference of the 72 observation points (corresponding to 72 days between 2010 and 2015) is only 0.01 • C.
In the SE borehole, the deepest measured temperatures seem less well reproduced by the model, but remain of reasonable accuracy.We observe a maximum difference of 0.7 • C between measured and modelled value at 10 m depth for the observation period, with a mean difference of 0.1 • C between the 60 measured points (two 3-4 months interruptions of the borehole records decrease the number of available data) and simulated temperatures at the same date.On the SE face, the almost continuous snow cover from autumn to spring/early summer (Magnin et al., 2017) and the high solar irradiation both affect the rock temperatures, which is not considered in the modelling approach.The solar radiation effect is visible on the measured profile from 1 July (Fig. 6a, b).

Annual features
On an annual average, differences between the measured and modelled temperature values are smaller than on a daily basis (Fig. 6c).At 10 m depth, in the worst case (2010), the modelled value is 0.2 • C higher than the measured one at the NW borehole, with a mean difference of 0.01 • C during 2010-2015.At the SE borehole, the maximal difference between measured and modelled value is 0.04 • C, but only two full years of observations (2010 and 2011) are available for model evaluation, the 2012-2015 time series being affected by significant data gaps.Therefore, the mean difference between observed and measured value was not calculated.The model satisfyingly reproduces the negative temperature gradient along the SE profiles, resulting from the heat sink effect of the opposite north face (Noetzli et al., 2007;Magnin et al., 2015c).On the NW face, the model shows a significantly lower temperature gradient than the measured one.Since 3-D effects seem well reproduced on the SE face, the difference between the measured and modelled temperature gradient on the NW face may result from the cooling effect of the fracture rather than from 3-D effects.
Borehole temperatures provide information at a very fine scale and are thus not appropriate for evaluating 2-D models with a multi-metric resolution.However, they are much more suitable than RST measurements since the temperature at 10 m depth results from the heat transfer of a multi-metric area at the surface.Further evaluation is possible using an electrical resistivity transect, which has been proven to be a relevant support in improving model evaluation (Noetzli et al., 2008).

Evaluation with ERT at GM
In October 2012, an ERT sounding was performed along a 160 m long profile of the GM NW face with an electrode spacing of 5 m.Field measurements were combined with labwww.the-cryosphere.net/11/1813/2017/The Cryosphere, 11, 1813-1834, 2017 Given the remarkable capability of the transient simulations to reproduce current temperature conditions at the AdM borehole locations and in the GM NW face, the model can be judged to be accurate enough to consider future long-term scenarios.

Future scenarios
In Fig. 8, time-dependent conditions for September 2099 in response to RCP4.5 and RCP8.5, which respectively result in +3 and +5 • C increases in air temperature from the present day to the end of the 21st century, are displayed.Future scenarios result in highly contrasted permafrost conditions, from an almost total disappearance (only relict body subsisting at the core) to the preservation of entire permafrost conditions, depending on both the RCP and the site.
At the GM, a relict body has subsisted in the internal part and below the topographical summit in both RCPs.Unlike the 20th century changes, with the isotherms retreating towards the NW face and in the interior of the summit, the permafrost body retreats downward in the core of the massif during the 21st century.Temperature gradients depend on the RCP, but are stronger than during the 20th century in both cases, with differences of 4 and 6 • C between the shallow layers and the deepest part for RCP4.5 and 8.5 respectively.
At the AdM, a 10-12 m-wide body of cold permafrost (−3 to −2 • C) still subsists under RCP4.5, located below the NW face and between 3710 and 3770 m a.s.l.Warm permafrost is present throughout, occurring down to 7-10 m below the top and at 20-25 m depth in the central part of the SE face.Thus, permafrost has disappeared in the AdM SE face.In the most pessimistic scenario, permafrost has totally disappeared from the AdM summit, but warm permafrost will still exist in the NW rock wall at 3750 m a.s.l.As with the GM, the permafrost body retreats downwards.
At the GPA, the entire geometry has been affected by the projected air temperature increases.In the case of RCP4.4140 m a.s.l.In that case, the coldest areas are retreating in the core of the NW half of the summit.

Discussion
In our study, we simulate the long-term temperature evolution in three rock wall permafrost sites with different topographical settings.We use a relatively simple approach since our transient simulations are only forced by air temperature changes applied to the RST and of heat conduction and latent heat exchange processes.The limitations in our modelling approach are examined prior to a discussion of the implications of our results for determining future permafrost changes in the steep rock walls of the Mont Blanc massif and for investigating rock wall destabilisation.

Model limitations
Model uncertainties arise from misrepresentations of the processes in the thermal system, unknown physical properties and errors in input data (Gupta et al., 2005;Gubler et al., 2013).Regarding the modelling approach adopted in this study, uncertainties mainly arise from five different sources.

Initialisation and forcing data
Assumptions and simplifications were necessary to generate the initial 2-D temperature field at the end of the LIA and run transient simulations.Starting with LIA equilibrium conditions ignores possible transient effects from Würm and Holocene climate variability.However a thermal perturbation in depth due to millennial timescale changes is unlikely given the geometry of the investigated sites (max.width of 350 m), and results from previous studies (Kohl, 1999;Noetzli and Gruber, 2009).
Simulated 2-D temperature fields seem to accurately reproduce real-world data as shown by the comparison of modelled temperature with the AdM borehole temperatures and the GM ERT transect.This underscores two strong points.Firstly, the similarity between measurements and model (Sect.4.2.2) in the uppermost layer of the rock walls (i.e. the 10 and 25 uppermost metres at the AdM and the GM respectively) suggests that the strategy to run transient simulations is relevant enough to simulate the permafrost changes since the LIA up to the current period.In the meantime, it emphasises the quality of the forcing data that have the advantage of being partly constructed from direct measurements of air temperature in the Mont Blanc area.This ensures a better representation of the local climate variation compared to air temperature time series extracted from kilometre-scale climate models.
The resolution of the input data is one of the most challenging issues when forcing permafrost models in highly heterogeneous land surfaces such as mountain terrains (see next section; Fiddes et al., 2015).Downscaling methods are under development for mountain terrains (Fiddes andGruber, 2012, 2014;Fiddes et al., 2015), but currently available homogenised air temperature data sets at a kilometre-scale remain better suited to drive numerical models on larger areas with coarser spatial resolution (Jafarov et al., 2012;Westermann et al., 2013).
More sophisticated approaches capable of simulating complex energy exchanges at the bedrock-atmosphere interface have been proposed to model surface temperature of steep mountain slopes using specific algorithms to compute solar radiations (Stocker-Mittaz et al., 2002;Gruber et al., 2004;Salzmann et al., 2007) and heterogeneous snow accumulation effects (Pogliotti, 2011;Haberkorn et al., 2015aHaberkorn et al., , b, 2017;;Magnin et al., 2017).However such complex approaches induce a high degree of uncertainty due to the numerous input variables owning to their respective sources of error.As a consequence, a certain proportion of the modelled RST variability between different model outputs is in the range of the model noise.Furthermore, relevant data and parameterisation techniques are still missing for long-term transient models of rock wall permafrost accounting for snow effects.

Future scenarios
In this study, possible climate evolutions over the future decades are obtained from one single climate model (namely, IPSL-CM5A-MR) and two scenarios for greenhouse gas emissions, out of four possible RCPs, and 20 to 40 models (depending on the model versions and the type of the climate simulations) that participated in the 5th Assessment Report of the IPCC (2013).The choice of the climate model was motivated by its realistic steady-state in present-day conditions when compared to observational time series and its median response to greenhouse effect evolutions.The retained radiative forcings for the 21st century consist of the RCP4.5 and RCP8.5, that is, a moderately optimistic and a pessimistic scenario (see Sect. 3.3.2, Moss et al., 2008).These RCPs can be considered currently as reasonable estimations of the lowest and highest air temperature changes likely to happen, and thus plausible contrasted boundaries within which the permafrost could evolve.According to the Paris agreement on climate change adopted on December 2015 during the COP21, climate change should be limited "well below 2 • C" above pre-industrial levels with more ambitious targets updated every 5 years.To date, however, the reductions in greenhouse gas emissions planned by the participating countries lead to an estimated global warming of roughly 2.5-3 • C, which lead us to discard the more optimistic RCP2.6 scenario, describing a rapid decrease in emissions as soon as the early 2020s, but makes RCP4.5 the most probable pathway for the 21st century.
Precipitation and snow were neglected in our modelling approach due to (i) the marginal and local effects of precipitation on long-term rock wall permafrost changes, (ii) the complexity of snow accumulation patterns on steep slopes, which do not only depend on precipitations but also on a variety of parameters such as slope roughness and steepness, aspect (melting), and wind patterns, and (iii) the large uncertainties in precipitation projections (Heinrich et al., 2013) as well as the persistent difficulties of current climate models to simulate them.

Model dimensions and resolution
Rock wall permafrost is highly sensitive to climate change since it is subject to multi-directional heat transfers from the different sides of a summit (Noetzli et al., 2007).In this study, we assess rock wall temperature changes in 2-D only.Therefore, the simulated changes in permafrost distribution are possibly slightly underestimated since the signals only penetrate through two sides (NW and SE).Sensitivity analyses have shown, however, that 2-D simulations show similar long-term transient temperature pathways than in 3-D situations (Noetzli and Gruber, 2009).Based on these previous findings, it can be assumed that the 2-D temperature fields are acceptable for drawing reasonable patterns of permafrost distribution and changes.
Additionally, the model resolution (4 m), defined by the initial RST resolution (Sect.3.2.1),and later refined for spatial discretisation of the model domain (Sect.3.2.2), is sufficient to represent the main topographical control on the rock wall temperature distribution, and to drive longterm changes.Permafrost distribution modelling often suffers from coarse topographical resolution which does not represent the natural variability of environmental parameters (Etzelmüller, 2013).Metric resolution is essential for a realistic simulation of rock wall permafrost at the summit scale (Magnin et al., 2015a).The satisfying quality of the 2-D model outputs for the current period (Sect.4.2.1)confirms that the model resolution is accurate enough to address longterm permafrost changes at the scale of the selected sites.Simulations at shorter spatial and temporal scales, especially in the uppermost layers, would certainly require higher spatial resolution of the model domain and consistent input data.

Thermal parameters
Subsurface thermal parameters have been defined by published values for hard rock.The thermal conductivity and heat capacity are assumed to be homogeneous in the model domain whereas they in fact vary with frozen/thawed conditions in the natural environment due to the changing properties of the interstitial ice/water.The variable state of the interstitial ice or water results in seasonal variation in the The Cryosphere, 11, 1813-1834, 2017 www.the-cryosphere.net/11/1813/2017/thermal conductivity of the porous and saturated rock media (Wegmann et al., 1998) and longer-term permafrost changes.This is accounted for in our modelling approach throughout Eqs. ( 3) and ( 4).In natural conditions, however, this changing conductivity is heterogeneous in the rock mass due to variable fracture density and porosity, which is not taken into account in our study.Based on the model evaluation (Sect.4.2), this lack of consideration of heterogeneous porosity and related latent heat processes is not an obvious limitation for simulating thermal fields in depth.However, to simulate near-surface thermal processes at finer spatial and temporal resolutions, heterogeneous porosity would have to be taken into account.
We selected conservative values for granitic rock but tested the model sensitivity to different thermal conductivity values (2.7 and 3 W m −1 K −1 ).Results confirmed findings from previous studies: a lower conductivity increased the geothermal heat flux control (Maréchal et al., 2002;Noetzli et al., 2007) but did not lead to substantial changes in the modelled temperature fields (Kukkonen and Safanda, 2001).Furthermore, the thermal conductivity is naturally anisotropic (Goy et al., 1996), whereas it is considered isotropic in our modelling approach.Increasing the conductivity in horizontal directions increases the topo-climatic control, whereas increasing the conductivity in a vertical direction gives more importance to the geothermal heat flux (Noetzli and Gruber, 2009).
The influence of the geothermal heat flux appears to be of relatively high importance for running steady-state conditions.Equilibrium conditions without the influence of the upwards and deep-seated flow only depend on climate control, which can lead to highly different conditions than when balancing the respective influences of the geothermal heat flow and the climate.In Fig. 9, an example of the geothermal heat flux effect is given for the GM site, which is the most sensitive to this parameter given its relatively low relief and wide geometry.Without geothermal heat flow, the equilibrium condition for the LIA shows almost entirely cold permafrost conditions with more vertical isotherms.However, when running a transient simulation, the influence of the geothermal heat flux is less significant on the temperature range.It only affects the shape of the isotherms leading to permafrost retreat in the core of the summit.Simulations without geothermal heat flux were also run with a lower thermal conductivity (2.7 W m −1 K −1 ) and a monthly time step, which had no impact on the results of this study (Fig. 9).Only the uppermost layers show different temperature ranges due to different time steps in the forcing data (daily versus monthly), but shallow thermal processes are beyond the scope of this study.Idealised test cases suggested that mountain summits are decoupled from the deep-seated geothermal flow influence (Noetzli et al., 2007), but such settings are not representative of most of the alpine study cases, which are more or less sharp and elevated.

Heat transfer processes
Energy transfer inside the rock mass is mainly driven by heat conduction processes, whereas fluid flows can be, in a first approximation, neglected to simulate long-term changes (Kukkonen and Safanda, 2001).Nevertheless, advective heat transport by water circulation along fractures may locally warm the bedrock at depth (Hasler et al., 2011a).Conversely, air circulation in open clefts would instead cool the bedrock (Hasler et al., 2011b;Magnin et al., 2015c).These nonconductive heat transfers are not accounted for in our modelling approach.The evaluation with borehole temperatures clearly shows their effect in the shallow layers (the first 6 m below the surface, Sect.4.2.1) while temperatures at deeper layers are accurately represented with consideration of 2-D heat conduction only.
The melting of ice bodies in the fractures may also be expected in the natural environments, as suggested by recent observations in rockfall scars (Ravanel et al., 2010b).Ice melting may significantly dampen and lower the rate of temperature changes by depth by the consumption of latent heat (Wegmann et al., 1998;Kukkonen and Safanda, 2001).Some ice-filled fractures can turn into thawing corridors during the thawing season (Krautblatter and Hauck, 2007), which can favour the melting of ice-filled fractures, and degrade rock wall permafrost in unexpected areas and depths (Hasler et al., 2011a).Such processes were approximated with a relatively high porosity value in the model domain (5 %), which fails to represent the anisotropic and heterogeneous character of such processes.Further developments to gain data on rock wall structure are highly encouraged.

Past and future permafrost degradation
Taking into consideration the model limitations, only the long-term and summit scale permafrost changes can be considered despite the fact that the shape of the isotherms remains somewhat uncertain due to limitations in the thermal parameters.In this section, we summarise the most probable permafrost changes since the end of the LIA to the current period, and examine possible changes by the end of the 21st century.The patterns of simulated 2-D permafrost degradation from the termination of the LIA to the end of the 21st century show a dependency on the topographical settings (summit width), the bedrock temperature (latent heat effects) and the intensity of the climate signal.

Permafrost degradation since the LIA termination
The rate of change between equilibrium conditions in 1850 and the early 1990s was in the same range as the one experienced during the past three decades.Indeed, for the first time period, the thermal perturbation has been detected down to 30 m below the surface (the depth at which the shape of the www.the-cryosphere.net/11/1813/2017/The Cryosphere, 11, 1813-1834, 2017 steady-state isotherms have changed), whereas it reaches at least twice that depth in 2015.Narrow peaks such as the AdM have been entirely affected by air temperature increases since the end of the LIA due to the short distance between the sides from which the signal penetrates, and the resulting intense lateral heat fluxes, especially at the top (Noetzli et al., 2007).
In wider geometries, where N and S-facing surfaces are both distant by more than 150 m, such as the GPA and the GM, the core has remained intact.Temperature changes during 1850-2015 were not as high at the GM as for the two other sites, and were not as high as the air temperature change, because the entire summit was in the temperature range within which latent heat exchanges occur.This pattern is aligned with the global trend, showing the delaying effect of latent heat uptake in boreholes with temperatures close to 0 • C (Romanovsky et al., 2010).The rate of temperature change strongly increases during the 21st century when the GM temperature becomes almost entirely positive.
Between the 1990s and 2010s, cold permafrost disappeared and warm permafrost largely penetrated in the shallow layers of the GM and AdM SE faces respectively.Quantitative interpretations about the permafrost lower boundary and its changes are limited by the discontinuous character of mountain permafrost mainly governed by local conditions (Etzelmüller, 2013), topography and transient controls; these different influences are difficult to distinguish (Noetzli and Gruber, 2009).Nevertheless, results of our 2-D simulations clearly suggest that climate change in between the LIA termination and the 2010s, especially since the 1990s, have led to permafrost disappearance below the S-exposed faces at least up to 3300 m a.s.l.(top of the GM), but not above 3700 m a.s.l.(foot of the AdM SE face).Thus, lower boundaries of snow-free and hard rock wall permafrost lie within this elevation interval, but Magnin et al. (2015a) have suggested that isolated permafrost bodies could exist down to 2800 m in favourable S-exposed slopes where conduction is not the prominent heat transfer process (due to high fracturing for example).Failing to account for snow and fracturing parameters could lead to a 3 • C underestimation of bedrock temperature in S-exposed rock walls (Hasler et al., 2011b).
Warm permafrost has thickened below the S-exposed face at least up to 3850 m a.s.l.during the past two decades while it was extending below the N-facing slopes up to 3300 m.Although warm permafrost reaches 20 m depth below the Sexposed AdM face in 2015, this depth cannot be extrapolated to other sites due to the site-specific effects of the opposite N face.

Permafrost degradation during the 21st century
By the end of the 21st century, even the core of the relatively wide summits such as the GPA will be strongly affected by the projected increase in air temperature.Warm permafrost in the S-exposed faces will extend with depth by at least up to 4300 m a.s.l., and in the N-facing faces according to the RCP4.5.Permafrost will certainly disappear in all aspects below 3300 m, and up to 3850 m a.s.l. at least (top of the AdM) in the S-exposed faces.Following the most pessimistic scenario, permafrost will disappear in the subsurface of the S-exposed rock walls at least up to 4300 m a.s.l.(top of the GPA), while cold permafrost will still exist at the same elevation in N-facing slopes if the S-face influence does not affect it.At lower elevation such as at the AdM, warm permafrost will still occur below the N-exposed faces, but could disappear in the narrowest sections due to the S-facing slope influence.When both mountain sides do not allow for permafrost conditions any more, such as at the GM, the permafrost body will retreat downwards, in the core of the summit.
These degradation patterns locally depend on the topographical control, especially the summit width.This either reinforces the intensity of the climatic signal in narrow geometries by mixing S-facing and N-facing slope influences, or maintains independency between the shallow layers of the S-and N-exposed slopes when it is wide enough (multiplehm).Site-specific patterns make the concept of "lower limit" not suitable to describe rock wall permafrost distribution and changes.

Conclusions
Former studies have described the 3-D processes and patterns of permafrost degradation in idealised high mountain geometries of the European Alps following former IPCC re-ports (Noetzli et al., 2007).In this study, we investigated permafrost degradation from the LIA steady-state until the end of the 21st century at three summits that are representative of the topographical and permafrost settings of the Mont Blanc massif rock walls, and with one of the most recent climate models.Simulations were performed in vertical cross sections with local air temperature forcing data and two possible air temperature scenarios, accounting for moderately optimistic and pessimistic 21st century pathways.They provide insights into the past and future changes experienced by rock wall permafrost in the Mont Blanc area, relevant for geomorphological applications.The main outcomes are as follows: 1. Thermal conditions for the current period (2010)(2011)(2012)(2013)(2014)(2015) are remarkably well represented when comparing simulated temperature fields to an independent data set.Our modelling approach is therefore well suited to run longterm transient simulations in rock walls.
2. Thermal perturbation induced by climate change since the end of the LIA is visible down to 60-70 m below the surface; i.e. the narrow Aiguille du Midi peak has been affected in its entirety.
3. Between the 1990s and the 2010s, (i) our simulations suggest that permafrost has disappeared in the shallowest 20 m of the S-exposed faces and warm permafrost has extended with depth within the N-exposed faces up to at least 3300 m a.s.l., (ii) warm permafrost has extended at least up to 3850 m and penetrated down to 20 m depth within the S-exposed faces at 3700-3800 m a.s.l.
4. At the end of the 21st century, only relict permafrost bodies will persist in the core of the wide summits below 3300 m a.s.l.
5. Considering a moderately optimistic scenario (RCP4.5),permafrost will disappear during the 21st century in any topographical setting lower than 3300 m a.s.l., and by depth within the S-exposed faces up to 3850 m a.s.l.Warm permafrost will extend below the N-exposed faces at similar elevations and at least up to 4300 m a.s.l.below the S-exposed faces.Cold permafrost will still exist below the N-exposed faces above 4000 m a.s.l.
6. Considering the most pessimistic scenario, permafrost will disappear from the S-exposed faces at least up to 4300 m a.s.l., and lead to permafrost disappearance in the narrow summits below 3850 m.Without the influence of a close S-exposed face, warm permafrost could persist at least down to 3700 m a.s.l. in the N-exposed faces.Cold permafrost will exist in spite of significant warming, at depths of the N-facing slopes higher than 4000 m a.s.l.
www.the-cryosphere.net/11/1813/2017/The Cryosphere, 11, 1813-1834, 2017 7. Locally, permafrost evolution patterns may be slightly different due to heterogeneous bedrock and nonconductive heat transfer occurring in fractures, as well as snow patch effects which are not taken into account in the modelling procedure.
8. The simulations show that the patterns of permafrost changes are mainly governed by the local topographical control, which emphasises the specificity of rock wall permafrost and restricts the extrapolation of the results.9. Transient simulations provide useful information for analysis of the thermal conditions at rockfall locations, but analysis of specific events would require the combination of hydrologic, mechanical and thermal models at shorter timescales with consistent input data and thermal parameters.

Perspectives
The main perspectives should include the improvement of the thermal modelling by downscaling climate models at the local scale in order to better account for topographical effects and by simulating 3-D and non-conductive heat transfer processes in the rock mass in order to approach short-term and shallow layers processes.These developments must be designed to bridge the gap between thermal and mechanical models in order to approach rockfall triggering mechanisms (Krautblatter et al., 2012).Climate downscaling exercises could be undertaken through different approaches.Dynamical downscaling consists of numerically resolving the regional atmospheric thermodynamics at high spatial resolutions in a limited area model, forced laterally by a coarser-resolution global model (e.g.Laprise, 2008).Computation costs are extremely high, and even kilometre-scale resolutions such as those achieved by the regional climate modelling community for climate change studies (e.g.Giorgi and Mearns, 1991;Feser et al., 2011) may strongly dampen small-scale topographic contrasts.As an alternative, statistical downscaling attempts to establish empirical statistical relationships between largescale variables, such as air temperature evolution throughout the century, and local-scale variables, mostly derived from the topography.This computationally efficient approach could allow for estimations at a hectometric scale.
FEFLOW has the potential to simulate hydrological processes but relevant data are required to achieve a realistic parameterisation of these processes.Geophysical measurements appear to be the most efficient approach to gain such in situ data (Mewes et al., 2016) and we consider that the use of geophysical results to parameterise complex thermal modelling approaches is one of the future challenges.This will constitute a first step in modelling the thermohydromechanical processes leading to bedrock failure.
Even though permafrost degradation appears to be the most prominent factor in the ever-increasing rockfall activity observed in the European Alps (Gruber and Haeberli, 2007;Deline et al., 2015;Luethi et al., 2015) and other alpine ranges such as in New Zealand (Allen et al., 2009), they are triggered by a complex combination of factors (such as the fracturing characteristics, the lithology, the glacial and periglacial influences, etc.) that thermal models currently do not represent.Future developments must account for hydrological and mechanical processes by including structural data such as the fracture characteristics and filling (air, ice, etc.), variable saturation and hydrostatic pressures, in combination with high spatial and temporal thermal simulations.Current developments attempt to link thermal and mechanical models (Mamot et al., 2016), and such an integrative approach will make major contributions towards understanding rock wall destabilisation.
Data availability.Data are available upon request to the first author.The borehole data are visible online at http://gtnp.arcticportal.org/ (Magnin, 2015b) and on http://www.permasense.ch/en.htmlwhere they can be downloaded as well.The ERT data will be uploaded on similar data repository in within the next year.The MARST map used for the model initialisation can be downloaded at https://hal-sde.archives-ouvertes.fr/hal-01120617(Magnin et al., 2015b).
Author contributions.FM conducted the study.She designed the modelling approach, prepared the input data, performed the model evaluation, interpreted the results, and wrote the manuscript with inputs from the co-authors.JYJ implemented the modelling approach in FEFLOW, performed the model runs and contributed in writing of the methods.LR is the responsible of the WP Permafrost in the VIP Mont Blanc project, and contributed in designing the study and writing of the manuscript.BP and JP extracted the climate forcing variables, contributed in the choice of the RCP forcings and in writing of the methods and the discussion.PD contributed to the writing of the manuscript.
Competing interests.The authors declare that they have no conflict of interest.Special issue statement.This article is part of the special issue "The evolution of permafrost in mountain regions".It is not associated with a conference.acknowledge the three anonymous reviewers and the two editors for their relevant comments, questions and suggestions, which have highly contributed in improving the paper quality.We thank Ben Snook for having reviewed and corrected the English grammar.
Edited by: Christian Hauck Reviewed by: three anonymous referees

Figure 2 .
Figure 2. Topographical profile locations on the three sites.

Figure 3 .
Figure 3. Boundary conditions: initial RST plotted over the models meshes for the three sites.The spatial resolution of the initial RST (coloured dots) has been refined at the mesh scale by linear interpolation.Below the topographies, a box of 5000 m elevation and a constant geothermal heat flux of 85 mW m −2 were set up.
, which participated within the 5th Coupled Model Intercomparison Project (CMIP5, Taylor et al., 2012)/AR5 of IPCC (2013).For this study and for climate projections in future decades, we used two contrasted radiative forcing scenarios, namely the Representative Concentration Pathway (Moss et al., 2008) RCP4.5 and RCP8.5.They correspond respectively to increases of +4.5 and +8.5 W m −2 for 2100 relative to pre-industrial values, resulting in air temperature increases of +3 and +5 • C by the end of the 21st century according to the comparison of the measured mean air temperature of the 1980-2010 period and the projected mean air temperature for the period 2070-2100.

Figure 4 .
Figure 4. Forcing air temperature data displayed at an annual time step.For running the transient models, the time step was refined for some periods as described on the figure.

4. 1
Permafrost evolution from the LIA to the current period 4.1.1Steady-state at the LIA termination

Figure 5 .
Figure 5.Initial steady-state (a) conditions and time-dependent conditions in September 1992 (b) and September 2015 (c) for the three studied sites (note differing figure scales).

Figure 6 .
Figure 6.Evaluation of modelled bedrock temperature against measured bedrock temperature in the AdM NW and SE boreholes at daily time step (a, b) and annual time step (c).

Figure 7 .
Figure 7. Model evaluation (a) against ERT transect (b) for October 2012 at the GM.The red line on panel (a) delineates the contour of the ERT transect.
5, cold permafrost is largely still present.Warm permafrost has penetrated into the SE face, reaching 40 m depth at 4100 m a.s.l.The −5 • C isotherm has kept a similar curving shape to in 2015, but has crept towards the NW face.A 50 m-wide colder body (< −6 • C) still persists at depths of the NW face between 4050 and 4250 m a.s.l.The simulations for the RCP8.5 show different results: permafrost has disappeared in the shallowest 20-30 m of the SE face and warm permafrost exists at to 40-50 m depth of this face, with isotherms roughly parallel to the rock surface.The −5 • C isotherm forms the coldest body which is 40 m wide, rounded and located at 40 m depth of the NW face between 4060 and

Figure 8 .
Figure 8. Time-dependent conditions in September 2099 after RCP4.5 and 8.5 forcings for the three investigated sites.

Figure 9 .
Figure 9. Bidimensional thermal fields at the GM site for steady-state conditions (1850) and transient conditions (2010).Models in the left panels (a, c) are run at daily time step with a geothermal heat flux (GHF) of 85 mW m −2 and a thermal conductivity of 3 W m −1 K −1 , whereas those on the right panels (b, d) are run at monthly time step without geothermal heat flux and with a thermal conductivity of 2.7 W m −1 K −1 .