Surface energy balance of the Sygyktinsky Glacier, south Eastern Siberia, during the ablation period and its sensitivity to meteorological fluctuations

The physically based melt of the low elevation Eastern Siberian glaciers is poorly understood due to the lack of direct micrometeorological studies. We used an automatic meteorological station to record the meteorological and energy characteristics of the Sygyktinsky Glacier, south Eastern Siberia (56.8° N, 117.4° E, 2,560 m a.s.l.), during two ablation seasons and computed the surface energy balance (SEB) for 30-min intervals. The glacier ablation was both modeled and measured by stakes and a thermistor cable. The net radiation (Rnet) was the main contributor (71–75 W m−2, 89–95%) to the SEB (79 W m−2, 100%), followed by sensible (2–4 W m−2, 3–5%) and latent (2–3 W m−2, 2–4%) heat fluxes. The net shortwave radiation was the main positive component of Rnet, while the net longwave radiation was weak and either negative (− 15 W m−2 in 2019) or positive (4 W m−2 in 2020). The small proportion of turbulent fluxes in the SEB is explained by the low wind speed (1.2 m s−1). The glacier ablation was found to be more sensitive to changes in shortwave radiation and wind speed, suggesting the need to consider the atmospheric conditions of the ablation period (summer snowfalls, cloudiness, wind speed) when analyzing long-term trends in glacial changes.


Results
The meteorological conditions averaged for two ablation seasons (for the same period from July 7 to August 23) are shown in Table 1, and daily mean values of the measured parameters are shown in Fig. 2. of cases above the ice melting point (0 °C) and ranged from − 3.2 to 20.4 °C (on average 7.3 °C, standard deviation 3.8 °C). In 2020, the summer period was cooler: the temperature was positive in 93% of cases and ranged from − 4.5 to 21.0 °C (on average 6.3 °C, standard deviation 3.6 °C). From July to August, the air temperature decreased in both seasons. The surface temperature was approximately 0 °C (melting snow and ice). The temperature gradient between the glacier surface and the 2 m level was positive, on average 3.1 °C m −1 . This indicates the predominance of stable conditions in the boundary atmospheric layer (temperature inversion).
The diurnal cycle of temperature fluctuates in the range of about 4 °C and has a pronounced maximum in the afternoon (Fig. 3). The minimum temperature is observed during the night (between 3:00 and 4:30). The daily cycle of relative humidity varies in the range of 17% and has a minimum in the afternoon. The specific humidity has a small daily cycle (within 0.4 g kg −1 ) with a maximum in the afternoon (between 15:00 and 18:00).
The average relative air humidity at 2 m above the glacier surface (AWS1) was 76 ± 20% in 2019 and 85 ± 14% in 2020 (Table 1). In 2019, the daily humidity varied from 41 to 100%, and in 2020, from 51 to 100%. Fluctuations in humidity with a frequency of several days are due to a change in synoptic processes over the region (maxima correspond to cyclones, and minima to anticyclones). The daily specific humidity in 2019 ranged from 4.1 to 8.6 g kg −1 , and in 2020, from 3.7 to 9.7 g kg −1 (Fig. 2). In both seasons, there is a tendency for the specific humidity to decrease from July to August.
Wind speed and direction. The wind regime on the glacier (AWS1) is characterized by the predominance of weak winds (Table 1). In 2020, daily average wind speeds at 2 m above the glacier surface ranged from 0 to 3.3 m s −1 , with an average of 1.3 m s −1 (standard deviation 0.8 m s −1 ). The observed half-hour wind speed reached 8.0 m s −1 (Fig. S2), and the maximum wind speed reached 19.1 m s −1 . The average and maximum wind speed on the glacier (AWS1) were, respectively, 20% and 10% higher than on the moraine (AWS2), due to the increased openness of AWS1. The correlation between the wind speed on the glacier and moraine is statistically significant but not high (R 2 = 0.35, n = 2247, p < 0.001). The correlation between the maximum wind speeds is higher (R 2 = 0.52, n = 2247, p < 0.001).
The prevailing wind direction on the glacier (AWS1) is defined by synoptic processes and local topography (Fig. S2). Westerly wind directions (W + WSW + WNW, 248-293°) have the highest occurrence (30%) and higher speeds, and they are associated with large-scale atmospheric circulation. Over the glacier, at the level of 700 hPa, southwestern winds prevailed (25%) during the observation period of 2020. There is a pronounced secondary maximum of the northerly wind (occurrence 10%), which is weak and flowing down from the surrounding mountain slopes. The low-frequency periodicity in the maximum wind speed (5-8 days) is associated with the passage of atmospheric fronts (Fig. 2).
The wind has a weakly expressed daily cycle with increased values in the first half of the day and decreased values in the second (Fig. 3). The maximum speed values are observed in the morning (8:30) and in the afternoon (13:00), and the minimum in the evening (21:00).

