Articles | Volume 18, issue 4
Research article
30 Apr 2024
Research article |  | 30 Apr 2024

Surface heat fluxes at coarse blocky Murtèl rock glacier (Engadine, eastern Swiss Alps)

Dominik Amschwand, Martin Scherler, Martin Hoelzle, Bernhard Krummenacher, Anna Haberkorn, Christian Kienholz, and Hansueli Gubler

We estimate the surface energy balance (SEB) of the Murtèl rock glacier, a seasonally snow-covered permafrost landform with a ventilated coarse blocky active layer (AL) located in the eastern Swiss Alps. We focus on the parameterisation of the turbulent heat fluxes. Seasonally contrasting atmospheric conditions occur in the Murtèl cirque, with downslope katabatic jets in winter and a strongly unstable atmosphere over the heated blocky surface in summer. We use a novel comprehensive sensor array both above the ground surface and in the coarse blocky AL to track the rapid coupling by convective heat and moisture fluxes between the atmosphere, the snow cover, and the AL for the time period September 2020–September 2022. The in situ sensor array includes a sonic anemometer for eddy-covariance flux above-ground and sub-surface long-wave radiation measurements in a natural cavity between the AL blocks. During the thaw seasons, the measurements suggest an efficient ( 90 %) export of the available net radiation by sensible and latent turbulent fluxes, thereby strongly limiting the heat available for melting ground ice. Turbulent export of heat and moisture drawn from the permeable AL contributes to the well-known insulating effect of the coarse blocky AL and partly explains the climate resiliency of rock glaciers. This self-cooling capacity is counteracted by an early snow melt-out date, exposing the low-albedo blocky surface to the intense June–July insolation and causing reduced evaporative cooling due to exacerbated moisture scarcity in the near-surface AL during dry spells. With climate change, earlier snowmelt and increased frequency, duration, and intensity of heat waves and droughts are projected. Regarding the parameterisation of the turbulent fluxes, we estimated the year-round turbulent fluxes using a modified Louis (1979) scheme. The monthly SEB is closed within 20 W m−2 except during the snowmelt months and under katabatic drainage winds in winter. Detected sensible turbulent fluxes from nocturnal ventilation processes, although a potentially important ground cooling mechanism, are within our 20 W m−2 uncertainty because nighttime wind speeds are low. Wintertime katabatic wind speeds needed to be scaled to close the SEB, which hints at the limits of parameterisations based on the Monin–Obukhov similarity theory in complex mountain terrain and katabatic drainage winds. The present work contributes to the process understanding of the SEB and climate sensitivity of coarse blocky landforms.

1 Introduction

Coarse blocky landforms such as rock glaciers, block fields, and talus slopes are covered by a thick clast-supported debris mantle typically  1–5 m thick. These “cold rocky landforms” (Brighenti et al.2021) exhibit lower ground temperatures compared with adjacent fine-grained or bedrock areas, a phenomenon referred to as “undercooling” or “algific” (Wakonigg1996; Delaloye and Lambiel2005; Millar et al.2014). This special ground thermal regime is owed to the interactive effects of energy exchange processes between the atmosphere and the ground arising from the debris mantle properties (Hoelzle et al.2001). Large blocks and typically sparse fine materials near the surface create a vast, highly connected and thus permeable pore space. Depending on the stability of the air column in the debris mantle controlled by temperature gradients, a variable part of the near-sub-surface debris mantle is ventilated and effectively participates in the turbulent heat exchange with the atmosphere. Heat and moisture stored in the ventilated near-surface sub-layer of the debris mantle are rapidly mobilised and contribute to turbulent fluxes at short ( hourly) timescales. These non-conductive sub-surface heat transfer processes related to storage and phase changes in heat, vapour, water, and ice produce the peculiar micro-climate observed in debris mantles. As a consequence, to calculate the surface energy balance, no unambiguous, clearly defined and fixed interface of energy conversion (“active surface”, Oke1987) separating radiative–convective processes in the atmosphere from (dominantly) diffusive processes in the sub-surface is available (Herz2006). These active surfaces do not necessarily coincide for different processes (radiation conversion, wind-forced and buoyancy-driven air convection, interception of precipitation, water flow) or quantities (solar radiation, momentum, water), and therefore, they vary in time and are not necessarily at ground surface. On seasonally snow-covered sites, snow controls ground–atmosphere heat exchange processes at the surface by its effects on albedo and thermal insulation (Hoelzle et al.2003; Zhang2005).

These complex interactions between debris mantle, snow cover, and atmosphere demand continuous high-resolution sub-surface measurements in addition to the “classical” above-surface weather station measurements and ground temperature. However, gathering the necessary measurements of sub-surface parameters is challenging in remote mountain terrain, and only a few comprehensive data sets beyond ground temperatures exist in mountain permafrost, exceptions being Rist and Phillips (2005) and Rist (2007). They deployed a heat flux plate, ultrasound probes, conductometer, vapour traps, and reflectometer probes to characterise the ground hydro-thermal regime of a steep, permafrost-underlain scree slope. Instead, non-conductive heat transfer processes in block fields, talus slopes, and debris cover of glaciers have been inferred from their effect on the ground thermal regime by means of (high-resolution) temperature measurements (Wakonigg1996; Wunder and Möseler1996; Humlum1997; Harris and Pedersen1998; Hoelzle et al.1999; Gorbunov et al.2004; Růžička et al.2012; Wagner et al.2019; Petersen et al.2022), in cases aided by gas/smoke tracer experiments (Popescu et al.2017) or thermal infrared imaging (Shimokawabe et al.2015).

These challenges have been described in energy balance studies conducted on the Murtèl rock glacier, situated in a cirque in the Upper Engadine (eastern Swiss Alps). In this permafrost landform, the debris mantle overlies a perennially frozen rock glacier core, is seasonally frozen, and roughly coincides with the thermally defined active layer (AL). We will use the term “coarse blocky AL” throughout the text to point to both material property and thermal state. In a micro-climatological study by Mittaz et al. (2000), significant deviations from a closed surface energy balance (SEB) of up to 78 W m−2 in winter and −130W m−2 in summer were found. These deviations exceed methodological uncertainties. The authors suggested that the apparent heat sink in summer and heat source in winter arose from processes unaccounted for in their calculations due to the lack of measurements, namely advective and convective heat transport in the coarse blocky AL. Scherler et al. (2014) addressed the seasonal SEB imbalance by adopting a volumetric energy balance approach consisting of adding a porous interfacial buffer layer able to store and release heat and including radiative and sensible turbulent heat transfer in the coarse blocky AL. The integration of additional sub-surface heat transfer mechanisms and storage components reduced the deviations to 26 W m−2 in summer and −29W m−2 in winter. The remaining seasonal deviations were interpreted as latent heat effects of freezing and thawing of ice in the AL and at the permafrost table. This agreed well with long-term melt rates derived from photogrammetric and three-dimensional borehole deformation leading to subsidence estimates ( 5 cm yr−1) (Kääb et al.1998; Müller et al.2014).

In this work, we estimate the SEB of the Murtèl rock glacier using in situ measurements on and within its ventilated coarse blocky AL. We revisit the micro-climatological studies by Mittaz et al. (2000) and Scherler et al. (2014) adding new measurements from a novel sensor array. Since the deviations present in these works were largely caused by uncertainties in the turbulent fluxes, here we focus on the parameter and parameterisations required for their calculation. We address two questions: a general one and a technical one. The general question is how large the individual surface heat fluxes on the Murtèl rock glacier are and how these are seasonally distributed. From a quantitative process understanding, we gain insight into how the insulating effect of the coarse blocky AL works and how the rock glacier, a thermally conditioned permafrost landform, responds to climate change. The technical question is which turbulent flux parameterisations and input parameters are appropriate for the study of a seasonally snow-covered ventilated landform in complex mountain terrain. Difficulties in calculating the turbulent fluxes on Murtèl arise at two spatial scales: the meso-scale relief (102103 m) and sub-landform scale (10−1101 m). The meso-scale relief (102103 m) modifies wind speeds and atmospheric stability in the Murtèl cirque, which must be accurately reflected by the flux parameterisations. In winter, persistent katabatic winds develop on the snow-covered rock faces and slopes and converge in the cirque. These katabatic jets determine the near-surface wind velocities and the vertical turbulent exchange. In summer, surrounding rock spurs weaken the regional valley wind in the sheltered Murtèl cirque. Thermals rising from the strongly heated debris surface create an unstable atmosphere with comparatively low horizontal wind speeds. We test different stability corrections commonly used for debris-covered glaciers. At sub-landform scale (10−1101 m), the coarse blocky AL is ventilated and participates in the convective exchange of momentum, heat, and moisture with the atmosphere, unless the ground is covered by a sufficiently thick snow cover. Heat and moisture can be drawn from an interfacial buffer layer coupled with the atmosphere. The active surface does not necessarily coincide with the ground or snow surface. Where should we appropriately measure the meteorological variables, namely surface temperature and surface humidity, needed to estimate the turbulent fluxes – on the ground/snow surface or at some depth? We describe these turbulent processes at the interface between atmosphere and uppermost coarse blocky AL ( 1.5 m depth) using in situ wind speed measurements and link these to measurements of eddy-covariance sensible fluxes. We then use the gained process understanding to define the appropriate “surface” temperature and humidity for the Bowen and bulk aerodynamic approaches to parameterise the turbulent fluxes.

Our work contributes to the quantitative process understanding of surface energy fluxes on a ventilated coarse blocky landform situated in a complex mountain terrain. The quantification of individual surface heat fluxes and near-surface storage terms will benefit the modelling of past and present mountain permafrost distribution and will help to anticipate the response of coarse blocky permafrost landforms to climate change.

2 Study site and past research

2.1 Murtèl rock glacier

The studied Murtèl rock glacier (WGS 84: 46°2547′′ N, 9°4915′′ E; CH1903+/LV95: 2'783'080, 1'144'820; 2620–2700 m a.s.l.; Fig. 1) is located in a north-facing periglacial area of Piz Corvatsch in the Upper Engadine (eastern Swiss Alps), a slightly continental rain-shadowed high valley (Fig. 2). Mean annual air temperature (MAAT) is −1.7 °C; mean annual precipitation is  900 mm (Scherler et al.2014). This tongue-shaped, single-unit (monomorphic, Frauenfelder and Kääb2000) active rock glacier is  250 m long and  150 m wide, surrounded by steep rock faces and directly connected to a talus slope (2700–2850 m a.s.l.). Crescent-shaped furrows ( 3–5 m deep) and ridges with steep and, in some places, near-vertical slopes dissect the slightly north-northwestward-dipping surface ( 10–12°) and create a pronounced furrow-and-ridge micro-topography in the lowermost part of the rock glacier. The snow cover is thicker and lasts longer in furrows than on ridges, influencing the ground thermal regime at small scale (Bernhard et al.1998; Keller and Gubler1993). The coarse-grained and clast-supported debris mantle is only 1–2 m thick in the colder furrows, while it is 3–5 m thick on the rest of the rock glacier. The ground ice table is accessible in a few places. Characteristic clast size ranges from 0.1 to 2 m edge length, with a few rockfall-deposited boulders of  3–5 m. Fine material ( sand) is virtually absent from most of the surface; its volume fraction increases with depth (inverse grading; Haeberli et al.2006). Rain and percolating meltwater quickly disappear, and the surface appears dry. Beneath the coarse blocky debris mantle, roughly coinciding with the thermally defined AL, lies the perennially frozen ice-supersaturated rock glacier core. Drill cores have revealed sand- and silt-bearing massive ice (3–28 m depth, ice content over 90 % by volume), although boreholes drilled within  30 m distance suggest some lateral small-scale heterogeneity (Vonder Mühll and Haeberli1990; Arenson et al.2010). Surface creep rates are  10 cm yr−1, show a coherent creep pattern, and have been slightly accelerating in the past decade (Noetzli and Pellet2023). Water emerges seasonally from several springs or seeps at the foot of the rock glacier front and flows to the rock glacier forefield of till-veneered bedrock not underlain by permafrost (lower boundary of discontinuous permafrost) (Schneider et al.2012).

The site is snow covered for 7–9 months annually with an up to 1–2 m thick snowpack. Persistent katabatic winds (Mittaz et al.2000) develop in the topographically shaded cirque (no direct insolation in November–February) and redistribute snow from the windswept ridges into the furrows, eroding the snow around large blocks (Bernhard et al.1998). These are preferential spots for snow funnels that form after the first snowfall in early winter. Oscillating airflow, resembling breathing, has been observed in these openings through the snow cover. This airflow allows for some vertical heat exchange to occur between the coarse blocky AL and the atmosphere, even during the winter months. As a result, the insulating effect of the snow cover is reduced to some degree, although the extent of this reduction has not been quantified (Bernhard et al.1998; Keller and Gubler1993; Keller1994).

Figure 1Location of Murtèl rock glacier in the Upper Engadine, a high valley in the eastern Swiss Alps. Inset map: location and extent (black rectangle) of regional map within Switzerland (source: Swiss Federal Office of Topography swisstopo).

2.2 Past micro-climatological research

The Murtèl rock glacier and the Murtèl–Chastelets periglacial debris slope have been intensely investigated since ca. 1970. One of the longest continuous mountain permafrost temperature time series worldwide (since 1987) and atmospheric measurements from an automatic weather station (AWS; since 1997) run by the Swiss Permafrost Monitoring Network (PERMOS) have turned this site into a “natural laboratory” for mountain permafrost research (summarised by Hoelzle et al.2002). Relevant statistical and process-oriented energy balance studies are the ones by Hoelzle and Haeberli (1995), Hoelzle (1996), Mittaz et al. (2000), Hoelzle et al. (2001), Stocker-Mittaz et al. (2002), Herz et al. (2003), Hoelzle et al. (2003), Hanson and Hoelzle (2004, 2005), Hoelzle and Gruber (2008), Schneider et al. (2012, 2013), Scherler et al. (2014), and Wicky and Hauck (2017, 2020), as well as at least 10 unpublished master theses.

Apart from the two studies on Murtèl by Mittaz et al. (2000) and Scherler et al. (2014) mentioned above, perhaps the most similar investigation to the present work in terms of ground properties is Isaksen et al. (2003), performed on a ventilated coarse blocky permafrost site in southern Norway (Juvvasshøe); its meso-scale landscape, however, is a windswept flat mountain top dissimilar to the Murtèl cirque.

3 Measurements and data processing

3.1 Sensor placement

A total of 50 m away from the existing Murtèl PERMOS cluster (Noetzli et al.2019; “AWS” in Fig. 2c), additional sensors were installed above the ground surface and within natural cavities of the porous coarse blocky AL in August 2020. This PERMA-XT sensor cluster comprised snow and atmospheric sensors above the ground surface, AL sensors distributed in natural cavities between the blocks (Table 1), and two automatic time-lapse cameras in the RGB colour and thermal infrared spectral range. We also used a four-component radiation sensor from the PERMOS cluster for our analysis; its specifications are reported in Scherler et al. (2014).

The above-ground sensors were located on a rock glacier ridge and included air temperature and humidity sensors, a barometer, a sonic ranger for snow height, and a sonic anemometer (CSAT; 3.87 m a.g.l.) for eddy-covariance measurements mounted on a custom-made sensor pylon. On ground level, an unheated tipping-bucket rain gauge measured liquid precipitation. Four unshielded snow thermistors at 0 (ground level), 25, 50, and 100 cm a.g.l. measured the vertical snow temperature profile. Sensor specifications are presented in Table 1.

The below-ground sensors were distributed in natural cavities at different micro-topographical positions (ridges, slopes, and furrows) around the meteo pylon within a 30 m distance and at different depths beneath ground surface. Five horizontally lying and vertically hanging thermistor strings and five thermo-anemometers (TP01/WS01) measured temperature and an airflow speed proxy at 5 and 30 min intervals, respectively. Most sub-surface sensors were concentrated in a 3 m deep and 0.5–1.5 m wide instrumented cavity (Fig. 2e). A thermistor string measured the vertical temperature profile of the cavity air (TK1/1–5), complemented by two hygrometers near the surface and in the mid-cavity (HV5; −0.7 and −2 m). Five thermistors (TK6/1–5) were drilled 5 cm into the blocks at depths corresponding to the TK1 thermistors. Three thermo-anemometers recorded wind speed at three levels: close to the surface (−0.35 m), mid-cavity (−1.5 m), and in a narrow extension at −2.1 m (not used in this study). Finally, a back-to-back pair of pyrgeometers mounted at mid-cavity level (CGR3; −1.55 m) measured the upward and downward long-wave radiation in the cavity. Detailed sensor specifications are presented in Table 1. Since accurate distances were required for the calculation of vertical gradients and fluxes, we triangulated the relative height of the sensors in the instrumented cavity with a laser distance meter and goniometer (Leica DISTO X310).

All sensors were solar-powered and wired to data loggers connected to the Internet via a mobile network. Hourly data transmission, when allowed by battery voltage, enabled timely detection of technical failures and intervention. To prevent uncontrolled power shortages during the no-insolation winter period, a protocol progressively sent power-demanding sensors into power-saving mode (50 % duty cycle or deactivated completely) as battery voltage decreased. This affected the most power-demanding sonic anemometer and, more rarely, the thermo-anemometer. Power-saving mode was active preferentially during nighttime. Such incomplete data sets are biased to daytime measurements, and thus daily average might not be representative (Sect. 3.2). For all other sensors, power from diffuse and snow-reflected radiation was sufficient for continuous operation.