Cloudiness.
The summer period at the Kodar is characterized by overcast conditions and rainy weather due to increased cyclonic activity 30 . In 2020 the cloudiness (data from the Chara weather station) was higher than in Table 1. Mean values of meteorological and radiation parameters on the glacier for the period 7 July to 23 August. Cloudiness data were obtained at the Chara weather station.  Fig. 2). In 2019, the cloudiness ranged from 3 to 100% (mean value is 68%, standard deviation 25%), and in 2020 from 18 to 100% (79% on average, standard deviation 23%). The highest cloud cover was associated with cyclones, and the lowest with anticyclones. The daily averaged cloud cover in Chara is statistically significantly correlated (2019: R 2 = 0.53, n = 48, p < 0.001; 2020: R 2 = 0.62, n = 48, p < 0.001) with the cloudiness over the glacier, calculated as the ratio between shortwave radiation incoming at AWS1 and at the top of the atmosphere S TOA 18 .
Energy fluxes. Shortwave radiation. The incoming summer shortwave radiation (S in ) at the glacier (AWS1) was low due to increased cloudiness related to cyclonic activity over the study area (Table 1). The atmospheric transmissivity is 0.34 on average, i.e., only one third of the solar radiation entering the top of the atmosphere reached the glacial surface. S in has a strong variability, both day-to-day and interannual (Fig. 4)  www.nature.com/scientificreports/ in solar radiation arriving at the top of the atmosphere. The lower amount of incoming solar radiation in 2020 compared to 2019 (21% less) is explained by more significant (17% more) cloud cover in 2020 (Table 1). S in has a distinct diurnal cycle (Fig. 3), with the daily range in 2019 being 1.2 times higher than in 2020. Reflected shortwave radiation (S ref , AWS1) is characterized by strong day-to-day variability and tends to decrease towards the end of the ablation season due to an increase in the glacial surface albedo (Fig. 4). At the same time, its interannual variability is relatively small (Table 1), since AWS1 is located near the long-term equilibrium line of the glacier, and the change of the snow surface from snow to ice occurs every year.
The net shortwave radiation ranged from 11 to 225 W m −2 in 2019 (on average 90 W m −2 ) and from 11 to 160 W m −2 in 2020 (on average 67 W m −2 ). The net shortwave radiation in August was higher than in July (in 2019 by 7 W m −2 , in 2020 by 20 W m −2 ) due to a decrease in the reflected radiation in August (in 2019, by 3.4 times; in 2020, by 1.5 times).
Albedo. From July to August 2019, albedo decreased linearly from 0.73 to 0.19 (Fig. 4). The decrease is associated with snow melting, dust deposition on the glacier surface from adjacent slopes, and the melting out of finegrained material. The snow cover completely melted on 4 August. The background albedo values in August were low (about 0.20); however, pronounced peaks in August (0.40 to 0.88) marked the summer snowfalls (e.g., 14-15 and 21-23 August). The duration of the periods with high albedo was no more than 2-3 days, indicating that the snow melted quickly. In general, the 2019 ablation season is characterized by two different albedo regimes, reflecting different physical conditions of the glacier surface: (i) a gradual linear-like decrease in albedo in July (snow) and (ii) low background values with pronounced peaks of short-term summer snowfalls in August (ice). The average albedo in the 2019 ablation season was 0.40 (0.57 in July, 0.29 in August). A stable (winter) snow cover on the glacier surface began to form on 6 September.
At the beginning of July 2020, due to the higher winter accumulation, the thickness of the snow cover on the glacier near AWS1 was 2.5 times (4.7 times in w.e.) more than at the beginning of July 2019. Furthermore, the snow cover completely melted here only at the end of August. Accordingly, the change in albedo in 2020 differed Longwave radiation. The sensors of longwave radiation were installed only on the moraine (AWS2) to prevent the risk of losing them due to snow avalanches failing down from surrounding slopes to glacier surface. Since during both ablation seasons the glacier surface was melting (about 0 °C) the emitted longwave radiation flux was taken constant (316 W m −2 ) 11 . Both weather stations (AWS1 and AWS2) are located close to each other (300 m in distance and 30 m in altitude) within the same boundary layer, therefore we assumed the incoming longwave radiation measured off the ice (L in , AWS2) was corresponded to that on the glacier (AWS1). Compared to shortwave radiation ( Fig. 4), L in had insignificant interannual and day-to-day variability and ranged from 236 to 344 W m −2 in 2019 (on average 301 W m −2 , standard deviation 28 W m −2 ) and from 222 to 361 W m −2 in 2020 (on average 319 W m −2 , standard deviation 31 W m −2 ). It changed slightly within the ablation season and was maximum in July 2020 (328 W m −2 ) and minimum in August 2019 (300 W m −2 ). However, the net longwave radiation in 2019 was slightly different from that in 2020. In 2019, it varied from − 80 to 29 W m −2 (average − 15 W m −2 ), while in 2020, it varied from − 94 to 46 W m −2 (on average 4 W m −2 ). That is, if in the 2019 observation period, the net longwave radiation was a sink of energy, then in 2020, on the contrary, the glacier surface received additional energy from the atmosphere. In 2019, the daily net longwave radiation was positive for 20 days (42% of the period), and in 2020, for 31 days (65%) due to increased cloudiness in 2020. The relationship between L in and cloudiness is confirmed by a good correlation between them (R 2 is 0.52 in 2019 and 0.56 in 2020, n = 48). www.nature.com/scientificreports/ L in is characterized by a weak diurnal variation with an afternoon maximum (Fig. 3). Some differences are observed between the diurnal cycles in 2019 and 2020. In 2019, the net longwave radiation was negative in the evening and at night, and close to zero in the daytime. In 2020, it was close to zero in the evening and at night, and positive in the daytime. This indicates a positive contribution of longwave radiation to the glacier melting in 2020 and a negative contribution in 2019. Accordingly, at night and in the evening, the net radiation (R net ) of the glacier surface in 2019 was negative (the glacier was losing energy), while in 2020, it was close to zero. In the diurnal cycle, the turbulent fluxes of both sensible and latent heat were always positive, with an afternoon maximum (Fig. 3). More pronounced maxima of sensible and latent heat were observed in 2020. In contrast to 2019, the night and evening turbulent fluxes in 2020 fully compensated the energy loss by longwave radiation. In general, during the ablation season, stable (inversion) air conditions prevailed both during the day and at night.
Surface energy balance and melting. The glacial surface was melting for almost the entire observation period, with the exception of night frosts, when all available heat was spent on heating the glacier surface to the melting point (4% of cases in 2019 and 7% in 2020). That is, almost all the heat supplied to the surface was spent on melting snow and ice. The main source of melting energy was R net (89-95%), followed by turbulent heat (H + LE, 5-9%) ( Table 2, Fig. 5). The energy input with rainfall Q r varied insignificantly (daily values from 0 to 6 W m −2 in 2019, and from 0 to 9 W m −2 in 2020) and was an insignificant positive component of the SEB (1-2%). The subsurface heat flux Q g contributed negatively to the SEB and was insignificant (< 1%). R net was dominated by the net shortwave radiation S net , but its contribution to the SEB differed in 2019 and 2020. If in 2019 S net was 115% and L net − 20%, then in 2020, S net was 85% and L net 4%. That is, if in 2019 the glacier lost the energy for melting due to longwave radiation (− 15 W m −2 ), then in 2020, it received it (4 W m −2 ). S in greatly influenced the day-to-day changes in S net , which is confirmed by their high correlation (in 2019, R 2 = 0.72, n = 48).
The daily values of the SEB were positive on all days, with the exception of 18 August 2020 (Fig. 5). They varied from 0 to 172 W m −2 in 2019 (on average 78 W m −2 ), and from − 20 to 188 W m −2 in 2020 (on average 79 W m −2 ). In 2019, the SEB in August was slightly higher than in July, while in 2020, on the contrary, it was lower. Despite the higher S net values in 2019 compared to 2020 (90 and 67 W m −2 , respectively), the mean SEB values in 2019 and 2020 were the same (79 W m −2 ). This is explained by the positive net longwave radiation (4 W m −2 ), increased turbulent heat fluxes (H + LE, 7 W m −2 ) on the background of increased cloudiness, and higher wind speeds in 2020 (Table 1).