Figure 2Sketch map of Murtèl–Corvatsch cirque (a) with two active rock glaciers, Murtèl and Marmugnun. (b) Photo of the above-surface PERMA-XT installations. (c) Wind patterns are seasonally varying in the sheltered cirque. Panels (d) and (e) show locations of sensors on Murtèl rock glacier and in the instrumented main cavity, respectively (Table 1). The term “cavity roof” used throughout the text refers to the uppermost  1 m of the instrumented cavity.

Table 1PERMA-XT sensor specifications.

Measurement range and accuracy by manufacturer/vendor. Specifications of PERMOS sensor available in Scherler et al. (2014) and Hoelzle et al. (2022). a CSI: Campbell Scientific, Inc. b Sensors in radiation shield RAD10E. c CSAT measures three-dimensional high-frequency sonic wind speed and derives sonic buoyancy flux. d Snow temperature setup and thermistor strings manufactured by Waljag GmbH. e Semi-quantitative airflow speed derived from measured heat flux.

Download Print Version | Download XLSX

3.2 Data processing

In the present work, we analysed data of 2 years from 1 September 2020 to 30 September 2022 except for the eddy-covariance data. The sonic anemometer (CSAT) was operational from 5 November 2020 onwards.

3.2.1 Surface radiation

All four components of the radiative heat fluxes at the surface – i.e. incoming and outgoing short-wave (S) and long-wave (L) radiation – were measured directly by the PERMOS micro-meteorological station with two back-to-back pairs of Kipp & Zonen CM3 pyranometers (0.3–3 µm) and CG3 pyrgeometers (5–50 µm) (Hoelzle et al.2022). We applied the same radiation corrections as in Hoelzle et al. (2022): (1) instrumental corrections of both pyranometers with nighttime short-wave radiation measurements (pyranometer offset correction), (2) correction for snow cover of the upward-looking pyranometer (S=S/0.86 when S>S), (3) Sicart et al. (2005) correction of the incoming long-wave radiation for interference with solar radiation, and (4) rejection of incoming short-wave radiation at low sun angles (if albedo exceeds unity).

3.2.2 Eddy-covariance data

We estimated the sensible turbulent flux QH using the eddy-covariance method, based on the covariance of vertical wind speed fluctuations w and temperature fluctuations T, measured by a CSI CSAT (CSI CSAT3B manual2020) and averaged over 30 min (Foken2017; see Sect. 4.2.2):

(1) Q H eddy = ρ a C p w T ,

where (ρaCp) refers to the volumetric heat capacity. We calculated buoyancy fluxes from the eddy raw data (10 Hz sampling rate) using Campbell Scientific's software EasyFlux™ and following the conventional pre-processing chain (Radić et al.2017; Fitzpatrick et al.2017; Steiner et al.2018). Processing steps were trend and outlier removal, despiking (Vickers and Mahrt1997; Foken et al.2012), coordinate rotation using the double rotation method (Tanner and Thurtell1969), and spectral corrections (Moncrieff et al.1997; Massman2000).

Due to the lack of fast-response moisture measurements (no high-frequency gas analyser), the eddy-covariance method yielded the sonic buoyancy flux wTsnc, which, for |QH|>|QLE|, is close to but not equal to the sensible heat flux wT (Liu and Foken2001). We assess the discrepancy using the Bowen ratio method modified by Liu and Foken (2001) based on Schotanus et al. (1983) (“SND correction”):

(2) w T = w T snc 1 + 0.51 T a C p L v Bo - 1 ,

where Tsnc is the sonic temperature, Ta is the air temperature (averaged over 30 min), and Bo refers to the Bowen ratio (Sect. 4.2.2). Cp and Lv are constants defined in Sect. 4.2.2.

3.2.3 Snow height and snow water equivalent

We measured snow height with a sonic distance sensor (Table 1) on a rock glacier ridge in addition to the nearby PERMOS snow height measurements. Raw measurements were compensated for variations in the speed of sound with air temperature, following the manufacturer's guidelines (CSI SR50A manual2020).

We used the semi-empirical ΔSNOW model (Winkler et al.2021) to convert the measured snow height to snow water equivalent (SWE). This parsimonious model requires snow height and its temporal changes as the sole input. Calibration data used to develop their model were gathered in the Swiss and Austrian Alps in climatic regions similar to the Engadine.

3.2.4 Precipitation

We measured liquid precipitation with an unheated tipping-bucket rain gauge (Table 1) mounted on a rock glacier ridge.

Snow or sleet (a mix of snow and rain) can fall in any month of the year. The distinction from rain is important because of the latent heat of melting Lm that ice carries and by far exceeds the sensible heat of rainwater (cwΔT). In the comparatively dry climate of the Engadine, wet-bulb temperature can be a better discriminator than dry-bulb temperature (Froidurot et al.2014). Wet-bulb temperature Twb [K] is defined in Stull (2017):

(3) e - e ( T wb ) = - γ ( T - T wb ) [ Pa ] ,

where γ=(cpP)/(0.622Lv) [Pa K−1] is the psychrometric “constant”. Equation (3) is solved iteratively. Based on hourly images from the automatic time-lapse camera, we defined a conservative wet-bulb temperature threshold of 2 °C that separates snow or sleet (Twb<2 °C) from rain.

3.2.5 Sub-surface long-wave radiation

The in-cavity net long-wave radiation QCGR3rad was calculated from the upwards Lal and downwards Lal long-wave radiation components measured with a back-to-back pyrgeometer pair installed in the instrumented cavity at a depth of 1.55 m beneath ground level (Table 1):

(4) Q CGR 3 rad = ( L al - L al ) ,

which is not to be confused with the long-wave radiation L measured at 2 m a.g.l. (Sect. 4.2.1).

We corrected raw outputs Lraw/ from the two pyrgeometers in the instrumented cavity by accounting for the long-wave radiation emitted by the instruments themselves (Kipp & Zonen CGR3 manual2014):

(5) L al / = L raw / + σ T CGR3 4 ,

where TCGR3 refers to the pyrgeometer housing temperature. Large (>0.5 °C) or rapid changes in housing temperature differences between the two back-to-back mounted pyrgeometers indicated dust or water deposition on the upward-facing pyrgeometer window. Such disturbed measurements appeared in the high-resolution (10 min) data but did not significantly affect the daily net long-wave radiation balance in the sheltered cavity.

3.2.6 Sub-surface airflow speed

We refer to the sub-surface “wind” in the cavity as “airflow” to differentiate it from the atmospheric wind. We used the Hukseflux WS01/TP01 sensor to perform airflow speed measurements in the AL cavities. This sensor included a heated foil that measures a cooling rate expressed as a convective heat transfer coefficient hWS [Wm-2K-1] related to airflow speed (Hukseflux WS01 manual2006). We did not convert hWS to airflow speed but process this variable as a semi-quantitative indicator assuming uhWS, sufficient for our purpose. WS01/TP01 do not resolve the direction of the airflow (hence the term “speed” instead of “velocity”). Deposition and evaporation of liquid water can disturb measurements (hWS increase), as revealed by wrapping the heated foils in moist tissues. WS01 measurements during precipitation events were filtered out. Repeated zero-point checks were performed throughout the snow-free season by enclosing the heated foil in small, dry plastic bags for a few hours, ensuring stagnant conditions with zero airflow speed. Neither drift nor temperature dependency beyond measurement uncertainty was detected.

4 Surface energy balance calculation

4.1 Energy balance of the Murtèl near-surface AL

The point-scale surface energy balance at seasonally snow-covered sites accounts for net radiation Q [W m−2], composed of short- and long-wave radiation components; turbulent fluxes, composed of sensible heat QH [W m−2] and latent heat QLE [W m−2]; melt energy of snow QM at the surface [W m−2]; energy from precipitation QP [W m−2]; and heat flux QG [W m−2] into and from the ground when snow-free, replaced by conductive heat flux across the snow cover QS [W m−2] when snow covered:

(6) S + S + L + L net radiation Q + Q H + Q LE turbulent fluxes + Q M + Q P + Q G / S = 0 [ W m - 2 ] .

Fluxes are counted as positive if these provide energy to the reference surface, i.e. the terrain or snow surface. Unlike Scherler et al. (2014), we consider the reference surface as an infinitely thin skin layer without storage. Fluxes must be balanced at all times (Eq. 6). Ground QG and snow QS heat fluxes measured beneath the surface are extrapolated to the reference surface by means of the calorimetric correction.

4.2 Flux parameterisations

We estimated the fluxes (terms in Eq. 6) as follows: all radiative fluxes were derived from on-site measurements, i.e. net radiation Q both above the surface (PERMOS data) and within the AL (PERMA-XT data); turbulent fluxes QH and QLE were estimated using the Bowen energy balance and the bulk aerodynamic methods and directly from eddy-covariance measurements. Snowmelt QM and precipitation QP heat fluxes were estimated using the calorimetric method from SWE estimates and on-site measured rainfall rates, respectively. The heat flux in the snowpack QS is the calorimetrically corrected conductive heat flux at the base of the snowpack. The ground heat flux QG (flux in the near-surface AL) was estimated analogously to the calorimetrically corrected net long-wave radiation measured in situ in the AL at 1.5 m depth.

4.2.1 Surface radiative fluxes

The net radiation Q is the sum of the measured and corrected radiation components (Sect. 3.2.1):

(7) Q = S + S S + L + L L ,

where S is the outgoing short-wave radiation, S is the incoming short-wave radiation, and net short-wave radiation S is the sum of both; correspondingly, L represents the outgoing long-wave radiation, L is the incoming long-wave radiation, and net long-wave radiation L is the sum of both.

4.2.2 Surface turbulent heat fluxes QH, QLE

We estimated the turbulent flux using three different methods: (1) the Bowen energy balance method (Bowen1926), (2) the bulk aerodynamic method (Mittaz et al.2000; Hoelzle et al.2022), and (3) directly with the eddy-covariance method from CSAT measurements (in the lack of a fast-response vapour analyser, only the sensible turbulent flux QH is estimated; Sect. 3.2.2).

Bowen energy balance method

The Bowen ratio Bo is defined as the ratio of sensible to latent heat flux (Bowen1926; Ohmura1982; Oke1987) and reflects the partitioning of the turbulent fluxes into sensible and latent components:

(8) Bo := Q H Q LE = K h K w C p Δ T L v Δ q C p ( T a - T s ) L v ( q a - q s ) ,

where Kh represents the eddy diffusivities for sensible heat and Kw is water vapour. Invoking the similarity principle and assuming KhKw (Ohmura1982), the Bowen ratio can be calculated from the isobaric specific heat capacity of air (Cp=Cpd(1+0.84q) with Cpd=1005Jkg-1K-1) (Reid and Brock2010; Hoelzle et al.2022), the latent heat of vaporisation Lv (if Ts≥0 °C) or sublimation Ls (if Ts<0 °C) (2.48×106J kg−1, 2.83×106J kg−1), and the gradients of temperature ΔT and specific humidity Δq (specified below). The sensible and latent turbulent fluxes are then expressed as a function of the available energy Q+QM+QG/S (QP considered negligible; Sect. 5.3.4):

(9) Q H Bo = Q + Q M + Q G / S 1 + 1 / Bo , Q LE Bo = Q + Q M + Q G / S 1 + Bo .

The ground and snow surface temperature Ts [°C] was calculated from the measured long-wave radiation components L and the Stefan–Boltzmann law via (Oke1987)

(10) ε σ T s 4 = L - ( 1 - ε ) L ,

where ε is the surface emissivity; σ represents the Stefan–Boltzmann constant (σ=5.670×10-8Wm-2K-4); and L and L represent the measured outgoing and incoming long-wave radiation, respectively. We took the emissivity value of ε=0.96 from Scherler et al. (2014).

The surface specific humidity qs [g g−1] is either the specific humidity of the air in the near-surface coarse blocky AL qsa (measured) or the saturated snow surface qss(Ts,Pa), depending on whether or not the snow cover is thick enough to suppress convective air exchange between the AL and the atmosphere (Sect. 6.3.1).

(11) q s = q sa , for h S SEB h ^ S crit (snow-free or open snow cover) q ss , for h S > SEB h ^ S crit (decoupling snow cover)

The critical snow height SEBh^Scrit is reached when the AL is so strongly decoupled from the atmosphere that the convective fluxes across the snow cover are no longer detectable by our SEB estimations. We determine the snow height SEBh^Scrit using the AL air temperature, specific humidity, and airflow speeds (Sect. 6.3.1). The specific humidity of the air qa was calculated from the measured air temperature Ta and relative humidity at the measurement level zm of 2 m a.g.l.

Bulk aerodynamic method

In the bulk parameterisation, sensible QHbulk and latent QLEbulk turbulent fluxes are driven by the gradients of temperature ΔT=Ta-Ts [K, °C] and specific humidity Δq=qa-qs [g g−1], respectively, and horizontal wind speed u [m s−1]:


where ρa [kg m−3] represents air density. The flux–gradient relationship is modified by atmospheric stability, accounting for enhanced turbulent fluxes in an unstable atmosphere and suppressed turbulent fluxes in a stable atmosphere, by means of the bulk exchange factors Cbt (bulk turbulent heat and vapour transfer coefficients) [unitless, –]. Different formulations of the stability functions exist. Here, we use the modified Louis (1979) scheme (Song1998), motivated by its use in the GEOtop model, a distributed hydrological model designed for complex terrain (Rigon et al.2006; Endrizzi et al.2014), and other studies (e.g. Essery and Etchevers2004). We compare the modified Louis (1979) scheme with the iterative Monin–Obukhov scheme and the widely used Businger–Dyer scheme (Oke1987). The latter has been used in the previous SEB studies on Murtèl by Mittaz et al. (2000) and Scherler et al. (2014) and in many SEB estimates on debris-covered glaciers (Brock et al.2010; Reid and Brock2010; Reid et al.2012; Steiner et al.2018, 2021), despite discrepancies (overestimated fluxes) at strongly unstable atmosphere noted by Steiner et al. (2018). The issue with the Businger–Dyer parameterisation is that near-surface wind speeds on topographically sheltered rough terrain tend to be lower for a given atmospheric stability compared with the conditions under which the empirical Businger–Dyer parameterisation was originally developed (flatland in Kansas, USA; Businger et al.1971; Haugen et al.1971). This approach is therefore problematic in complex terrain, including our study site. Finally, we test the parameterisation without stability correction (“bulk c0”), which corresponds to the special case of a neutral atmosphere. Detailed explanations of the parameterisation schemes of the turbulent heat transfers are described in the appendix (Appendix A).

Accurate estimates of turbulent fluxes rely on representative values for the roughness lengths for momentum z0m, heat z0h, and moisture z0q (Steiner et al.2018; Smeets et al.1999; Smith2014). We calculated the roughness length for momentum z0m from the eddy-covariance data following Conway and Cullen (2013) and Fitzpatrick et al. (2017):

(14) z 0 m = ( z u - d ) exp - k u h u - ψ m z u - d L ,

where uh [m s−1] represents the mean horizontal wind speed, u [m s−1] is the friction velocity, k [0.40] is the von Kármán constant, zu [m] is the measurement height of wind speed, L [m] is the Obukhov length, ψ is the integrated stability function (Appendix A), and d is the zero-plane displacement height. We set the (unknown) zero-plane displacement height d to zero (Miles et al.2017). The sensor height zu(t) is the variable distance above the ground or snow surface: zu(t):=zCSAT-hS(t). The recommended filtering is to use only uh>3m s−1 under near-neutral conditions |zu/L|<0.05 (Radić et al.2017; Nicholson and Stiperski2020). We compared our eddy-covariance-derived roughness length for momentum z0m with the aerodynamically derived values reported in Mittaz et al. (2000). For simplicity, the scalar lengths for heat z0h and humidity z0q are often considered equal to the momentum roughness length z0m (Sicart et al.2005; Reid and Brock2010) or 1–3 orders of magnitudes smaller (Smeets and van den Broeke2008). Here, we could not independently estimate the scalar roughness lengths z0h and z0q due to the lack of humidity-corrected sonic temperature or high-frequency humidity measurements. We assumed an equal roughness length for heat and moisture, z0h=z0q, and used the unknown ratio R^z0:=z0h/z0m=z0q/z0m as a calibration parameter.

Finally, we compared the sensible turbulent flux with the katabatic model of Oerlemans and Grisogono (2002) specifically developed for conditions of katabatic or nocturnal drainage winds. This is different from all the above-mentioned parameterisations as this predicts a quadratic rather than near-linear relation between QH and the driving temperature difference ΔT. Since we lacked vertical profile observations of potential temperature required for this parameterisation, we could not test this bulk method. Instead, we checked the validity of the quadratic QHΔT relation on Murtèl (cf. Radić et al.2017) and show wind profile measurements collected by Mittaz et al. (2000) in the years 1997–2000 (Appendix D).

4.2.3 Snowmelt heat flux QM, snow heat flux QS, and snowpack sensible heat storage HS