Ablation.
In 2019 (48-day period, 7 July to 23 August), the total measured mass balance (ablation) was − 1.38 m w.e. (29 mm w.e. day −1 ), with July and August accounting for − 0.85 and − 0.53 m w.e. (61 and 39%), respectively. The mean daily ablation in July and August differed by 1.4 times (34 and 24 mm w.e. day −1 , respectively). By the time the snow cover completely melted (4 August 2019), the total ablation was − 1.05 m w.e. (76% of the total ablation). For the same period in 2020, the total ablation was − 1.47 m w.e. (31 mm w.e. day −1 ). The total and daily average ablation in July were slightly higher than in August (− 0.80 and − 0.63 m w.e. and 32 and 29 mm w.e. day −1 , respectively).
The modeled and measured 48-day mass balances were in good agreement both in 2019 and 2020 (Fig. 6). In general, the ablation period was characterized by a more or less uniform decrease in the glacier mass balance. The root-mean square error (RMSE) of the model is 0.01 m w.e., which is within the accuracy of the ablation

Discussion
To test the sensitivity of the modeled SEB and mass balance to changes in meteorological parameters, we repeatedly ran the model with the changed input meteorological parameters: S in and S ref (± 50 W m −2 ), air temperature (± 1 °C), relative humidity (± 10%), and wind speed (± 0.3 and + 1.0 m s −1 ). The tests showed (Fig. 7) that with an increase in air temperature by 1 °C, the total ablation will increase by only 0.014 m (1.0%), and with a decrease in temperature by 1 °C, ablation will decrease by 0.029 m (1.9%). With an increase in relative humidity of 10%, ablation will increase by 0.01 m (0.7%), and with a decrease, it will decrease by 0.031 m (2.1%). With an increase in wind speed by 0.3 m s −1 , ablation will increase by 0.042 m (2.8%), and with a decrease by 0.3 m s −1 , it will decrease by 0.041 m (2.7%). At the same time, with an increase in wind speed by 1 m s −1 , ablation will noticeably increase by 0.195 m (13%). In our calculations, changes in the air temperature, humidity, and wind speed only affected the turbulent component of the SEB. Sensible heat is sensitive to changes in the wind speed and air temperature, while latent heat is sensitive to changes in the temperature, humidity, and wind speed. An increase in wind speed has the strongest effect on the magnitude of turbulent fluxes and their proportion in the melt energy. For example, with an increase in wind speed by 1 m s −1 , the turbulent heat fluxes (H + LE) will increase by 3 times, and their proportion will increase to 21%. However, changes in the shortwave radiation have the greatest influence on glacier ablation. With an increase in S in by 50 W m −2 , the total ablation would increase by 0.664 m (44%), and with a reduction, it would decrease by 0.536 m (36%). With an increase in reflected radiation by 50 W m −2 , melting would decrease by 0.65 m (43%), and with a decrease, accordingly, it would increase by 0.453 m (30%). S ref has the greatest effect on changes in albedo (from − 0.17 to 0.25, respectively). In general, the summer mass balance of the glacier is most sensitive to changes in the radiation components of SEB (S in and S ref , albedo). We compared the SEB and summer meteorological conditions of the Sygyktinsky Glacier with those of some other mid-latitude (40-65°N) glaciers of Eurasia (Table 3). As the data were obtained in different ablation seasons and using different approaches, the analysis is only very preliminary. However, we can see some similarities and differences in the SEB components of the glaciers and their relations with meteorological conditions. In general, S in decreases with increasing latitude. The relatively low average value of S in at the Kodar (< 200 W m −2 ) is close to that on the glaciers of Southern and Eastern Siberia (Altai, Suntar-Khayata), as well as on the Storbreen glacier in Norway. This is probably due to significant cloudiness, which has maximum values in Kodar, Altai, and Norway (> 70%). Accordingly, the net shortwave radiation on these glaciers (with the exception of Altai) has rather low values (< 100 W m −2 ). All glaciers lost radiative heat by outgoing longwave radiation (L out ), while due to high cloudiness, the smallest losses were observed on the Storbreen and Sygyktinsky glaciers (< 10 W m −2 ). The contribution of R net to the SEB at Kodar (92%) was slightly higher than on the glaciers in the Alps, Scandinavia, Altai, and Suntar-Khayata (76-88%). Accordingly, the contribution of turbulent heat to the melt energy at the Sygyktinsky Glacier was minimal (8%), while the glacier received additional heat both in sensible and latent (condensation) form. Such an insignificant proportion of turbulent heat on the Sygyktinsky Glacier, despite the high air temperature, is explained by extremely low wind speeds (< 1.5 m s −1 ), the lowest among all compared glaciers (Table 3). In turn, low wind speeds on the glacier are probably explained by the predominance of the low gradient baric field over the Eastern Siberia in summer 30 . It is known that under stable air stratification (dT > 0) at low wind speeds, turbulent heat transfer is suppressed as a result of the action of hydrostatic stability 31 . Studies www.nature.com/scientificreports/ in Greenland have shown that aerodynamic stability reduces the sensible heat flux over the melting ice surface compared to that predicted for a stable boundary layer 21 . Thus, changes in such meteorological characteristics as cloudiness and wind speed most strongly affect the structure of the SEB and the summer mass balance of the Sygyktinsky Glacier.

Conclusions
Using meteorological data, we computed all energy fluxes on the Sygyktinsky Glacier, Kodar ridge, Eastern Siberia, during two ablation seasons (July-August 2019 and 2020). The data obtained made it possible to calculate the physically based energy-balance model of the glacier. The model includes both the direct measurement of radiation fluxes and calculated turbulent fluxes (aerodynamic method with stability correction). The model is in good agreement with the measured summer mass balance. Sensitivity tests have shown that the SEB is most sensitive to changes in shortwave radiation and weakly sensitive to changes in temperature and humidity. Meanwhile, the turbulent heat fluxes are largely controlled by wind speed. With an increase in wind speed, the proportion of turbulent heat in the SEB increases sharply. Net shortwave radiation is a dominant source (89-95%) of melt energy, and it is strongly controlled by cloudiness and short-term summer snowfalls (in August). The net longwave radiation is weakly negative/positive, depending on the prevailing weather conditions (cloudiness). Sensible and latent turbulent fluxes are positive components of the SEB. However, their contribution is insignificant (≤ 10% balance) due to the low wind speed (on average 1.2 m s −1 ). Nevertheless, on some windy days, turbulent fluxes can reach quite large values. The heat of precipitation and subsurface fluxes are insignificant over the ablation period and compensate for each other. The SEB on the Sygyktinsky Glacier is similar to that on the other mid-latitude glaciers in South and Eastern Siberia, as well as on the Storbreen glacier in Norway. We explain it by the dominant overcast weather conditions during the ablation period. The prevailing influence of R net in the SEB indicates the need to take into account changes in atmospheric processes (cloudiness, summer precipitation, wind speed) when analyzing long-term glacial trends.