The snowmelt heat flux QM is the heat flux consumed by the melting snowpack. We estimated QM from the daily change in SWE [kg m−2] derived from the snow height data hS (sonic ranger data) with the semi-empirical parsimonious ΔSNOW model (Winkler et al.2021):

(15) Q M L m Δ ( SWE ) Δ t ,  if  T b = 0 ° C ,

where Lm [336 kJ kg−1] represents the latent heat of fusion for ice. We considered runoff-generating snowmelt to occur in spring when temperature at the base of the snowpack Tb reaches the melting point and ignored melting–refreezing events within the snowpack.

The snow heat flux QS is the conductive heat flux across the single-layer snowpack. We calculated Q̃S in the lowermost 25 cm of the snowpack, ignoring transient effects (Mittaz et al.2000):

(16) Q ̃ S = - k S d T d z - k S T s - T b h S - k S T 25 - T 0 Δ z ,  if  h S > 30 cm and 0 otherwise ,

where dT/dz is a linearised temperature gradient in the lowermost 25 cm of the snowpack. To ensure that both thermistors are snow covered, the minimum snow height hS required is 30 cm. We related the snow thermal conductivity kS to the snow density via the empirical equation kS:=kS(ρS)=2.93(10-6ρS2+10-2) developed for Murtèl (Keller and Gubler1993; Hoelzle et al.2022), where we estimated the bulk density of the snowpack with the ΔSNOW model via ρS=SWE/hS (Winkler et al.2021).

Q̃S was calculated near the base of the snowpack instead of at the snow reference surface to which Eq. (6) refers. With this approach, we used the more stable snow density and thermal conductivity near the snowpack base, which is less affected by compaction than the near-surface layers that receive fresh snow. The heat flux Q̃S was extrapolated to the snow surface by adding the sensible heat storage changes in the snowpack ΔHS above zS=25 cm (“calorimetric correction”) (Foken2008; Liebethal et al.2005; Boike et al.2003). We estimated the changes in cold content ΔHS as

(17) Δ H S c S ( m S / A ) Δ T S / Δ t c S SWE Δ T S / Δ t ,

where (mS/A) represents the mass of the snowpack per area, i.e. the SWE [kg m−2]; cS [2.050kJkg-1K-1] is the specific heat capacity of snow or ice at 0°C; and 〈ΔTS [°C] represents the layer-averaged snow temperature changes. The calorimetrically corrected snow heat flux QS was then calculated as

(18) Q S = Q ̃ S + Δ H S .

4.2.4 Precipitation heat flux QP

The rainfall heat flux QPr was estimated via (Sakai et al.2004; Reid and Brock2010)

(19) Q Pr = C w r ( T P - T s ) ,

where Cw=ρwcw (4.18 MJm-3K-1) is the water volumetric heat capacity and r [m3m-2s-1] is the rainfall rate intercepted at the surface. Precipitation temperature TP was approximated using the wet-bulb temperature Twb, calculated from air temperature and relative humidity (Eq. 3). We assumed precipitation in the form of rain if Twb≥2 °C. Water contributions from upslope flowing onto the rock glacier and liquid precipitation falling into the snowpack were not accounted for.

The heat flux QPs due to rapid melting of sleet or shallow summertime snow (Twb<2 °C, hS<10 cm) was roughly estimated as

(20) Q Ps L m ρ w r s ,

where rs [m3m-2s-1] represents the amount of solid precipitation that melts in the pluviometer per time period Δt [s]. Our measurement setup was not designed to accurately record the precipitation rate during mixed rain- and snowfall or the fraction of ice crystals and liquid water in the total precipitation. The heat flux QPs is intended as an order-of-magnitude estimate.

4.2.5 Near-surface ground heat flux QG and sensible heat storage Halθ

Analogously to the calorimetrically corrected snow heat flux QS, the ground heat flux QG is the sum of the measured net long-wave radiation in the instrumented cavity QCGR3rad (Eq. 4) and the sensible heat storage changes ΔHalθ of the virtual layer between the ground surface and the depth of the pyrgeometer:

(21) Q G = Q CGR 3 rad + Δ H al θ .

Rock mass subject to changing temperatures constitutes a heat source or sink (Scherler et al.2014). The sensible heat storage change in the (dry) blocks ΔHalθ is proportional to the rate of change in rock temperature, assuming a constant volumetric heat capacity {(1-ϕal)ρrcr}/t=0. In the lack of rock temperatures over the 2-year time period analysed in this work, we use in-cavity air temperatures Ta averaged over a thermal adjustment timescale Δtr (1 d) as a surrogate (Sect. 6.1.2; Appendix C). The rate of sensible heat storage and release of the blocks was then approximately calculated as (Liebethal and Foken2007; Ochsner et al.2007; Brock et al.2010)

(22) Δ H al θ = - z s 0 t { ( 1 - ϕ al ) ρ r c r T r ( z ) } d z ( 1 - ϕ al ) ρ r c r Δ t r i { T a ( z i , t + Δ t r ) - T a ( z i , t ) } Δ z i ,

where zs=155 cm is the distance from the pyrgeometer pair to the ground surface (reference level); ρr is the rock density (2690 kg m−3) (Corvatsch granodiorite; Schneider2014); cr is the specific heat capacity (790 Jkg-1K-1); ϕal=0.4 (Scherler et al.2014) is the AL porosity; and Tr(z) and Ta(z) are the vertical rock and in-cavity air temperature profile [°C], respectively. In the discretised formulation, temperatures T(zi) are layer-wise averages in the ith layer with thickness Δzi (denoted by “〈⋅〉”) derived from the thermistor string TK1/1 and the radiometric surface temperature Ts.

5 Results

5.1 Meteorological conditions

The weather in each season differed markedly between the 2 years analysed in the present work (2020–2022; Fig. 3). The winter 2020–2021 was colder than the 2021–2022 one (November–April: average temperature: −6.2 °C vs. −5.3 °C; minimum daily average temperature: −16.5 °C vs. −15.1 °C) and richer in terms of snow amount (November–April: average snow height measured on a wind-swept ridge: 76 cm vs. 54 cm) and duration (early onset of snow cover: 5 October vs. 3 November; later melt-out: mid-June vs. mid-May). Summer 2021 was cool and wet compared with the hot and dry summer 2022; temperatures were lower (July–August: average: 6.9 °C vs. 9.3 °C) with frequent passage of synoptic fronts, often bringing cold air (≤3 °C; minimum daily average temperature: 0.7 °C vs. 5.6 °C) and mixed precipitation (sleet). Snowfall occurred in a few days throughout the summer and melted within hours. A few snow patches survived over the summer after melt-out of the winter snowpack in mid-June. In contrast, the hot and dry summer 2022 was marked by three heat waves (in June, July, and August) and daily minimum temperatures not below 5 °C. Several dry spells occurred during this season; the longest one was an 11 d long dry spell within the 5–19 July heat wave. Almost no precipitation was recorded between 20 June and 1 August despite some convective precipitation events recorded on the nearby Piz Corvatsch cable car station (MétéoSuisse station at 3294 m a.s.l.). Discharge data of the rock glacier outflow (own measurements, not shown), camera images, and field observations (fresh debris flow deposits, flooding of furrows) revealed rainwater funnelled onto the rock glacier. Data gaps for this period were filled with MétéoSuisse precipitation data from the nearby station “Piz Corvatsch”.

Wind speeds (measured by PERMOS) in the sheltered Murtèl cirque were generally low (hourly means: 1–3 m s−1; peak wind speed ∼5m s−1; Fig. 2c). The wind pattern was often marked by a strong diurnal cycle. Peak wind speed was reached during the night in winter (strong katabatic wind blowing downslope from southeast) and in the afternoon in summer (regional valley wind known as the “Maloja wind” blowing from west-northwest–west-southwest overruled a local anabatic wind). Summer nights were calm or with weak katabatic winds (wind speed: ≤2m s−1 from west–southeast).

Figure 3Meteorological conditions. (a) Air temperature (daily mean) and precipitation (daily sum). (b) Snow height and SWE. Rain and sleet (mixed precipitation) separated based on a wet-bulb temperature threshold of 2 °C. Precipitation data at Piz Corvatsch from MétéoSuisse.


5.2 Snow and ground thermo-hygric conditions

During the snow-rich winter 2020–2021, a thick and insulating snow cover (exceeding 70 cm) sealed the cavity from the atmosphere. The cavity air was kept isothermal (within ±0.2 °C; Fig. B1; Appendix B) and isohume at saturation at a much higher moisture level than that in the (colder) atmosphere (Fig. B2; Appendix B). Virtually no dry and cold air from the atmosphere was mixed into the closed cavity system. A stable winter equilibrium temperature (WEqT; bottom temperature of snow cover (BTS, Haeberli1973; Haeberli and Patzelt1982) was reached in March 2021 (-4.59±0.05 °C). In contrast, during the snow-poor winter 2021–2022, the rapidly fluctuating temperature and humidity (relative and specific) indicated a connection across the thin snow cover and an exchange with the dry, cold air from the atmosphere. Following snow melt-out (end of zero curtain) and re-connection with the atmosphere, the cavity air began to warm and “desaturate” from the surface downwards (July 2021; June 2022).

In summer, the near-surface cavity roof was generally warmer and more humid compared with both the atmosphere (despite slightly lower relative humidity) and the deeper cavity (specific humidity profiles are shown in Appendix Fig. B3). Daily average specific humidity in the cavity roof and that in the atmosphere were strongly correlated (R2=0.868). The near-surface cavity consistently presented a moisture surplus in relation to that of the atmosphere throughout the year. This surplus persisted even during dry spells. Summertime in-cavity temperature gradients were more stable the higher the surface temperature was. In summer 2021, frequent passages of cold fronts with rapid atmospheric cooling destabilised the in-cavity air column by reducing the temperature gradients. In summer 2022, dry spells impacted the sub-surface moisture conditions down to 2 m within 3–6 d after the last precipitation event: water infiltrated within minutes, near-surface relative humidity started to decrease within hours, and the cavity air drying front (evaporation front at the isohume of rH = 100 %) receded to greater depths within days. At a depth of 2 m, saturation was lost  5 d after the last precipitation event. Consequently, the in-cavity humidity gradients reversed, indicating a switch from downwards to upwards humidity transport. This was most pronounced in July 2022 during the intense mid-July dry spell accompanying a heat wave (Appendix Figs. B3, B2). In contrast, near the surface, the measured gradient in specific humidity between the near-surface AL and the atmosphere and, therefore, QLE remained largely constant regardless of the different weather conditions and the humidity gradient within the AL (Sect. 6.3.3).

Figure 4Indicators of AL–atmosphere coupling. (a) Relationship between the daily temperature amplitude in the cavity roof (normalised by the 2 m temperature amplitude) and the corresponding snow depth hS measured on a wind-swept rock glacier ridge (PERMA-XT station) and a broad flat area (PERMOS station). (b) Normalised airflow speed (u/umax) at different locations as a proxy for sub-surface ventilation and convective coupling. Measurements below the level of detection (LoD) were considered zero. Airflow speed decreases with increasing snow height relative to the maximum speeds under snow-free conditions. Onset of decoupling varies with sensor location. WS/5 is beneath a large wind-exposed block on a ridge where snow funnels remain open longer than those on flat terrain (WS/3, WS/4).


We determine the snow height necessary to close the snow cover and to decouple the AL from the atmosphere with the sub-surface AL air temperature, specific humidity (not shown, correlated with air temperature), and airflow speeds (Fig. 4). With increasing snow depth, the differences in air temperature and specific humidity between the cavity and the atmosphere increased (Appendix Figs. B1, B2), the correlation between in-cavity and atmospheric signal was lost, and rapid (hourly–daily) fluctuations in the in-cavity temperature and specific humidity weakened (Fig. 4a). Normalised near-surface AL temperature amplitudes indicate the degree of convective coupling: amplitudes similar to that in the atmosphere mean coupled and strongly attenuated amplitudes mean decoupled (Fig. 4a; additional winter 2022–2023 data shown). The decoupling proceeded gradually with snow depth (sketched by the schematic envelope), and the thresholds are only approximate values. Also, the sub-surface airflow speed was attenuated gradually according to the micro-topographic setting that controls the local snow depth (Fig. 4b): on gently sloping and less rough terrain, the vertical connection between the coarse blocky AL and the atmosphere was reduced at snow heights of 5–20 cm and lost at ∼50–60 cm of snow (WS/3 and WS/4); on wind-swept ridges with wind erosion, a much thicker snow cover was necessary to shut down the last snow funnel (Bernhard et al.1998; Keller and Gubler1993) (∼80 cm at WS/5 at the upwind side of a big block on a ridge).

5.3 Surface energy balance

Monthly SEB (Eq. 6) is dominated by the short-wave S and long-wave L radiation components, followed by the sensible QH and latent QLE turbulent heat fluxes and the ground heat flux QG (Fig. 5; Table 2). In winter, turbulent fluxes compensate for the energy lost by net radiation. In summer, turbulent fluxes export roughly 90 % of the available net radiation. Snowmelt QM absorbs practically the entire net radiation (QM-Q) in the respective melt-out month (June 2021; May 2022), and the sensible and latent turbulent fluxes either are then small or roughly cancel each other (Fig. 6). The ground heat flux QG is downwards (negative) during snowmelt (infiltration of meltwater into the frozen coarse blocky AL releases latent heat) and net negative during the thaw season. The sensible rain heat flux QP is a negligible SEB component. In short, SL|QH|>|QLE|>|QG|QP. The parameters and constants used in the calculations are tabulated in Appendix E (Table E1).

Figure 5(a) Monthly energy balance components. Turbulent fluxes calculated using the modified Louis (1979) parameterisation (cL). (b) Daily and monthly surface albedo A as a snow-cover indicator.


Table 2Season-averaged heat fluxes and heat-flux ratios.

Turbulent fluxes calculated using the modified Louis (1979) (cL) bulk parameterisation.

Download Print Version | Download XLSX

5.3.1 Surface radiation

Short- and long-wave radiation are by far the largest SEB components (Fig. 5a). Consequently, the well-known effect of the snow cover on albedo is the single largest control of net radiation and hence of the entire SEB (Fig. 5b).

In mid-winter, the meso-scale relief controls the solar radiation budget. The shaded north-facing Murtèl cirque receives no direct insolation from November to February. Depending on cloud cover and amount of incoming diffuse or terrain-reflected short-wave radiation, the net radiation is negative and dominated by the long-wave radiation budget (LS<0). Net radiation Q was less negative in the cloudy and precipitation-rich winter 2020–2021 (December–February; average: −19W m−2) than in the sunny and dry winter 2021–2022 (−40W m−2) due to greater incoming long-wave radiation L emitted by the clouds (233 vs. 215 W m−2). A case in point is the exceptionally warm and sunny November 2020, which received in total 33 W m−2 less short- and long-wave radiation than November 2021. The monthly mean radiation deficits are up to −53W m−2, which, together with the steep snow-covered slopes, is a favourable setting for strong katabatic winds.

5.3.2 Turbulent heat fluxes

The Bowen ratio partitions the available energy from net radiation, the ground heat flux, and the snowmelt heat flux into sensible and latent turbulent fluxes. Throughout the seasons of both years, the heat fluxes were as follows (Fig. 6): during wintertime, from November to March, without direct insolation in the shaded cirque, the energy loss of 20–50 W m−2 due to long-wave emission and terrain shading was largely compensated by the sensible heat flux QH. Latent heat flux QLE at the cold snow surface was small, with mostly resublimation (moist winter 2020–2021) and sublimation fluxes (dry winter 2021–2022) within -6±10W m−2 (Fig. 6a). With the latent heat flux being small and in a variable direction, the computed Bowen ratio showed a large scatter (|Bo|>1; Fig. 6b). In the drier winter 2021–2022, QH compensated for the heat lost by sublimation in addition to the radiative heat loss. During spring, from April to May, before the snowmelt period, the magnitude of the fluxes remained similar but the direction reversed. In summer, after complete snow melt-out, from June/July to October, the turbulent fluxes were larger to compensate for the large radiation surplus. In the more rainy summer 2021, sensible and latent heat exports were similarly large (Bo=1.25±2.71), whereas in the drier summer 2022, heat export was more dominated by sensible flux (Bo=1.48±2.47). Towards early autumn, with decreasing temperatures, the sensible flux lost importance relative to the latent flux (Bo decreased and approached zero). In September/October, with mixed precipitation (sleet) and first snow falls at still warm conditions near 0 °C, the latent heat export by sublimation from the snow cover intermittently dominated over the sensible heat uptake (Bo between 0 and −1). With further cooling, the moisture supply became limited again, and sensible heat uptake took over to offset the increasing radiation deficit (Bo<-1).

We compare the estimates of turbulent fluxes (monthly averages) with the net radiation (Fig. 6) and among each other (Fig. 7). In the snow-free summer months, the Louis (1979) bulk fluxes tend to slightly underestimate the fluxes, while the Businger–Dyer fluxes (cB&D) clearly overestimate them. The measured eddy-covariance flux (buoyancy flux) falls short of closing the SEB, roughly by 50 %–75 % depending on how large the (unmeasured) eddy latent flux would be. The imbalance (“missing” flux to close the SEB; Foken2008) is, in any case, substantial.

Figure 6(a) Sensible and latent (dotted boxes) turbulent surface fluxes. Available energy from net radiation Q, snowmelt QM, and ground heat flux QG (hatched boxes) is partitioned into sensible QHBo and latent QLEBo turbulent fluxes according to the Bowen ratio. The Bowen SEB is closed by design. Bulk fluxes differ according to the stability function and do not necessarily close the SEB. (b) Bowen ratio (Eq. 8). November–March: large scatter in winter because QLEBo is small and changes direction (unstable calculation). June/July–October: Bo 0.5–2.5, reflecting the wet and cool summer 2021 and the dry and hot summer 2022 (Fig. 3). Monthly average based on 16–31 valid values per month.


Figure 7Comparison of turbulent flux estimates from different parameterisations (daily and monthly averages; see Fig. 6 for time series). Comparison metrics: root mean square error (RMSE [W m−2]), mean bias error (MBE [W m−2]), and Pearson correlation coefficient (r [–]) on daily estimates.


5.3.3 Snowmelt heat flux

The latent heat of the melting snowpack QM largely consumes the available net radiation Q in the respective snowmelt months; the snowmelt heat flux is roughly as large as typical summer sensible turbulent heat fluxes QH (Fig. 5).

5.3.4 Precipitation heat flux

The sensible rain heat flux QPr exerts a negligible cooling effect of -2W m−2 at the daily timescale (not shown). Short but intense rainfall (thunderstorms) with fluxes up to −100W m−2 in 10 min is averaged out because such high-precipitation events are short. The heat fluxes QPs arising from mixed precipitation or a shallow snow cover that melts within hours (“summer snow”; 10–30 W m−2 daily average) are similar to QG and not negligible. Such events typically occur in early autumn (September–October) but also occurred throughout the wet and cool summer 2021 (June–October). Their timing often coincides with episodes of rapid ground cooling (Fig. 8).

5.3.5 Snow and ground heat fluxes

Snow heat flux QS is in the range of 2 to −4W m−2 and ground heat flux QG in the range of 30 to −20W m−2 at the daily timescale (Fig. 8). The snowpack is a thermal insulator. Within-snowpack conductive fluxes are 2–3 orders of magnitude smaller than typical SEB components. The conductive snow heat flux is negligible compared with the overall SEB and across-snowpack convective fluxes (snow funnels; Figs. 8, 4). The simplifications from the single-layer snowpack that ignore vertically varying snow density do not detract from this finding.

Heat storage changes in the snowpack (Fig. 8a) and uppermost coarse blocky AL (Fig. 8b) are asymmetric (Guodong et al.2007) in opposite directions. In the snowpack, downwards heat transfer or warming (storage gain) is more intense (larger maxima) but less frequent, likely due to non-conductive fluxes from refreezing meltwater (synchronised with warm spells Ta>0 °C or rapid warming to the zero curtain). In the near-surface AL, upwards heat transfer or cooling (storage loss) is more intense but less frequent. Different weather conditions in the two summers analysed are reflected by QG: the passage of cold fronts in summer 2021 led to more frequent convective cooling, often accompanied and enhanced by sleet; conversely, ground warming is pronounced during dry spells in summer 2022. Similarly for the two winters, virtually no storage changes (QG≈0) occurred in the snow-rich winter 2020–2021, whereas temperature and sensible heat storage fluctuations continued in the snow-poor winter 2021–2022, albeit strongly attenuated compared with those in the snow-free season (Fig. 4; Sect. 6.3.1). The non-zero QG in winter 2021–2022 beneath the thin snow cover reflects convective heat exchange across the snow funnels since the conductive snow heat flux QS is too small to account for QG.

Figure 8(a) Snow heat flux QS and changes in cold content of the snowpack ΔHS. (b) Ground heat flux QG, net long-wave radiation in the instrumented cavity QCGR3rad, and sensible heat storage changes ΔHalθ (uppermost 1.5 m of the coarse blocky AL). Beneath the insulating snow cover in winter 2020–2021, the ground thermal regime is stable, and QG≈0W m−2. In contrast, the rapid fluctuations in QG during the snow-poor winter 2021–2022 are faster than those during the conduction time, indicating convective processes as a driver. Heat input was larger in summer 2022 than in summer 2021.


5.4 Aerodynamic roughness length

Valid values of roughness length for momentum z0m (Eq. 14) scatter over 2 orders of magnitude in the range of 10−2 to 100 m. Few (2.1 %) data points met the quality criteria u>0.1m s−1 under near-neutral conditions |ζ|<0.05 or u>3.0m s−1 (Radić et al.2017; Nicholson and Stiperski2020) (Fig. 9a). The bin-wise average (median) is 0.19 m (0.23 m) for snow-free conditions (hS<10 cm) and slightly decreases with increasing snow height up to 0.07 m (0.08 m) at 90–100 cm of snow (Fig. 9b). Averages were calculated from the average of the logarithmised values (Radić et al.2017). Snow heights between 10 and 50 cm or exceeding 100 cm are rarely observed in the 2-year data set, hence the observation gap. z0m is itself a function of sensor distance and hence of snow height hS (Eq. 14). To control for possible spurious correlation in the z0mhS relation, we eliminated the confounding variable hS and tested with the constant measurement height zCSAT=3.78 m. This did not significantly affect the z0mhS relation.

Figure 9Calculated momentum roughness length z0m vs. (a) horizontal wind speed uh and (b) snow height hS. Of all values (grey dots), 2.1 % met the quality criteria (blue dots). Momentum roughness length decreases from 0.19 to 0.07 m (×10-0.43) with snow height increasing from 0 to 100 cm.


6 Discussion

6.1 Surface fluxes and uncertainties

6.1.1 Fluxes

Monthly SEBs are closed within the calculation uncertainties of ≤20W m−2 (Fig. 10, blue bars; Table 3). This represents a substantial improvement compared with Mittaz et al. (2000) and Scherler et al. (2014) and is due to the novel sensor array used in the present work. The largest SEB imbalances (deviation from closure) occur during the transition between seasons, when ground thermal conditions linger around 0 °C, during extreme meteorological conditions, and in mid-winter (December–March). These are the snowmelt months (June 2021, May 2022), early autumn (September), and the July 2022 heat wave and accompanying dry spell that strongly impacted the ground thermal and moisture regime (Figs. B1B2). Larger deviations that occurred in the snow-rich mid-winter 2020–2021 are reduced by a “katabatic wind correction” for the variable anemometer height above the snow surface (Sect. 6.2).

Figure 10SEB imbalance (monthly averages). Turbulent fluxes calculated using modified Louis (1979) (cL) bulk parameterisation. Correcting excessively high wind speed measured in the low-level katabatic jet at thick snow cover improves the SEB in the snow-rich winter 2020–2021 (katabatic wind correction; Sect. 6.2, Eq. 23).


6.1.2 Parameter sensitivity

Uncertainties arise from terrain or snow cover variability across the rock glacier, spatially variable properties of the coarse blocky AL (e.g. emissivity, porosity, intrinsic permeability), and instrumental measurement errors, among other factors. We assessed the impact of the largest sources of parameter uncertainty of the heat fluxes (Table 3; sensor accuracy from Scherler et al.2014; Hoelzle et al.2022). We estimate our overall accuracy as ∼20W m−2 at the daily timescale due to uncertainties in QG (uncertainties in sensible heat storage changes ΔHalθ) and QH (uncertainty in surface temperature Ts propagated from the emissivity ε; Eq. 10). These uncertainties play an important role in this study. The SEB uncertainties determine which processes are included in the SEB estimates and which are not, complementing measurement-driven criteria. We consider insignificant (i.e. within the uncertainty) the processes with associated fluxes smaller than our 20 W m−2 uncertainty threshold. The processes considered insignificant in the context of the SEB are not necessarily insignificant at depth. In fact, in the AL, all daily-averaged sub-surface heat fluxes are within 20 W m−2 (Fig. 8b). We apply this criterion on the nocturnal Balch ventilation (Sect. 6.3.2) and the decoupling snow height SEBh^Scrit (Sect. 6.3.3).

Insolation differences due to the local shading effect of the coarse terrain surface, the terrain slope effects, or instrumental effects (cosine response) might lead to differences of up to 30 W m−2 (daily average) in summer (Table 3) (Otto et al.2012; Scherler et al.2014; Fitzpatrick et al.2017; Kraaijenbrink et al.2018). Uncertainties in the outgoing long-wave radiation L of ±10 % corresponding to a radiometric surface temperature difference of 7 °C might arise from patchy snow cover. The lateral advective heat transport altering the boundary layer characteristics (Essery et al.2006; Mott et al.2013, 2015) could explain the large SEB imbalance (Fig. 5) of roughly 100–150 W m−2 during the melt-out phase. Smaller but still significant Ts differences of 2 °C can occur due to local shading (micro-topography) or uncertainty in emissivity (e.g. variable snow emissivity; Steiner et al.2018) and lead to considerable QH uncertainties of up to 50 W m−2. Surface temperature differences from different emissivity values (ε=1 instead of 0.96 from Scherler et al. (2014); Eq. 10) are within −1.2 and +0.25 °C. The deviations tend to be largest on hot clear-sky days with little incoming long-wave radiation L (Eq. 10) and hence might translate into considerable uncertainties in QH precisely when these are largest.

Temperature and specific humidity uncertainties driving bulk fluxes (Eqs. 8, 1213) lead to considerable QH and QLE uncertainties. The Essery and Etchevers (2004) parameterisation is moderately sensitive to the spatially variable wind field which might arise from the micro-topography, for example wind sheltering in the furrows. The momentum roughness length z0m is varied by a factor of 100.5 (between 0.09 and 0.45 m), which reflects the standard deviation of the measurements (Fig. 9). The roughness lengths are perhaps among the most critical parameters for the estimation of turbulent fluxes when using the bulk approach (besides the emissivity ε): sensitive yet hard to constrain on rough and complex mountainous terrain (Steiner et al.2018; Rounce et al.2015; Giese et al.2020). Equivalence between momentum z0m and scalar roughness lengths for heat z0h or moisture z0q would lead to prohibitively large deviations and can be excluded. As the rough terrain turns sensor height above the surface into a somewhat vague parameter, we vary it by a typical block edge length of 0.5 m. The arising QH deviations of 10–20 W m−2 are similar to other parameter uncertainties. Air temperature and specific humidity measurements from the PERMOS and PERMA-XT stations, located within 50 m, differ by 0.6±1.0 °C and 0.004±0.24g kg−1, respectively. The calculated fluxes are similar: QH is within ±5W m−2 and QLE within 2 W m−2. Maximum deviations are temporarily up to 10–15 W m−2 in winter and during the snowmelt period. The PERMOS and PERMA-XT station data yield flux estimates indistinguishable within their uncertainty.

Ground heat flux QG primarily reflects the rate of temperature change (RTC) of the near-surface AL and is weakly sensitive to the pyrgeometer flux measurements QCGR3rad. The reason behind is that QG (Eq. 21) is dominated by the sensible heat storage changes ΔHalθ of the 1.5 m thick blocky layer above the long-wave radiation measurement (Liebethal et al.2005). Consequently, rather large uncertainties in QG (20 W m−2) come from uncertainties in the storage changes δ(ΔHalθ) from two factors, namely the porosity ϕal and the thermal adjustment time Δtr (time lags). First, a high porosity limits the heat storage capacity that scales with 1−ϕal. We tested a plausible range of δϕal±0.2 and obtained uncertainties up to 20 W m−2. Porosity might laterally vary closest to the surface, precisely where daily temperature amplitudes are largest. Second, the assumption of similar air and rock temperature profiles Ta(z)Tr(z) – local thermal equilibrium (LTE) assumption (Nield and Bejan2017)) – becomes problematic in the roof of the ventilated cavity, where conditions are highly transient (convective heat transfer and effect of insolation). The entire rock mass might not adjust to the rapid temperature fluctuations of +2 to −3 °C from day to day. The chosen thermal relaxation time Δtr of 1 d is a minimum duration (Appendix C). We assess the influence of longer adjustment times Δtr by comparing 1 d with 5 d running averages. The differences are similar to the uncertainties related to porosity. We conclude that, due to the large and rapid heat turnover in the ventilated near-surface coarse blocky AL under highly transient conditions, the ground heat flux QG is arguably the least-constrained flux. This might not be a surprising finding on a landform that does not present a clearly defined surface.

Finally, we let SEBh^Scrit=80 cm to simulate an open snow cover in early–mid-winter 2020–2021 (October–February; Fig. 3) and throughout the snow-poor 2021–2022. Due to the large humidity differences across the snow cover (Fig. B2), this resulted in a massive (up to 102 fold) increase in QLE. The closing of the snow cover is a key control factor for the SEB and ground thermal regimes (Eq. 11; Sect. 6.3.1) (Hanson and Hoelzle2004).

Table 3Instrumental sensitivity: uncertainty due to parameter values and meteorological variables.

Uncertainty of daily average fluxes. QH and QLE using the modified Louis (1979) parameterisation (cL).

Download Print Version | Download XLSX

In the sensitivity analysis (Table 3), input parameters and meteorological variables are varied independently from each other and based on likely maximum measurement errors. However, the meso-scale relief and local sub-surface processes like ventilation might add some systematic biases that exceed instrumental errors. We will explore parameterisation uncertainties (stability corrections) in Sect. 6.2 and uncertainties in the meteorological input variables in Sect. 6.3.

6.2 Meso-scale terrain effects on turbulent flux parameterisations

6.2.1 Katabatic wind

The interaction of the steep snow-covered wintertime terrain with a negative radiation balance induces downslope katabatic winds that govern the near-surface wind field. This, in turn, affects the calculation of the bulk turbulent fluxes that require a representative wind speed as input. Our initially calculated wintertime turbulent fluxes were “too large” by 10–35 W m−2 on monthly average (Fig. 10), especially during the overcast snow-rich winter 2020–2021, while the radiative cooling and the forcing temperature deficit ΔT were weaker than in the snow-poor winter 2021–2022 (Fig. B1; Appendix B; Eqs. 1213). This result suggests that measured wind speeds were off in relation to the snow height and not to the temperature deficit. In fact, wind tower measurements on Murtèl performed by Stocker-Mittaz (2002) showed strong and persistent katabatic winds in winter, with a wind speed maximum a few metres above the surface weakly correlated to snow height (Appendix D; see Fig. 2c). The growing snow cover “shifted” the high-velocity region of the low-level katabatic jet to the level of the wind sensor, causing an apparent wind speed increase. We compensate the wind speed us measured at variable height above the snow cover (zshS) for the snow height with a simple power-law relation:

(23) u = u s ( z u / ( z s - h S ) ) - m ^ ,

where m^ is the exponent. Note that the usual correction of the sensor height above the variable snow surface further increased rather than decreased the turbulent fluxes. Since our measurements cannot resolve the near-surface wind speed profile of the low-level katabatic jet, the intention of this ad hoc “katabatic wind correction” is not to accurately describe the wintertime wind profile but rather to pragmatically render our input data amenable to the flux–gradient parameterisations based on the Monin–Obukhov theory. We use m^ as a calibration parameter (denoted by the circumflex). A literature value of m^=0.6 for stable conditions (Arya2001) reduced the wind speed sufficiently to yield turbulent fluxes that close the SEB within our uncertainty threshold (Fig. 10). Although the Monin–Obukhov theory is not strictly valid under such conditions (Grisogono et al.2007), our consistent flux estimates based on a scaled wind speed are in line with Endrizzi et al. (2014) and Denby and Greuell (2000), who argue that the bulk method still provides reasonable flux estimates when measured close to the surface, as was the case here (measurements within 2 m above the variable snow cover). This illustrates the importance of accurate wind speeds for the calculation of turbulent fluxes. We suggest an alternative solution in Sect. 6.2.2.

6.2.2 Comparison of turbulent flux parameterisations

From the comparison of the measured eddy-covariance fluxes and the calculated Bowen and bulk fluxes (daily averages), we reach the following conclusions.

Overall, the modified Louis (1979) parameterisation (cL) seems the best choice to parameterise the turbulent fluxes on Murtèl for three reasons. First, as expected, the modified Louis (1979) cL bulk fluxes were near-identical to the iteratively calculated Monin–Obukhov fluxes (Fig. 7d; r=0.996) but at less computational cost. Second, neither filtering nor post-processing was necessary because virtually all estimates met the quality criteria. Third, these show the least deviations from the Bowen fluxes (Fig. 7e), which we consider the benchmark as the Bowen SEB is closed by design. The Bowen fluxes offer reasonable daily estimates based on minimal measurement requirements, namely temperature and specific humidity profiles. Wind speed is not required. The minimum time resolution is daily scale. Hourly values either are numerically unstable (small ΔT/Δq, especially over a “warm” spring–early-summer snow cover) or do not account for the systematic diurnal wind speed variations. In summer, for example, hourly sensible Bowen fluxes are overestimated during the calm nights and compensated by excessive daytime fluxes.