Methods
Automatic weather stations. Two automatic weather stations (AWSs) were installed on 6 July 2019 on the glacier and its terminal moraine (Fig. 1). AWS1 was installed on a relatively flat surface (≤ 10) near the ice divide of the glacier at 2561 m a.s.l. (56°51.02′ N, 117°25.09′ E), near the long-term average equilibrium line altitude (Table S1). Temperature, humidity, and shortwave radiation (incoming and reflected) sensors were mounted on a vertical mast. As the glacier surface descended, the mast was thrice redrilled; accordingly, the height of the sensors above the glacier surface was corrected and the orientation of the solar radiation sensors was checked. The glacier is quite narrow (about 300 m, Fig. 1) and avalanches falling down to its surface in winter and spring are common. Therefore, to prevent a risk of losing the expensive sensors of incoming and outgoing longwave radiation we permanently installed them on the moraine station (AWS2), at a vertical mast at height 2.5 m above the ground surface (Table S2). A thermistor cable with temperature sensors (a distance of 10 cm www.nature.com/scientificreports/ away from each other) was installed in the 2.2 m borehole drilled 3 m from AWS1 to measure the snow/ice temperature and daily ablation. In the winter of 2019/20, AWS1 failed due to an avalanche and thermistor cable was destroyed and the mast with sensors was buried under the snow. On 6 July 2020, a tripod was installed at the same place with the same sensors for air temperature and humidity (at heights of 0.5 and 2.0 m above the glacier surface), wind speed and direction (1.0 and 2.0 m), and S in and S ref (2 m). Measurements of air temperature, humidity, and wind speed at two levels were carried out in order to more accurately estimate the turbulent heat fluxes and the roughness length parameter. The tripod was descending with the glacier surface, while the height of the sensors above the glacier surface always remained constant. AWS2 was installed on the terminal moraine, as a rigid tripod, at an elevation of 2529 m a.s.l. (56°50.84′ N, 117°25.06′ E), in order to calibrate AWS1 records, and to be a backup in the event of an unexpected AWS1 failure (e.g., due to avalanches). AWS2 included sensors for air temperature and relative humidity, atmospheric pressure, wind speed and direction, shortwave and longwave radiation (incoming and outgoing), and ground temperature (at the surface and at a depth of 10 cm). Here, we assumed that the incoming longwave radiation on the moraine (AWS2) was similar to that on the glacier (AWS2). This assumption is based on the following premise. AWS1 and AWS2 are located close to each other in space (300 m in distance and 30 m in altitude). This means that they are located within the same boundary layer of the atmosphere, which is often estimated at values of the order of 50-100 m. Accordingly, this assumption is confirmed by similar meteorological parameters, temperature and humidity, measured at both stations (R 2 for temperature and relative humidity are 0.89 and 0.91, respectively). On average, in the summer of 2019, the difference in temperature and humidity between AWS1 and AWS2 was only 0.3 °C and 1%, respectively, which is within the measurement error.
The temperature and humidity sensors were naturally ventilated and protected from solar radiation. The sensor readings were recorded in the memory with a 30 min frequency. The instantaneous wind speed was measured every 30 min, while the maximum wind speed was recorded within each 30-min interval. Irkutsk time (GMT + 8) was used as the closest to local time. The AWS data were loaded onto removable memory cards during field observations on 24 August 2019 and 2020. Due to a save error, wind direction data for 2019 was lost. Thus, continuous time series were obtained for the period from 7 July to 23 August during two seasons of ablation (this period is analyzed in the article). At the beginning and end of the observation period, the density of snow and ice was measured in shallow pits near AWS1. Ablation was also measured in the vicinity of AWS1 using 11 plastic rods. The characteristics of the sensors used in both AWSs are listed in Table S2.
Data treatment. The raw AWS data were thoroughly checked for errors using basic quality control procedures recommended by the World Meteorological Organization (WMO) 34 . In total, ≤ 1% of the original data (solar radiation, temperature, humidity, and wind speed) were rejected. Erroneous sporadic data were removed, and gaps were filled using linear interpolation. The main errors were typical for solar radiation data and caused by the effect of atmospheric precipitation on the upper sensors (especially in mixed or solid form). Since the rain gauge measured only liquid precipitation correctly, 13% of the precipitation data was rejected without filling.
Albedo was calculated as the "accumulative albedo" α acc with a 30-min resolution as the ratio of the sums of reflected and incoming shortwave radiation over a 24-h time window 20 : The use of "accumulative albedo" instead of the conventional one makes it possible to neutralize the measurement errors of the incoming shortwave radiation 20 . The specific air humidity q (g kg −1 ) was calculated using formulas adopted by the WMO 35 with the measured values of the atmospheric pressure p (hPa), air temperature T (°C), and relative humidity RH (%) and the calculated value of the water vapor pressure in humid air e (hPa) as inputs: Surface energy balance model. The surface energy balance (SEB) was computed for 30-min intervals as where S in and S ref are the incoming and reflected shortwave radiation, L in and L out are the incoming and outgoing longwave radiation, H and LE are the sensible and latent heat, Q r is the heat supplied with rain, and Q g is the subsurface heat flux. S in and S ref were directly measured at the glacier (AWS1) while L in was measured at the moraine (AWS2). Since L out was not measured at AWS1 we assumed a constant value (316 W m −2 ) for a melting glacial surface (0 °C) 11 . All terms are taken to be positive when directed towards the surface and expressed in W m −2 .
The surface melt M (mm day −1 ) was calculated as