The eddy-covariance flux systematically underestimates the sensible turbulent flux by a factor of ∼2, leading to a large imbalance (underclosure) in the eddy SEB. The error from the lacking SND correction (Sect. 3.2.2) cannot account for the eddy SEB imbalance. The ratio between the eddy sensible flux and the (sonic) buoyancy flux is between 0.93 and 1.05 (10 % and 90 % quantile, respectively; Eq. 2) on a daily average. Since both the air and the ground surface are most often far from saturation, the eddy sensible flux error due to the lacking high-frequency gas analyser is less than 10 % on 617 out of 692 d (89.2 %) with valid Bowen ratios. Larger deviations occur on the remaining 75 d with |QLE|>|QH| (Bowen ratio within −0.6 and 0.5), typically during snowmelt periods (“warm” saturated snow surface), cloudy and rainy summer days (more frequently in summer 2021 than 2022), or during the first snowfall in autumn. Applying the SND correction showed little added value. Nonetheless, the summertime eddy-covariance fluxes correlate reasonably well with the Bowen sensible turbulent fluxes (Fig. 7a) and the modified Louis (1979) bulk fluxes (Fig. 7b) or the driving temperature gradient (Fig. 11a), and therefore random instrumental or data processing errors represent an unlikely explanation. We hypothesise that some hidden systematic reason causes the flux underestimation and the large eddy SEB imbalance, e.g. secondary circulation (Foken2008) of the anabatic afternoon winds.

The empirical Businger–Dyer (cB&D) formulation, developed over flat terrain, overestimates the turbulent fluxes and deviates more strongly for larger negative fluxes (“banana” shape; Fig. 7f), despite extensive filtering. Even the bulk parameterisation without stability correction (c0) outperformed cB&D for the summer fluxes under unstable atmospheric conditions (Fig. 7c). This finding agrees well with Steiner et al. (2018) on the Lirung debris-covered glacier that shares many topo-climatic features with the Murtèl rock glacier (unstable atmosphere over a strongly heated debris surface with anabatic valley winds).

Finally, we describe the relation between the wintertime sensible Bowen fluxes with the driving temperature gradient ΔT (“surface temperature deficit”) as quadratic (Fig. 11b). A quadratic functional QHΔT relation is unique to the katabatic model of Oerlemans and Grisogono (2002). Note that the Bowen fluxes are independent of wind speed (Eq. 8). Together with the ad hoc “katabatic correction”, our data show that the turbulent fluxes in a snow-covered cirque with katabatic winds might be better parameterised by a katabatic model than by the common Monin–Obukhov bulk method, as also found by Radić et al. (2017) on a sloped glacier surface and discussed by Grisogono et al. (2007) in the context of atmospheric boundary layer modelling.

Figure 11Comparison of the turbulent sensible flux estimates with the driving temperature difference ΔT (a–c) and atmospheric wind speed u (d–f) motivated by Eq. (12) (daily and monthly averages). Measured winter eddy fluxes are weakly sensitive to local ΔT (a) and u (d). Sensible Bowen flux shows a seasonally differing QHΔT relation: quasi-linear in summer and quadratic in winter (b) (with best-fit lines). Relation to wind speed u also differs seasonally (e), with a clearer increase in QH with wind speed in summer. Note that wind speed is not an input variable for the Bowen parameterisation. The seasonally differing QHΔT relation reflects the seasonally differing atmospheric stability and wind conditions: stronger dependency on ΔT in the unstable summer atmosphere. The quadratic relation found for the wintertime fluxes suggests that katabatic nocturnal drainage winds govern the turbulent heat transfer (Oerlemans and Grisogono2002). Different parameterisations of the bulk fluxes (c) show different sensitivities to ΔT: c0 and cB&D are over-sensitive to ΔT in winter and summer, respectively.


6.3 Micro-scale landform effects on turbulent flux parameters

The “surface” temperature Ts and “surface” specific humidity qs are key inputs for the Bowen ratio (Eq. 8) and bulk methods to estimate turbulent fluxes (Eqs. 1213). However, heat and moisture can be drawn from the ventilated AL beneath the surface, provided the snow cover is sufficiently thin to allow convective exchange (hS<SEBh^Scrit). Although an “open” snow cover that allows vertical convective exchange between the AL and the atmosphere (Sect. 6.3.1) shapes the wintertime ground thermal regime up to a snow height of ∼70 cm, we will argue in Sect. 6.3.2 and 6.3.3 that already at more than 20 cm of snow convective exchange across the snow cover can be ignored in the context of our SEB.

6.3.1 Snow cover and AL–atmosphere coupling

The snow cover controls the coupling between the coarse blocky AL and the atmosphere. The turbulent fluxes draw heat and moisture from the AL as long as the snow cover is “open” via snow funnels, a phenomenon widely observed on coarse blocky landforms (Sawada et al.2003; Delaloye et al.2003; Delaloye and Lambiel2005; Morard et al.2008; Schneider2014; Kellerer‐Pirklbauer et al.2015; Popescu et al.2017). A thickening snow cover gradually suppresses the vertical convective coupling between the AL and the atmosphere, as more snow funnels close (sketched by the schematic envelopes in Fig. 4). A snow height of ∼10 cm begins to decouple the coarse blocky AL from the atmosphere (“semi-closed” in Fig. 4), but much more snow (∼70 cm) is necessary to achieve the insulating effect of the snow cover on the ground thermal regime (“closed/insulating”). Then, the snow cover is thick enough to suppress convective air exchange, rapid sub-daily fluctuations are strongly attenuated, and large temperature and moisture gradients to the outside air build up (Haberkorn et al.2015), indicating that the AL–atmosphere coupling is weak. Such a high value is typical for a terrain as rough and blocky as on Murtèl, agreeing with, for example, Hanson and Hoelzle (2004) and Herz (2006). The insulating effect of a thick snow cover was already observed decades ago and used to indirectly map the permafrost distribution using the bottom temperature of snow cover (BTS) method (Haeberli1973; Haeberli and Patzelt1982).

As regards the SEB, a snow cover as thin as SEBh^Scrit=20 cm causes a decoupling that is strong enough for our SEB to be significant. This is shown in the snow-poor winter 2021–2022 with snow depths between 40 and 70 cm, when the ground heat flux QG kept fluctuating but remained below the 20 W m−2 uncertainty threshold (Table 3). The convective flux across the snow cover cannot deviate strongly from QG because it was the single largest heat flux to supply/extract the heat for the AL sensible heat storage changes (small conductive snow heat flux, QS±2W m−2, Ts<0 °C, and no latent effects from snowmelt in that period). We take this threshold as the critical snow thickness (SEBh^Scrit; Eq. 11). We emphasise that the 20 cm threshold is an operational definition for the “decoupled” snow cover in the context of our SEB and its uncertainties.

6.3.2 “Surface” temperature and near-surface ventilation

From Sect. 6.3.1 it follows that for the decoupling or insulating snow cover (hSSEBh^Scrit), the reference surface coincides with the snow surface. The radiometric surface temperature represents the surface temperature Ts (Eqs. 8, 12). Whenever the snow cover is open or under snow-free conditions (hS<SEBh^Scrit), the radiometric surface temperature is the input meteorological variable Ts (Eq. 12). The surface is where the solar radiation is intercepted and transformed to thermal energy, hence the main source of sensible heat for QH. The large gradient (Ta-Ts)/zT drives the sensible turbulent flux QH (Eqs. 8, 12; Fig. 12a) and the wind-forced ventilation (Fig. 12b) (Panz2008; Stocker-Mittaz2002). This is shown by the measured eddy-covariance flux (Fig. 12c) that largely follows the difference between air and radiometric surface temperatures (Fig. 12a, red area), as do the local (anabatic) wind in the atmosphere above the surface and the airflow speed in the near-surface AL (Fig. 12b). During daytime with strong heating on the low-albedo surface, the blocky surface exceeds air temperatures by >15 °C.

However, we found evidence of near-surface ventilation that cannot be parameterised by the 2 m air temperature and the radiometric ground surface temperature: nocturnal Balch ventilation, a nighttime cooling process owed to the interplay of the air permeability, and the thermal inertia (heat retention) of the coarse blocky AL that is most apparent during fair-weather summer days with clear nights. Upwards turbulent heat export is protracted into the evening hours, long after sunset (S; Fig. 12c), as shown by the measured in-cavity airflow speed (WS/3; Fig. 12b) and the eddy-covariance flux (Fig. 12c). Furthermore, the measured eddy-covariance flux remained upwards-directed in the early morning hours, when the terrain surface has radiatively cooled to air temperature TsTa or even during reversals Ts<Ta in clear-sky nights (that frequently occur in the Engadine: e.g. 24 out of 31 nights in August 2022). The calculated bulk fluxes driven by 2 m air temperature and radiometric ground surface temperature (rGST; red area in Fig. 12a; Eq. 12) are then close to zero (or even downwards) and do not fully capture this nocturnal drawing of heat from the near-surface AL. The large thermal inertia of the rock mass and the protection from long-wave radiative cooling stabilise the sub-surface air temperatures and keep the air in the roof of the ventilated cavity during the nights warmer than that in the atmosphere outside and that on the ground surface (Fig. 12a; TK1/1TK6/1>TaTs by ∼4 °C). In the nighttime, the locally stably stratified air in the uppermost cavity roof is hence warmer and more unstable than that in the atmosphere (non-local static stability, Stull1991). The high permeability of the coarse blocky AL allows this air to escape upwards into the atmosphere or, equivalently, allows the colder outside air to sink into the cavity, thereby exporting heat (Balch ventilation; Balch1900; Millar et al.2014). Panz (2008) describes an occasional nocturnal sub-surface ventilation on Murtèl when Ts<Ta. Also Yao et al. (2014) interpreted nocturnal near-surface ventilation on the debris-covered Koxkar glacier (Xinjiang, China) from eddy-covariance and temperature data.

We assess this uncertainty by using the cavity roof temperature TK1/1 instead of the conventional 2 m air temperature as input Ta for the bulk fluxes (Eq. 12) and compare the nighttime fluxes (when TK1/1 >Ts). Due to the low wind speeds at night, the estimated nocturnal Balch fluxes are small despite the appreciable temperature gradients (10–20 W m−2), maximally 10 W m−2 larger than the “conventional” bulk flux and within our SEB 20 W m−2 uncertainty (Table 3). Since the short-wave radiative forcing S and the wind speed covary in phase, the “conventional” radiometric ground surface temperature and 2 m air temperature measurements are sufficiently adequate to parameterise the turbulent fluxes when sub-daily resolution is not required – especially considering that the chosen parameterisation (stability function) shows a much larger effect on the calculated turbulent flux (Fig. 12c; compare cL with cB&D). We reach a conclusion equivalent to that in Sect. 6.3.1 for the critical snow height: the nocturnal convective processes in summer are within the uncertainties of the daily-averaged SEB but might exert an important cooling effect on the sub-surface energy balance.

Figure 12Hourly averages of (a) temperatures, (b) atmospheric wind speed, (c) near-surface AL airflow speed, and (d) sensible turbulent fluxes and net short-wave radiation S during a summer fair-weather period (15–19 July 2022). The cavity air is stably stratified (TK1/1 > TK1/3), but the air in the cavity roof is warmer and unstable compared with the outside air. During clear summer nights, the atmospheric air and the blocky surface cool down more strongly than the air in the near-surface coarse blocky AL that receives heat from the blocks (air temperature TK1/1 approaches the rock temperature TK6/1 in the night). Conversely to the calculated bulk fluxes, the measured eddy flux decays slowly in the evening–midnight (18–24 h) and remains upwards (negative) despite the small temperature gradient in the early morning (3–7 h, a). The turbulent sensible flux draws heat from the near-surface AL that is nocturnally warmer than Ta and thus unstable. (b, c) Atmospheric wind and AL airflow speeds co-vary in phase (but at 10× smaller magnitudes), follow the net radiation with some delay, and are closely related to the turbulent fluxes (d). The specific humidity shows small diurnal oscillations; the associated QLE varies less than QH in absolute terms.


6.3.3 “Surface” humidity and evaporation

From Sect. 6.3.1 it follows that for the decoupling or insulating snow cover (hSSEBh^Scrit), the snow surface humidity at saturation qss represents the surface specific humidity qs (Eqs. 8, 13) . Whenever the snow cover is open or under snow-free conditions (hS<SEBh^Scrit), the humidity measurement in the cavity roof qsa is the input meteorological variable qs (Eq. 11).

Perhaps contrary to the impression of a dry-looking blocky surface, the estimated summertime latent turbulent flux QLE is  30 W m−2, corresponding to an evaporation rate of  1 mmw.e.d-1, even during the severe July 2022 dry spell. This value is in the range of evaporation rates reported for debris-covered glaciers of 0.6–2.8 mmw.e.d-1 (Collier et al.2014; Yao et al.2014; Steiner et al.2018). However, note that our evaporative flux estimate QLE might be an upper bound, in particular during the July 2022 dry spell. The moisture source for evaporative QLE is in the AL, not at the surface (except during precipitation events). However, our bulk parameterisation (Cbt; Eq. 13) ignores the additional resistance to vapour transport imposed by the blocky layer between the specific humidity measurement in the cavity roof and the atmosphere (qsa is measured 0.7 m beneath the terrain surface; Fig. 2e). The neglected resistance to vapour transport in the coarse blocky AL during moisture-limited evaporation stages, when moisture is drawn from deeper levels and governed by vapour diffusion through the porous AL (Pérez1998; Collier et al.2014), might lead to an overestimation of QLE. The longer the dry spell lasts, the deeper in the AL the moisture is drawn from (Sect. 5.2), and the more important this effect of vapour transport resistance becomes. This likely explains why the specific humidity in the cavity roof is almost always (down to sub-hourly timescales) higher than that in the atmosphere despite the relentless mixing by ventilation (moisture surplus relative to the atmosphere; Sect. 5.2), and it possibly accounts for the negative July 2022 SEB imbalance (Fig. 10). The upwards vapour transfer in the comparatively large and strongly ventilated cavity during the dry spell might be overly efficient compared with that in the surrounding AL. This might have led to a non-representatively high specific humidity qsa in the cavity roof and an overestimated evaporative flux QLE.

6.3.4 Aerodynamic roughness lengths

Our mean (median) values of the filtered aerodynamic momentum roughness lengths z0m of 19 to 7 cm (23–8 cm) agree with the values of Mittaz et al. (2000) and Stocker-Mittaz (2002) of 18 and 7 cm for the snow-free and snow-covered Murtèl surface, respectively (Fig. 9). The scatter range of  2 orders of magnitude is similar to that in other studies (Brock et al.2006; Sicart et al.2014; Radić et al.2017; Fitzpatrick et al.2017, 2019). The absolute values are at the upper limit of what has been typically found on debris-covered glaciers (Miles et al.2017) or on Juvvasshøe by Isaksen et al. (2003) (5 cm), which is plausible given the rough terrain of the Murtèl rock glacier.

The calculated roughness length decreases slightly with increasing snow height by 0.1 cm per cm of snow (Fig. 9b). A similar but much stronger relation between snow height and roughness lengths was found on the Haut Glacier d'Arolla, where Brock et al. (2006) found z0m to decrease by 2 orders of magnitude with snow heights up to 3 m. This represents a much stronger relation than that on Murtèl, where z0m decreases by 0.5× over the observed range of snow thickness. With a maximum snow height of 120 cm (PERMA-XT measurement; Fig. 3) or 200 cm (PERMOS) during the measurement period, the snow cover is thin compared with the terrain roughness (edge length of blocks  50 cm) and undulations. Both the comparatively high roughness length and the weak sensitivity to snow height might suggest that the aerodynamic roughness is largely controlled by the furrow-and-ridge micro-topography that is not smoothed out by the snow cover or by the few largest blocks that stick out of the snow cover rather than by the average block size. In forests, for comparison, the tallest trees of the canopy have a disproportionally large influence on the aerodynamic roughness (Foken2017).

For the bulk parameterisation, we linearly interpolate log 10{z0m} with snow height. Due to the lack of accurate humidity-corrected sonic temperature measurements, we cannot calculate z0h or z0q. We used the unknown ratio R^z0:=z0h/z0m=z0q/z0m as a calibration parameter and found R^z0=2×10-2. While the approximation z0hz0q seems applicable, z0mz0h is not compatible with our parameterisation (Table 3). This agrees with previous studies on debris-covered glaciers (Steiner et al.2018, 2021) and hummocky ice surfaces (Smeets and van den Broeke2008; Sicart et al.2014).

6.4 Synthesis

The thick coarse blocky AL strongly insulates the underlying permafrost body. We quantified this well-known effect on Murtèl: during the thaw seasons 2021 and 2022, roughly 90 % of the net radiation Q was exported by the turbulent fluxes and not available to melt ground ice (Fig. 13). The ratio between received surface net radiation Q to downwards-transmitted heat flux Qr=QCGR3rad is Q:Qr9:1. This follows from the ratio of the measured in-cavity net long-wave radiation QCGR3rad to the surface net radiation Q and the relation QCGR3rad=-0.12Q+5.6W m−2 (R2=0.680). Our measurements corroborate the relation between surface net radiation and ground heat flux that has been found by Hoelzle et al. (2022) using PERMOS data from 1997–2019. Export of the received net radiation is predominantly by sensible fluxes QH (50 %–70 %) and secondarily by latent fluxes QLE (30 %–50 %; Table 2). Surface albedo (spring–early summer snow cover) and sub-surface thermal and moisture regime of the thick coarse blocky AL control the energy partitioning at/near the surface and the rock glacier's efficiency/ability to export the heat supplied by the surface net radiation Q, the primary heat input (Fig. 5; Table 2).

More heat was transferred into the ground (and more ground ice observed to melt) during the hot and dry summer 2022 than during the cool and wet summer 2021. Two sets of conditions enhanced the downwards heat transfer QG in the hot and dry summer 2022 as measured by the sub-surface pyrgeometer, in particular during dry spells and heat waves: first, less snow in winter 2021–2022 and a warm 2022 spring resulted in an early snowmelt in May, 1 month earlier than in summer 2021, and an early start of the thaw season. This exposed the dark low-albedo blocky surface to the strong insolation and the June heat wave, resulting in a rapid rise in ground temperatures and hence an increase in sensible heat storage and downwards ground heat flux (Fig. 8; see Hoelzle et al.2022, for a 20-year perspective). Second, the 11 d dry spell in July 2022 exhausted the near-surface moisture stores in the coarse blocky AL as indicated by the specific humidity deficit (qsa-qsa) in the cavity roof (Appendix Fig. B2). Lack of moisture limited the evaporative cooling QLE and countered the rock glacier's ability to export the heat (Fig. 8). The moisture storage capacity of the coarse blocky Murtèl AL is limited because the near-surface AL presents a rapid drainage and little fine material (silt) to hold water. During dry spells, the evaporation front recedes quickly to greater depths.

Figure 13Relation between measured 10 d averaged surface net radiation Q and measured 10 d averaged in-cavity long-wave radiation QCGR3rad (of the consecutive 10 d window) summarises the relation between Q, QH, and QLE. The in-cavity net long-wave radiation QCGR3rad represents the radiative downward heat transfer Qr towards the ground ice table (note that upward heat transfer is primarily by convection and not represented by the radiation measurements; QCGR3rad=Qr0W m−2 during cooling episodes). Of the available net radiation Q, 90 % is exported into the atmosphere by QH+QLE. Only 10 % is transmitted deeper into the ground and available to warm the AL and to melt ground ice. Latent turbulent heat export is less efficient during dry spells, when the coarse blocky AL dries out and evaporation becomes increasingly moisture-limited (outliers from June, July 2022). The relations between in-cavity long-wave radiation and air temperature, surface temperature (L), or QH are qualitatively similar.


Efficient evaporative cooling relies on frequent precipitation events to resupply, as occurred in summer 2021. With climate change, both factors – early start of the thaw season and hot and dry weather spells in summer – are projected to worsen. A shift to earlier snowmelt and an increase in the number of snow-free days has already been observed in the Swiss Alps (Hoelzle et al.2022) and is projected to continue with climate change (Scherler et al.2013; Marmy et al.2016). Also, the frequency, duration, and intensity of heat waves and dry spells are likely to increase. Changes in the SEB of the thermally conditioned rock glaciers and other mountain permafrost landforms entail changes in the ground thermal regime and ground ice content first and, ultimately, morphological changes: thawing, melting of ground ice, and degradation.

7 Conclusions

We estimated the year-round surface energy balance (SEB) of the seasonally snow-covered ventilated coarse blocky active layer (AL) of the active Murtèl rock glacier that is situated in a cirque in the Upper Engadine (eastern Swiss Alps). The meso-scale landscape produces seasonally contrasting atmospheric conditions of downslope katabatic jets in winter and a strongly unstable atmosphere in summer. At landform scale, the ventilated near-surface AL acts as a buffer layer where heat and moisture transfer is coupled to the atmosphere, unless sealed by a thick snow cover. Based on a novel sensor array located above the ground surface and in the AL that expands an on-site automatic weather station from the Swiss Permafrost Monitoring Network (PERMOS), we were able to improve previous SEB calculations on Murtèl by Mittaz et al. (2000) and Scherler et al. (2014). The measurement period is from September 2020 to September 2022. Our monthly SEB imbalances are within 20 W m−2 except during the snowmelt months. The main findings concern (i) the climate resilience of the Murtèl rock glacier and (ii) technical aspects of the turbulent flux parameterisations.

  • i.

    Two crucial factors for climate resilience of the Murtèl rock glacier are the insulating high-albedo snow cover and a near-complete energy turnover during the snow-free thaw season. The two meteorologically contrasting years studied in this work with the cool and wet summer 2021 and the hot and dry summer 2022 were traced into the SEB and ground thermal and moisture regime. First, an early snowmelt in 2022 prolonged the thaw season, exposed the dark blocky surface to the intense July insolation, and increased AL temperature gradients and the downward heat flux. Second, about 90 % of the received surface net radiation was exported and only 10 % effectively transferred towards the ground ice table and available to melt ground ice. Heat export occurs predominantly via sensible turbulent flux and secondarily via the latent turbulent flux from evaporation. The degree of energy turnover and the turbulent flux partitioning is co-controlled by the availability of moisture for evaporation. Dry spells and heat waves counter the rock glacier's ability to export heat by limiting evaporative fluxes leaving the rapidly drying coarse blocky AL, since the moisture storage capacity is limited. Both trends – earlier snowmelt and more heat waves and dry spells – are projected to continue with climate change. This potentially renders coarse blocky landforms vulnerable to heat waves and dry spells.

  • ii.

    With our in-mountain permafrost unprecedentedly comprehensive in situ measurements of eddy-covariance flux, liquid precipitation, snow height, AL temperature and humidity profiles, sub-surface airflow speeds, and sub-surface long-wave radiation, we tested different bulk parameterisations and constrained their input parameters.

    • We parameterised the year-round turbulent fluxes using the modified Louis (1979) scheme (Song1998) despite seasonally contrasting atmospheric stability and wind profiles. Wintertime katabatic wind speeds needed to be scaled to close the SEB, which hints at the limits of parameterisations based on the Monin–Obukhov similarity theory in complex mountain terrain and katabatic drainage winds.

    • Sensible and latent partitioning of the turbulent fluxes using the simple Bowen ratio approach agreed with the bulk fluxes on monthly to daily resolution. Given its parsimony and independence from atmospheric stability and wind speed, the Bowen approach represents a valuable and robust tool to estimate the turbulent fluxes in complex terrain with a strongly heterogeneous wind field, provided sub-daily resolution is not required.

    • Since daily oscillations of wind speed and insolation are nearly in phase, the sensible turbulent flux QH is driven by the gradient between radiometric ground surface temperature Ts and air temperature Ta. The heat flux QH from nocturnal ventilation of the permeable coarse blocky AL is small due to low wind speeds at night. Using a radiometric surface temperature is convenient for modelling QH with remotely sensed data.

    • During the thaw season, the evaporative turbulent flux QLE is 30 %–50 % of the sensible turbulent flux as a function of moisture availability in the near-surface AL. During dry spells, the near-surface moisture stores are exhausted within days, and moisture supply to the surface is limited by the convective vapour transport from the deeper AL.

    • Measured eddy-covariance fluxes were systematically too small to close the SEB. Still, the eddy-derived momentum roughness length z0m corroborates previous aerodynamic estimates from Mittaz et al. (2000). z0m varies seasonally between 7 and 20 cm as a function of snow height that smooths the landscape.

Appendix A: Turbulent flux parameterisations

Different formulations of the stability functions exist. Here, we compare the widely used Businger–Dyer relations with the formulation from Louis (1979) and the iterative Monin–Obukhov scheme.

A1 Businger and Dyer parameterisation

In the Businger and Dyer parameterisation (denoted by “cB&D”), the bulk transfer coefficient Cbt is expressed as

(A1) C bt = k 2 ln z u z 0 m ln z T / q z 0 h / q ( ϕ m ϕ h / q ) - 1 ,

where the Businger–Dyer non-dimensional stability functions for heat and water vapour are expressed as a function of the bulk Richardson number Rib (Eq. A3) (Brutsaert1982; Oke1987),

(A2) ( ϕ m ϕ h ) - 1 = ( ϕ m ϕ q ) - 1 = ( 1 - 5 Ri b ) 2 , for 0.0 Ri b 0.2 (stable case) ( 1 - 16 Ri b ) 0.75 , for Ri b < 0.0 (unstable case)

that corrects for non-neutral stability of the near-surface atmosphere (atmospheric stratification). These are unity in the case of neutrally stable atmosphere (Rib=0), denoted by “c0”. Used constants and parameters are the von Kármán constant k (0.40); the roughness lengths for momentum z0m, heat z0h, and water vapour z0q; and the (constant) measurement heights for wind speed zu [m], temperature zT, and humidity zq above the ground surface. In this work, the variable sensor height above the snow surface is accounted for by correcting the wind speed u rather than simply subtracting the snow height from the sensor distance above the ground for site-specific reasons discussed in Sect. 6.2.

The Businger–Dyer and the modified Louis (1979) parameterisations characterise atmospheric stability with the bulk Richardson number Rib defined by

(A3) Ri b := g T T a - T s z T - z 0 h u z u - z 0 m 2 ,

with the average temperature T:=(Ta+Ts)/2 [K].

However, this widely used stability function (Reid and Brock2010; Brock et al.2010; Mittaz et al.2000; Favier et al.2004; Steiner et al.2018; Scherler et al.2014) yields implausible fluxes when the atmosphere strongly deviates from near-neutral stability conditions, |Rib|>0.2. Equation (A2) is invalid for Rib>+0.2 and tends to diverge for highly unstable atmospheres at low wind speeds where Rib-. Filtering out these fluxes as proposed by Steiner et al. (2018) potentially deletes a sizeable portion of the calculated fluxes at the wind-sheltered Murtèl cirque.

A2 Modified Louis (1979) scheme

The Louis (1979) scheme is an analytical approximation of the iterative Monin–Obukhov scheme. The bulk exchange factor is expressed as Cbt=fhCHn (notation from Essery and Etchevers2004), where CHn=k2ln{zu/z0}-2 (the neutral exchange coefficient) and

(A4) f h = ( 1 + 10 Ri b ) - 1 , for Ri b 0 (stable case) 1 - 10 Ri b ( 1 + 10 C Hn - Ri b / f z ) - 1 , for Ri b < 0 (unstable case)


(A5) f z = 1 4 z 0 z u .

However, one shortcoming of the Louis (1979) scheme critical on rough surfaces is the assumption of equal momentum roughness length (z0m) and sensible heat roughness length (z0h), i.e. z0m=z0h=z0. Hence, we use the Louis (1979) scheme modified by Song (1998) (denoted by “cL”). The implementation is described in Rigon et al. (2006) and Endrizzi et al. (2014).

A3 Iterative Monin–Obukhov scheme

In the Monin–Obukhov scheme (denoted by “cMO”), the atmospheric stability is characterised by the Obukhov length L defined by

(A6) L = u 3 k g T a Q H ρ a C p ,

where u [m s−1] is the friction velocity. The bulk exchange factor for heat Cbt=CMO is given by (Conway and Cullen2013; Fitzpatrick et al.2017; Steiner et al.2018)

(A7) C MO = k 2 ln z u z 0 m - ψ m z u L + ψ m z 0 m L ln z T z 0 h - ψ h z T L + ψ h z 0 h L ,

where ψm and ψh are the vertically integrated stability functions for momentum and heat, respectively. The latent flux QLE is calculated analogously to scalar roughness length for water vapour z0q and the integrated stability function for water vapour ψq, both assumed equal to the corresponding quantity and function for heat. Since L depends itself on QH, CMO is calculated iteratively.

Different formulations are available for the integrated stability function ψ as used by, for example, Fitzpatrick et al. (2017), Steiner et al. (2018), and Sauter et al. (2020). We use the following momentum ψm and heat ψh stability functions (again assuming ψh=ψq) (Fitzpatrick et al.2017):

(A8) ψ m ( ζ ) = - a ζ + b ζ - c d exp - d ζ + b c d , ζ > 1 (stable (Beljaars and Holtslag, 1991)) ln 1 + χ 2 2 1 + χ 2 2 - 2 arctan χ + π 2 , ζ < 0 (unstable (Dyer, 1974))

and the heat stability function ψh

(A9) ψ h ( ζ ) = - 1 + 2 a 3 ζ 1.5 + b ζ - c d exp - d ζ + b c d - 1 , ζ > 1 (stable (Beljaars and Holtslag, 1991)) 2 ln 1 + χ 2 2 , ζ < 0 (unstable (Dyer, 1974))

with the dimensionless Obukhov parameter ζ:=zu/L (for ψm) or zT/L (for ψh), χ=(1-16ζ)1/4, a=1, b=2/3, c=5, and d=0.35.

Appendix B: Ground thermal and hygric regimes

The ground thermal and hygric regimes are shown in Figs. B1 and B2. These are the meteorological input variables for the calculation of the turbulent fluxes.

Figure B1Thermal regime. (a) Radiometric ground surface temperature Ts (“rGST”) and air temperatures in the atmosphere Ta and in the instrumented cavity (cavity roof Tsa at 0.67 m and mid-cavity level at 1.95 m depth). (b) Vertical temperature gradients zT between surface and 2 m air temperature (drives turbulent sensible fluxes; Eq. 12), across the snow cover (TaTsa), and within the cavity (indicates the stability of the in-cavity air column). The summer 2022 dry spells are referred to in the text.


Figure B2Hygric regime. (a) Relative and specific humidity of the saturated snow surface qss, in the atmosphere qa, and in the instrumented cavity (qsa cavity roof at 0.67 m and mid-cavity level at 1.95 m depth). (b) Vertical specific humidity gradients (zq) between the atmospheric air and the surface (either qss or qsa according to Eq. 11), across the snow cover (qaqsa), and within the cavity (indicates the moisture transport direction). The moisture transport within the cavity is generally downwards (zq>0), and condensation rather than evaporation occurs in the deep AL. Severe dry spells that last long enough to exhaust the near-surface moisture storage ( 5 d) strongly impact the ground moisture regime, reverse the moisture gradients, and lead to upwards moisture transport.


Figure B3(a) Vertical temperature and (b) specific humidity profile (thaw-season average). Temperature and humidity are highest in the near-surface AL. The near-surface AL in contact with the atmosphere responds to dry spells and dries out within a few days after the last precipitation event. The deep AL remains close to saturation.


Appendix C: Estimate of the thermal relaxation time

An estimate for the thermal relaxation time Δtr of a packed bed is

(C1) Δ t r ρ r c r h eff V A s = ρ r c r h eff ( d p / 2 ) 3 ( 1 - ϕ al ) ,

where heff is the effective heat transfer coefficient calculated via 1/heff=1/h+dp/(10kr) (Esence et al.2017). This correction accounts for the additional resistance arising from the temperature gradients within large blocks. “Large” means blocks whose Biot number exceeds 0.1 (Esence et al.2017), with the Biot number defined by

(C2) Bi := h k r L .

The characteristic length L is given by the ratio of the block volume to surface area: L:=V/As=dp/6 for spheres. A minimum value for the convective heat transfer coefficient h under forced convection at low wind speed (u<1m s−1) is  10 Wm-2K-1. The value is derived from inverting Eq. (12) with h:=QH/(Ta-Ts) (nighttime values) and agrees with an estimate in Bozhinskiy et al. (1986). With a thermal conductivity of the rock of kr=2.5Wm-1K-1, the critical diameter is 0.15 m. With dp=60 cm (typical dimension of the blocks enclosing the instrumented cavity), Eq. (C1) yields an estimate of Δtr=12 h. Thermal adjustment (within >95 %) is attained after tr or ∼1.5 d for large blocks with 60 cm edge length. Smaller blocks reach thermal equilibrium much faster, e.g. 20 cm blocks within 10 h.

Appendix D: 1997–2000 wind speed profiles

Here, we justify our “katabatic correction” of wind speed measurements with wind profile data collected by Stocker-Mittaz et al. (2002) and Stocker-Mittaz (2002). A 10 m tower installed on the Murtèl rock glacier at the PERMOS monitoring site measured the wind velocity, air temperature, and relative humidity at 1.5, 2.0, 6.5, and 9.1 m a.g.l. between January 1997 and March 2000.

Average vertical wind speed profiles (Fig. D1a) show the winter-time low-level jet between December and March with maximum wind speeds close to the snow-covered surface and upwards decreasing wind speed. In summer, wind speeds increase with height and approximately show the log wind profile (Stocker-Mittaz2002). The wind speed gradient is predominantly negative for more than  75 cm of snow measured at the PERMOS station (located on a plateau) (Fig. D1b), corresponding to  50 cm at the PERMA-XT station (located on a more wind-swept ridge). This finding justifies the compensation for snow height (Eq. 23) and possibly explains why the “katabatic correction” has less effect on the SEB of the snow-poor winter 2021–2022 and a negative effect on the November 2020 SEB (Fig. 10).

Figure D1(a) Wind speed profiles for winter (December–March average) and summer (July–September) 1997–2000 (Stocker-Mittaz et al.2002). Winter storms where wind speed at 9.1 m exceeded 2.5 m s−1 are filtered out (threshold from Stocker-Mittaz2002). (b) The average wind speed gradient between the 2.0 and 6.5 m levels switches from a (dominantly) positive log wind profile to a negative low-level jet profile at a (PERMOS) snow height of  60 cm (R2=0.271). Own figure based on data from Stocker-Mittaz et al. (2002).


Appendix E: Nomenclature

Parameters and constants used in this study are tabulated in Table E1.