Turbulent fluxes.
A bulk aerodynamic approach based on the Monin-Obukhov similarity theory, including stability correction, was used to calculate turbulent fluxes. This approach is used when calculating turbulent flows from micrometeorological measurements on glaciers in different climatic settings 8,14,17 . This method shows a good correlation of calculated turbulent flows with those measured with eddy-covariance systems 7 . The sensible and latent turbulent fluxes, H and LE, were calculated as where c p is the specific heat capacity for air at constant pressure (1010 J kg −1 K −1 ); ρ 0 is the air density at the standard sea level (kg m −3 ); p and p 0 are the atmospheric pressure at glacier and standard sea levels (hPa); k is the von Karman constant (0.38); u , T, and e are the wind speed (m s −1 ), air temperature (K), and water vapor pressure in humid air (hPa) at a measurement level z above the glacial surface (2 m); T s and e s are the glacier surface temperature (273.15 K for a melting surface) and water vapor pressure (6.11 hPa) at the glacial surface level z 0 at 0 °C; L v is the latent heat of evaporation of snow/ice for the melting surface (2514 × 10 3 J kg −1 ); and z 0 m , z 0 t , and z 0 h are the surface roughness lengths (m) for moment, temperature, and humidity, respectively. Dimensionless stability functions for the moment ( m ), temperature ( t ), and humidity ( h ) were calculated using the bulk Richardson number Ri b 17 .
Ri b is a dimensionless characteristic of the relationship between the thermal and dynamic factors of turbulence and was calculated between the measurement level z (2 m) and the glacier surface z 0 as 17 where g is the acceleration due to gravity (9.8 m s −2 ).
Taking into account that at very low wind speeds Ri b takes unrealistically high values (> 1), we used only wind speeds > 0.5 m s −1 . In addition, we adopted a critical value of Ri b equal to 0.4, assuming that when this value is exceeded, ( Ri b > 0.4), turbulence ceases completely (H = LE = 0) and the flow becomes laminar 36 . Thus, the influence of the thermal stability of the near-glacial air layer on the turbulence coefficient was taken into account. In unstable air (a decrease in temperature with height), Ri b < 0 and the turbulence coefficient is greater than that in stable air, while in stable (inversion) air, Ri b > 0 and the turbulence coefficient is less than that in unstable air. Stable (inversion) conditions dominated on the glacier during the entire observation period. For example, in summer 2020, stable air conditions ( Ri b > 0) were observed in 91% of cases, while unstable conditions ( Ri b < 0) were observed in only 9% of cases.
The roughness lengths for the momentum z 0m were calculated from the wind speed measurements at upper (2 m) and lower (1 m) levels at AWS1 during the 2020 ablation season. For this, we used data under the nearneutral air conditions (Ri ≈ 0). When using the criterion |Ri|< 0.1 (17% of observations), the median value of z 0m was 0.61 mm (6.1 × 10 −4 m). The roughness lengths for the temperature z 0t (4.3 × 10 −4 m) and humidity z 0h (4.9 × 10 −4 m) were calculated using the calculated z 0m and the roughness Reynolds number Re * (here, we used Re * = 2.5) in accordance with the frequently used parameterization 37