Table E1Nomenclature: measurement variables, site-specific calibration parameters, dimensionless numbers, and constants.

Download Print Version | Download XLSX

Data availability

The PERMOS data can be obtained from the PERMOS network (, Hoelzle et al.2021, 2022) and the PERMA-XT measurement data from (last access: 5 April 2024; DOI:, Amschwand et al.2023).

Author contributions

DA performed the fieldwork, model development, and analyses for the study and wrote the manuscript. MS, MH, and BK supervised the study, provided financial and field support, and contributed to the manuscript preparation. AH and CK provided logistical support and editorial suggestions on the manuscript. HG designed the novel sensor array, regularly checked data quality, contributed to the analyses, and provided editorial suggestions on the manuscript.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.


This work is a collaboration between the University of Fribourg and GEOTEST. The authors wish to thank Walter Jäger (Waljag GmbH, Malans) and Thomas Sarbach (Sarbach Mechanik, St. Niklaus) for the technical support and the Corvatsch cable car company for logistical support. This publication is dedicated to Martin Scherler who laid the conceptual foundation.

Financial support

This research has been supported by Innosuisse – Schweizerische Agentur für Innovationsförderung (grant no. 36242.1 IP-EE).

Review statement

This paper was edited by Emily Collier and reviewed by two anonymous referees.


Amschwand, D., Scherler, M., Hoelzle, M., Krummenacher, B., Tschan, S., Aschwanden, L., and Gubler, H.: Murtèl rock glacier PERMA-XT data set (2020–2023), University of Fribourg/GEOTEST Zollikofen [data set],, 2023. a

Arenson, L. U., Hauck, C., Hilbich, C., Seward, L., Yamamoto, Y., and Springman, S. M.: Sub-surface heterogeneities in the Murtèl-Corvatsch rock glacier, Switzerland, in: Proceedings of the joint 63rd Canadian Geotechnical Conference and the 6th Canadian Permafrost Conference (GEO2010), Calgary (Alberta), Canada, 12–15 September 2010, ISBN 978-0-920505-46-5 (last access: 4 April 2024), 2010. a

Arya, S. P.: Introduction to Micrometeorology, no. 79 in International Geophysics Series, 2nd edn., Academic Press, ISBN 9780120593545, 2001. a

Balch, E. S.: Glacières or Freezing Caverns, Allen, Lane & Scott, Philadelphia, (last access: 4 April 2024), 1900. a

Beljaars, A. C. M. and Holtslag, A. A. M.: Flux Parameterization over Land Surfaces for Atmospheric Models, J. Appl. Meteorol., 30, 327–341,<0327:FPOLSF>2.0.CO;2, 1991. 

Bernhard, L., Sutter, F., Haeberli, W., and Keller, F.: Processes of snow/permafrost-interactions at a high mountain site, Murtèl/Corvatsch, eastern Swiss Alps, in: Proceedings of the 7th International Conference on Permafrost, Yellowknife, Northwest Territories, Canada, 23–27 June 1998, Université Laval (Québec), Canada, 35–41, (last access: 4 April 2024), 1998. a, b, c, d

Boike, J., Roth, K., and Ippisch, O.: Seasonal snow cover on frozen ground: Energy balance calculations of a permafrost site near Ny-Ålesund, Spitsbergen, J. Geophys. Res.-Atmos., 108, ALT 4-1–ALT 4-11,, 2003. a

Bowen, I. S.: The ratio of heat losses by conduction and by evaporation from any water surface, Phys. Rev., 27, 779–787,, 1926. a, b

Bozhinskiy, A. N., Krass, M. S., and Popovnin, V. V.: Role of debris cover in the thermal physics of glaciers, J. Glaciol., 32, 255–266,, 1986. a

Brighenti, S., Hotaling, S., Finn, D. S., Fountain, A. G., Hayashi, M., Herbst, D., Saros, J. E., Tronstad, L. M., and Millar, C. I.: Rock glaciers and related cold rocky landforms: Overlooked climate refugia for mountain biodiversity, Glob. Change Biol., 27, 1504–1517,, 2021. a

Brock, B. W., Willis, I. C., and Sharp, M. J.: Measurement and parameterization of aerodynamic roughness length variations at Haut Glacier d’Arolla, Switzerland, J. Glaciol., 52, 281–297,, 2006. a, b

Brock, B. W., Mihalcea, C., Kirkbride, M. P., Diolaiuti, G., Cutler, M. E. J., and Smiraglia, C.: Meteorology and surface energy fluxes in the 2005–2007 ablation seasons at the Miage debris-covered glacier, Mont Blanc Massif, Italian Alps, J. Geophys. Res.-Atmos., 115, D09106,, 2010. a, b, c

Brutsaert, W.: Evaporation into the Atmosphere, Environmental Fluid Mechanics, vol. 1, Springer Dordrecht, Netherlands,, 1982. a

Businger, J. A., Wyngaard, J. C., Izumi, Y., and Bradley, E. F.: Flux–profile relationships in the atmospheric surface layer, J. Atmos. Sci., 28, 181–189,<0181:FPRITA>2.0.CO;2, 1971. a

Collier, E., Nicholson, L. I., Brock, B. W., Maussion, F., Essery, R., and Bush, A. B. G.: Representing moisture fluxes and phase changes in glacier debris cover using a reservoir approach, The Cryosphere, 8, 1429–1444,, 2014. a, b

Conway, J. and Cullen, N.: Constraining turbulent heat flux parameterization over a temperate maritime glacier in New Zealand, Ann. Glaciol., 54, 41–51,, 2013. a, b

CSI CSAT3B manual: Manual for the Campbell Scientific CSAT3B three-dimensional sonic anemometer (v09/2020), Tech. Rep., Logan: Campbell Scientific Inc., (last access: 4 April 2024), 2020. a

CSI SR50A manual: Manual for the Campbell Scientific SR50A-series sonic ranging sensors (v02/2020), Tech. Rep., Logan: Campbell Scientific Inc., (last access: 4 April 2024), 2020. a

Delaloye, R. and Lambiel, C.: Evidence of winter ascending air circulation throughout talus slopes and rock glaciers situated in the lower belt of alpine discontinuous permafrost (Swiss Alps), Norsk Geogr. Tidsskr., 59, 194–203,, 2005. a, b

Delaloye, R., Reynard, E., Lambiel, C., Marescot, L., and Monnet, R.: Thermal anomaly in a cold scree slope (Creux du Van, Switzerland), in: Proceedings of the 8th International Conference on Permafrost, 21–25 July 2003, Zurich, Switzerland, edited by: Phillips, M., Springman, S. M., and Arenson, L. U., 175–180, Swets & Zeitlinger, Lisse, Zürich, ISBN 90 5809 582 7, (last access: 4 April 2024), 2003. a

Denby, B. and Greuell, W.: The use of bulk and profile methods for determining surface heat fluxes in the presence of glacier winds, J. Glaciol., 46, 445–452,, 2000. a

Dyer, A. J.: A review of flux–profile relationships, Bound.-Lay. Meteorol., 7, 363–372,, 1974. 

Endrizzi, S., Gruber, S., Dall'Amico, M., and Rigon, R.: GEOtop 2.0: simulating the combined energy and water balance at and below the land surface accounting for soil freezing, snow cover and terrain effects, Geosci. Model Dev., 7, 2831–2857,, 2014. a, b, c

Esence, T., Bruch, A., Molina, S., Stutz, B., and Fourmigué, J.-F.: A review on experience feedback and numerical modeling of packed-bed thermal energy storage systems, Sol. Energy, 153, 628–654,, 2017. a, b

Essery, R. and Etchevers, P.: Parameter sensitivity in simulations of snowmelt, J. Geophys. Res.-Atmos., 109, D20111,, 2004. a, b, c

Essery, R., Granger, R., and Pomeroy, J.: Boundary-layer growth and advection of heat over snow and soil patches: modelling and parameterization, Hydrol. Process., 20, 953–967,, 2006. a

Favier, V., Wagnon, P., Chazarin, J.-P., Maisincho, L., and Coudrain, A.: One-year measurements of surface heat budget on the ablation zone of Antizana Glacier 15, Ecuadorian Andes, J. Geophys. Res.-Atmos., 109, D18105,, 2004. a

Fitzpatrick, N., Radić, V., and Menounos, B.: Surface Energy Balance Closure and Turbulent Flux Parameterization on a Mid-Latitude Mountain Glacier, Purcell Mountains, Canada, Front. Earth Sci., 5, 67,, 2017. a, b, c, d, e, f, g

Fitzpatrick, N., Radić, V., and Menounos, B.: A multi-season investigation of glacier surface roughness lengths through in situ and remote observation, The Cryosphere, 13, 1051–1071,, 2019. a

Foken, T.: The energy balance closure problem: an overview, Ecol. Appl., 18, 1351–1367,, 2008. a, b, c

Foken, T.: Micrometeorology, Springer Berlin, Heidelberg, 2 edn., ISBN 978-3-642-25440-6,, 2017. a, b

Foken, T., Leuning, R., Oncley, S. R., Mauder, M., and Aubinet, M.: Corrections and Data Quality Control, in: Eddy Covariance. A Practical Guide to Measurement and Data Analysis, edited by: Aubinet, M., Vesala, T., and Papale, D., Springer Netherlands, 85–131,, 2012. a

Frauenfelder, R. and Kääb, A.: Towards a palaeoclimatic model of rock-glacier formation in the Swiss Alps, Ann. Glaciol., 31, 281–286,, 2000. a

Froidurot, S., Zin, I., Hingray, B., and Gautheron, A.: Sensitivity of precipitation phase over the Swiss Alps to different meteorological variables, J. Hydrometeorol., 15, 685–696,, 2014. a

Giese, A., Boone, A., Wagnon, P., and Hawley, R.: Incorporating moisture content in surface energy balance modeling of a debris-covered glacier, The Cryosphere, 14, 1555–1577,, 2020. a

Gorbunov, A. P., Marchenko, S. S., and Seversky, E. V.: The thermal environment of blocky materials in the mountains of Central Asia, Permafrost Periglac., 15, 95–98,, 2004. a

Grisogono, B., Kraljević, L., and Jeričević, A.: The low-level katabatic jet height versus Monin–Obukhov height, Q. J. Roy. Meteor. Soc., 133, 2133–2136,, 2007. a, b

Guodong, C., Yuanming, L., Zhizhong, S., and Fan, J.: The “thermal semi-conductor” effect of crushed rocks, Permafrost Periglac., 18, 151–160,, 2007. a

Haberkorn, A., Hoelzle, M., Phillips, M., and Kenner, R.: Snow as a driving factor of rock surface temperatures in steep rough rock walls, Cold Reg. Sci. Technol., 118, 64–75,, 2015. a

Haeberli, W.: Die Basis-Temperatur der winterlichen Schneedecke als moglicher Indikator fur die Verbreitung von Permafrost in den Alpen, Zeitschrift für Gletscherkunde und Glazialgeologie, 9, 221–227, 1973. a, b

Haeberli, W. and Patzelt, G.: Permafrostkartierung im Gebiet der Hochebenkar-Blockgletscher, Obergurgl, Ötztaler Alpen, Zeitschrift für Gletscherkunde und Glazialgeologie, 18, 127–150, 1982. a, b

Haeberli, W., Hallet, B., Arenson, L., Elconin, R., Humlum, O., Kääb, A., Kaufmann, V., Ladanyi, B., Matsuoka, N., Springman, S., and Mühll, D. V.: Permafrost creep and rock glacier dynamics, Permafrost Periglac., 17, 189–214,, 2006. a

Hanson, S. and Hoelzle, M.: The thermal regime of the active layer at the Murtèl rock glacier based on data from 2002, Permafrost Periglac., 15, 273–282,, 2004. a, b, c

Hanson, S. and Hoelzle, M.: Installation of a shallow borehole network and monitoring of the ground thermal regime of a high alpine discontinuous permafrost environment, Eastern Swiss Alps, Norsk Geogr. Tidsskr., 59, 84–93,, 2005. a

Harris, S. A. and Pedersen, D. E.: Thermal regimes beneath coarse blocky materials, Permafrost Periglac., 9, 107–120,<107::AID-PPP277>3.0.CO;2-G, 1998. a

Haugen, D. A., Kaimal, J. C., and Bradley, E. F.: An experimental study of Reynolds stress and heat flux in the atmospheric surface layer, Q. J. Roy. Meteor. Soc., 97, 168–180,, 1971. a

Herz, T.: Das Mikroklima grobblockiger Schutthalden der alpinen Periglazialstufe und seine Auswirkungen auf Energieaustauschprozesse zwischen Atmosphäre und Lithosphäre [The microclimate of coarse debris covers in the periglacial belt of high mountains and its effects on the energy exchange between atmosphere and lithosphere], PhD thesis, Justus-Liebig-Universität Gießen,, 2006. a, b

Herz, T., King, L., and Gubler, H.: Microclimate within coarse debris of talus slopes in the alpine periglacial belt and its effect on permafrost, in: Proceedings of the 8th International Conference on Permafrost, Zurich, Switzerland, 21–25 July 2003, edited by: Phillips, M., Springman, S. M., and Arenson, L. U., Swets & Zeitlinger, Lisse, Zürich, 383–387, 2003. a

Hoelzle, M.: Mapping and modelling of mountain permafrost distribution in the Alps, Norsk Geogr. Tidsskr., 50, 11–15,, 1996. a

Hoelzle, M. and Gruber, S.: Borehole and ground surface temperatures and their relationship to meteorological conditions in the Swiss Alps, in: Proceedings of the 9th International Conference on Permafrost, Fairbanks, Alaska, 29 June–3 July 2008, edited by: Kane, D. and Hinkel, K., 723–728, Institute of Northern Engineering, University of Alaska, Fairbanks, Alaska, ISBN 978-0-9800179-2-2,, 2008. a

Hoelzle, M. and Haeberli, W.: Simulating the effects of mean annual air-temperature changes on permafrost distribution and glacier size: an example from the Upper Engadin, Swiss Alps, Ann. Glaciol., 21, 399–405,, 1995. a

Hoelzle, M., Wegmann, M., and Krummenacher, B.: Miniature temperature dataloggers for mapping and monitoring of permafrost in high mountain areas: first experience from the Swiss Alps, Permafrost Periglac., 10, 113–124,<113::AID-PPP317>3.0.CO;2-A, 1999. a

Hoelzle, M., Mittaz, C., Etzelmüller, B., and Haeberli, W.: Surface energy fluxes and distribution models of permafrost in European mountain areas: an overview of current developments, Permafrost Periglac., 12, 53–68,, 2001. a, b

Hoelzle, M., Mühll, D. V., and Haeberli, W.: Thirty years of permafrost research in the Corvatsch–Furtschellas area, Eastern Swiss Alps: A review, Norsk Geogr. Tidsskr., 56, 137–145,, 2002. a

Hoelzle, M., Haeberli, W., and Stocker-Mittaz, C.: Miniature ground temperature data logger measurements 2000–2002 in the Murtèl-Corvatsch area, Eastern Swiss Alps, in: Proceedings of the 8th International Conference on Permafrost, Zurich, Switzerland, 21–25 July 2003, edited by: Phillips, M., Springman, S. M., and Arenson, L. U., Swets & Zeitlinger, Lisse, Zürich, 419–424, ISBN 90 5809 582 7, (last access: 4 April 2024), 2003. a, b

Hoelzle, M., Hauck, C., Mathys, T., Noetzli, J., Pellet, C., and Scherler, M.: Energy balance measurements at three PERMOS sites: Corvatsch, Schilthorn, Stockhorn, PERMOS [data set],, 2021. a

Hoelzle, M., Hauck, C., Mathys, T., Noetzli, J., Pellet, C., and Scherler, M.: Long-term energy balance measurements at three different mountain permafrost sites in the Swiss Alps, Earth Syst. Sci. Data, 14, 1531–1547,, 2022. a, b, c, d, e, f, g, h, i, j, k

Hukseflux WS01 manual: Manual for the Hukseflux WS01 sensor for ultra low wind speeds and boundary layer conductance (WS01 manual v0608), Tech. Rep., Delft: Hukseflux Thermal Sensors, (last access: 4 April 2024), 2006. a

Humlum, O.: Active layer thermal regime at three rock glaciers in Greenland, Permafrost Periglac., 8, 383–408,<383::AID-PPP265>3.0.CO;2-V, 1997. a

Isaksen, K., Heggem, E., Bakkehøi, S., Ødegård, R., Eiken, T., Etzelmüller, B., and Sollid, J.: Mountain permafrost and energy balance on Juvvasshøe, southern Norway, in: Proceedings of the 8th International Conference on Permafrost, Zurich, Switzerland, 21–25 July 2003, edited by: Phillips, M., Springman, S., and Arenson, L., Swets & Zeitlinger, Lisse, Zürich, 467–472, ISBN 90 5809 582 7, (last access: 4 April 2024), 2003. a, b

Kääb, A., Gudmundsson, G. H., and Hoelzle, M.: Surface deformation of creeping mountain permafrost. Photogrammetric investigations on Murtèl rock glacier, Swiss Alps, in: Proceedings of the 7th International Conference on Permafrost, Yellowknife, Northwest Territories, Canada, 23–27 June 1998, Université Laval (Qébéc), Canada, 531–537, (last access: 4 April 2024), 1998. a