Scientific Reports
| (2021) 11:21260 | https://doi.org/10.1038/s41598-021-00749-x www.nature.com/scientificreports/ where k τ is the thermal conductivity (0.4 W m −1 K −1 for old snow and 2.2 W m −1 K −1 for ice), and T g is the temperature (K) of the upper glacier layer at a depth z g (m) from the surface. For the calculation, we used the data of the temperature-depth profile obtained only in the 2019 ablation season (in 2019/20 winter, the thermistor cable broke down due to the avalanche and did not work). According to these data, the temperature of the glacier at a depth of 2 m was − 1.7 °C.
Ablation measurements. Continuous ablation measurements were made using the thermistor cable installed in a 2.2 m borehole next to AWS1 with a daily resolution for the period from 6 July to 22 August 2019. The technique of automatic measurement of ablation using a thermistor cable (vertical arrays of temperature sensors) was recently successfully tested on an Italian glacier 39 . The technique is based on the ability to differentiate between sensors buried in ice/snow and sensors exposed to the atmosphere. Sensors located beneath snow cover and within ice show much lower temperature variability compared to exposed sensors due to solar radiation. Thus, the analysis of temperature data makes it possible to identify the location of a sensor (below, above or on the ice/snow surface). Here we used the daily temperature variance calculations for all sensors to estimate their positions relatively glacier surface. The distance between the temperature sensors (10 cm) made it possible to measure the lowering of the glacier surface with a standard error of ± 5 cm 39 . In addition, ablation was measured in the vicinity of AWS1 by 11 stakes (3 readings per season). In 2020, due to the failure of the thermistor cable, aperiodic ablation measurements were made only by stakes (6 readings per season).