Keller, F.: Interaktionen zwischen Schnee und Permafrost, Mitteilungen der VAW/ETH Zürich 127, VAW/ETH Zürich, (last access: 4 April 2024), 1994. a

Keller, F. and Gubler, H.: Interaction between snow cover and high mountain permafrost: Murtèl/Corvatsch, Swiss Alps, in: Proceedings of the 6th International Conference on Permafrost, Beijing, China, July 5–9 1993, edited by: Guodong, C., vol. 1, Lanzhou Institute of Glaciology and Geocryology, Chinese Academy of Sciences & Chinese Society of Glaciology and Geocryology, South China University of Technology Press, 332–337, (last access: 4 April 2024), 1993. a, b, c, d

Kellerer‐Pirklbauer, A., Pauritsch, M., and Winkler, G.: Widespread occurrence of ephemeral funnel hoarfrost and related air ventilation in coarse‐grained sediments of a relict rock glacier in the Seckauer Tauern Range, Austria, Geogr. Ann. A, 97, 453–471,, 2015. a

Kipp & Zonen CGR3 manual: Manual for the Kipp & Zonen CGR 4 pyrgeometer (v1401), Tech. Rep., Delft: Kipp & Zonen, (last access: 4 April 2024), 2014. a

Kraaijenbrink, P. D. A., Shea, J. M., Litt, M., Steiner, J. F., Treichler, D., Koch, I., and Immerzeel, W. W.: Mapping surface temperatures on a debris-covered glacier with an unmanned aerial vehicle, Front. Earth Sci., 6, 1–19,, 2018. a

Liebethal, C. and Foken, T.: Evaluation of six parameterization approaches for the ground heat flux, Theor. Appl. Climatol., 88, 43–56,, 2007. a

Liebethal, C., Huwe, B., and Foken, T.: Sensitivity analysis for two ground heat flux calculation approaches, Agr. Forest Meteorol., 132, 253–262,, 2005. a, b

Liu, H. and Foken, T.: A modified Bowen ratio method to determine sensible and latent heat fluxes, Meteorol. Z., 10, 71–80,, 2001. a, b

Louis, J.-F.: A parametric model of vertical eddy fluxes in the atmosphere, Bound.-Lay. Meteorol., 17, 187–202,, 1979. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r

Marmy, A., Rajczak, J., Delaloye, R., Hilbich, C., Hoelzle, M., Kotlarski, S., Lambiel, C., Noetzli, J., Phillips, M., Salzmann, N., Staub, B., and Hauck, C.: Semi-automated calibration method for modelling of mountain permafrost evolution in Switzerland, The Cryosphere, 10, 2693–2719,, 2016. a

Massman, W.: A simple method for estimating frequency response corrections for eddy covariance systems, Agr. Forest Meteorol., 104, 185–198,, 2000. a

Miles, E. S., Steiner, J. F., and Brun, F.: Highly variable aerodynamic roughness length (z0) for a hummocky debris-covered glacier, J. Geophys. Res.-Atmos., 122, 8447–8466,, 2017. a, b

Millar, C. I., Westfall, R. D., and Delany, D. L.: Thermal regimes and snowpack relations of periglacial talus slopes, Sierra Nevada, California, U.S.A., Arct. Antarct. Alp. Res., 46, 483–504,, 2014. a, b

Mittaz, C., Hoelzle, M., and Haeberli, W.: First results and interpretation of energy-flux measurements over Alpine permafrost, Ann. Glaciol., 31, 275–280,, 2000. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o

Moncrieff, J., Massheder, J., de Bruin, H., Elbers, J., Friborg, T., Heusinkveld, B., Kabat, P., Scott, S., Soegaard, H., and Verhoef, A.: A system to measure surface fluxes of momentum, sensible heat, water vapour and carbon dioxide, J. Hydrol., 188-189, 589–611,, hAPEX-Sahel, 1997. a

Morard, S., Delaloye, R., and Dorthe, J.: Seasonal thermal regime of a mid-latitude ventilated debris accumulation, in: Proceedings of the 9th International Conference on Permafrost, Fairbanks, Alaska, 29 June–3 July 2008, edited by: Kane, D. L. and Hinkel, K. M., vol. 2, 1233–1238, ISBN 978-0-9800179-2-2, 2008. a

Mott, R., Gromke, C., Grünewald, T., and Lehning, M.: Relative importance of advective heat transport and boundary layer decoupling in the melt dynamics of a patchy snow cover, Adv. Water Resour., 55, 88–97,, 2013. a

Mott, R., Daniels, M., and Lehning, M.: Atmospheric flow development and associated changes in turbulent sensible heat flux over a patchy mountain snow cover, J. Hydrometeorol., 16, 1315–1340,, 2015. a

Müller, J., Gärtner-Roer, I., Kenner, R., Thee, P., and Morche, D.: Sediment storage and transfer on a periglacial mountain slope (Corvatsch, Switzerland), Geomorphology, 218, 35–44,, 2014. a

Nicholson, L. and Stiperski, I.: Comparison of turbulent structures and energy fluxes over exposed and debris-covered glacier ice, J. Glaciol., 66, 543–555,, 2020. a, b

Nield, D. A. and Bejan, A.: Convection in Porous Media, 5 edn., Springer, New York,, 2017. a

Noetzli, J. and Pellet, C., eds.: PERMOS 2023. Swiss Permafrost Bulletin 2022 (Annual report No. 4 on permafrost observation in the Swiss Alps), Cryospheric Commission of the Swiss Academy of Sciences,, 2023. a

Noetzli, J., Pellet, C., and Staub, B. (Eds.): Permafrost in Switzerland 2014/2015 to 2017/2018, Glaciological Report Permafrost No. 16–19 (PERMOS Report 2019), Cryospheric Commission of the Swiss Academy of Sciences, Fribourg,, 2019. a

Ochsner, T. E., Sauer, T. J., and Horton, R.: Soil heat storage measurements in energy balance studies, Agron. J., 99, 311–319,, 2007. a

Oerlemans, J. and Grisogono, B.: Glacier winds and parameterisation of the related surface heat fluxes, Tellus A, 54, 440,, 2002. a, b, c

Ohmura, A.: Objective criteria for rejecting data for Bowen ratio flux calculations, J. Appl. Meteorol. Clim., 21, 595–598,<0595:OCFRDF>2.0.CO;2, 1982. a, b

Oke, T. R.: Boundary Layer Climates, 2 edn., Routledge, London,, 1987. a, b, c, d, e

Otto, J., Keuschnig, M., Götz, J., Marbach, M., and Schrott, L.: Detection of mountain permafrost by combining high resolution surface and subsurface information – an example from the Glatzbach catchment, Austrian Alps, Geogr. Ann. A, 94, 43–57,, 2012. a

Panz, M.: Analyse von Austauschprozessen in der Auftauschicht des Blockgletschers Murtèl, Corvatsch, Oberengadin, Master's thesis, Geographisches Institut der Ruhr-Universität Bochum, 2008. a, b

Pérez, F. L.: Conservation of soil moisture by different stone covers on alpine talus slopes (Lassen, California), CATENA, 33, 155–177,, 1998. a

Petersen, E., Hock, R., Fochesatto, G. J., and Anderson, L. S.: The significance of convection in supraglacial debris revealed through novel analysis of thermistor profiles, J. Geophys. Res.-Earth, 127, e2021JF006520,, 2022. a

Popescu, R., Vespremeanu-Stroe, A., Onaca, A., Vasile, M., Cruceru, N., and Pop, O.: Low-altitude permafrost research in an overcooled talus slope–rock glacier system in the Romanian Carpathians (Detunata Goală, Apuseni Mountains), Geomorphology, 295, 840–854,, 2017. a, b

Radić, V., Menounos, B., Shea, J., Fitzpatrick, N., Tessema, M. A., and Déry, S. J.: Evaluation of different methods to model near-surface turbulent fluxes for a mountain glacier in the Cariboo Mountains, BC, Canada, The Cryosphere, 11, 2897–2918,, 2017. a, b, c, d, e, f, g

Reid, T. D. and Brock, B. W.: An energy-balance model for debris-covered glaciers including heat conduction through the debris layer, J. Glaciol., 56, 903–916,, 2010. a, b, c, d, e

Reid, T. D., Carenzo, M., Pellicciotti, F., and Brock, B. W.: Including debris cover effects in a distributed model of glacier ablation, J. Geophys. Res.-Atmos., 117, D18105,, 2012. a

Rigon, R., Bertoldi, G., and Over, T. M.: GEOtop: a distributed hydrological model with coupled water and energy budgets, J. Hydrometeorol., 7, 371–388,, 2006. a, b

Rist, A.: Hydrothermal processes within the active layer above alpine permafrost in steep scree slopes and their influence on slope stability, PhD thesis, University of Zurich, Zürich,, 2007. a

Rist, A. and Phillips, M.: First results of investigations on hydrothermal processes within the active layer above alpine permafrost in steep terrain, Norsk Geogr. Tidsskr., 59, 177–183,, 2005. a

Rounce, D. R., Quincey, D. J., and McKinney, D. C.: Debris-covered glacier energy balance model for Imja–Lhotse Shar Glacier in the Everest region of Nepal, The Cryosphere, 9, 2295–2310,, 2015. a

Růžička, V., Zacharda, M., Němcová, L., Šmilauer, P., and Nekola, J. C.: Periglacial microclimate in low-altitude scree slopes supports relict biodiversity, J. Nat. Hist., 46, 2145–2157,, 2012. a

Sakai, A., Fujita, K., and Kubota, J.: Evaporation and percolation effect on melting at debris-covered Lirung Glacier, Nepal Himalayas, 1996, Bulletin of Glaciological Research, 21, 9–15, 2004. a

Sauter, T., Arndt, A., and Schneider, C.: COSIPY v1.3 – an open-source coupled snowpack and ice surface energy and mass balance model, Geosci. Model Dev., 13, 5645–5662,, 2020. a

Sawada, Y., Ishikawa, M., and Ono, Y.: Thermal regime of sporadic permafrost in a block slope on Mt. Nishi-Nupukaushinupuri, Hokkaido Island, Northern Japan, Geomorphology, 52, 121–130,, 2003. a

Scherler, M., Hauck, C., Hoelzle, M., and Salzmann, N.: Modeled sensitivity of two alpine permafrost sites to RCM-based climate scenarios, J. Geophys. Res.-Earth, 118, 780–794,, 2013. a

Scherler, M., Schneider, S., Hoelzle, M., and Hauck, C.: A two-sided approach to estimate heat transfer processes within the active layer of the Murtèl–Corvatsch rock glacier, Earth Surf. Dynam., 2, 141–154,, 2014. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r

Schneider, S.: The heterogeneity of mountain permafrost – A field-based analysis of different periglacial materials, PhD thesis, University of Fribourg (Switzerland), 2014. a, b

Schneider, S., Hoelzle, M., and Hauck, C.: Influence of surface and subsurface heterogeneity on observed borehole temperatures at a mountain permafrost site in the Upper Engadine, Swiss Alps, The Cryosphere, 6, 517–531,, 2012. a, b

Schneider, S., Daengeli, S., Hauck, C., and Hoelzle, M.: A spatial and temporal analysis of different periglacial materials by using geoelectrical, seismic and borehole temperature data at Murtèl–Corvatsch, Upper Engadin, Swiss Alps, Geogr. Helv., 68, 265–280,, 2013. a

Schotanus, P., Nieuwstadt, F., and Bruin, H. D.: Temperature measurement with a sonic anemometer and its application to heat and moisture fluxes, Bound.-Lay. Meteorol., 26, 81–93,, 1983. a

Shimokawabe, A., Yamaura, Y., Akasaka, T., Sato, T., Shida, Y., Yamanaka, S., and Nakamura, F.: The distribution of cool spots as microrefugia in a mountainous area, PLoS ONE, 10, e0135732,, 2015. a

Sicart, J. E., Wagnon, P., and Ribstein, P.: Atmospheric controls of the heat balance of Zongo Glacier (16° S, Bolivia), J. Geophys. Res.-Atmos., 110, D12106,, 2005. a, b

Sicart, J. E., Litt, M., Helgason, W., Tahar, V. B., and Chaperon, T.: A study of the atmospheric surface layer and roughness lengths on the high-altitude tropical Zongo glacier, Bolivia, J. Geophys. Res.-Atmos., 119, 3793–3808,, 2014. a, b

Smeets, C. J. P. P. and van den Broeke, M. R.: The parameterisation of scalar transfer over rough ice, Bound.-Lay. Meteorol., 128, 339–355,, 2008. a, b

Smeets, C. J. P. P., Duynkerke, P. G., and Vugts, H. F.: Observed wind profiles and turbulence fluxes over an ice surface with changing surface roughness, Bound.-Lay. Meteorol., 92, 99–121,, 1999. a

Smith, M. W.: Roughness in the Earth Sciences, Earth-Sci. Rev., 136, 202–225,, 2014. a

Song, Y.: An improvement of the Louis scheme for the surface layer in an atmospheric modelling system, Bound.-Lay. Meteorol., 88, 239–254,, 1998. a, b, c

Steiner, J. F., Litt, M., Stigter, E. E., Shea, J., Bierkens, M. F. P., and Immerzeel, W. W.: The importance of turbulent fluxes in the surface energy balance of a debris-covered glacier in the Himalayas, Front. Earth Sci., 6, 144,, 2018. a, b, c, d, e, f, g, h, i, j, k, l, m

Steiner, J. F., Kraaijenbrink, P. D. A., and Immerzeel, W. W.: Distributed melt on a debris-covered glacier: Field observations and melt modeling on the Lirung Glacier in the Himalaya, Front. Earth Sci., 9,, 2021. a, b

Stocker-Mittaz, C.: Permafrost Distribution Modeling Based on Energy Balance Data, PhD thesis, University of Zurich, 2002. a, b, c, d, e, f

Stocker-Mittaz, C., Hoelzle, M., and Haeberli, W.: Modelling alpine permafrost distribution based on energy-balance data: a first step, Permafrost Periglac., 13, 271–282,, 2002. a, b, c, d

Stull, R.: Practical Meteorology: An Algebra-based Survey of Atmospheric Science, University of British Columbia, version 1.02b, 940 pp., (last access: 2 April 2024), 2017. a

Stull, R. B.: Static stability–an update, B. Am. Meteorol. Soc., 72, 1521–1529,<1521:SSU>2.0.CO;2, 1991. a

Tanner, C. B. and Thurtell, G. W.: Anemoclinometer measurements of Reynolds stress and heat transport in the atmospheric surface layer, techreport, US Army Electronics Command, Atmospheric Sciences Laboratory TR ECOM 66-G22-F, Wisconsin Univ-Madison Dept of Soil Science, (last access: 2 April 2024), 1969. a

Vickers, D. and Mahrt, L.: Quality control and flux sampling problems for tower and aircraft data, J. Atmos. Ocean. Tech., 14, 512–526,<0512:QCAFSP>2.0.CO;2, 1997. a

Vonder Mühll, D. and Haeberli, W.: Thermal characteristics of the permafrost within an active rock glacier (Murtèl/Corvatsch, Grisons, Swiss Alps), J. Glaciol., 36, 151–158,, 1990. a

Wagner, T., Pauritsch, M., Mayaud, C., Kellerer-Pirklbauer, A., Thalheim, F., and Winkler, G.: Controlling factors of microclimate in blocky surface layers of two nearby relict rock glaciers (Niedere Tauern Range, Austria), Geogr. Ann. A, 101, 310–333,, 2019. a

Wakonigg, H.: Unterkühlte Schutthalden [Undercooled talus], Arbeiten aus dem Institut für Geographie der Karl-Franzens Universität Graz, 33, 209–223, 1996.  a, b

Wicky, J. and Hauck, C.: Numerical modelling of convective heat transport by air flow in permafrost talus slopes, The Cryosphere, 11, 1311–1325,, 2017. a

Wicky, J. and Hauck, C.: Air convection in the active layer of rock glaciers, Front. Earth Sci., 8, 335,, 2020. a

Winkler, M., Schellander, H., and Gruber, S.: Snow water equivalents exclusively from snow depths and their temporal changes: the Δsnow model, Hydrol. Earth Syst. Sci., 25, 1165–1187,, 2021. a, b, c

Wunder, J. and Möseler, B. M.: Kaltluftströme auf Basaltblockhalden und ihre Auswirkung auf Mikroklima und Vegetation, Flora, 191, 335–344,, 1996. a

Yao, J., Gu, L., Han, H., Wang, Y., and Liu, S.: The surface energy budget on the debris-covered Koxkar Glacier in China, Environ. Earth Sci., 72, 4503–4510,, 2014. a, b

Zhang, T.: Influence of the seasonal snow cover on the ground thermal regime: An overview, Rev. Geophys., 43, RG4002,, 2005. a

Short summary
Rock glaciers are coarse-debris permafrost landforms that are comparatively climate resilient. We estimate the surface energy balance of rock glacier Murtèl (Swiss Alps) based on a large surface and sub-surface sensor array. During the thaw seasons 2021 and 2022, 90 % of the net radiation was exported via turbulent heat fluxes and only 10 % was transmitted towards the ground ice table. However, early snowmelt and droughts make these permafrost landforms vulnerable to climate warming.