Regional snow-depth estimates for avalanche calculations using a two-dimensional model with snow entrainment

Abstract

The currently adopted approach for avalanche-hazard mapping in Switzerland includes avalanche-dynamics modelling coupled with estimation of the greatest annual 3 day snowfall depth, H₇₂, for high-return periods, used as the release depth. New advances in avalanche dynamics show that this approach can be improved using models with mass entrainment, requiring in turn a statistical definition of the erodible snow cover. We propose a regional approach, based on index value, to evaluate release depth and erodible snow cover for large-return periods. The territory of Switzerland is divided into seven climatologically homogeneous regions. Generalized extreme value (GEV) distributions for the growth factors coupled with index-value estimation based on altitude provide an accurate estimate of snow depths, also for large-return periods. RAMMS, a two-dimensional avalanche-dynamics model including snow entrainment, is used for hazard mapping for a site used as an example of the Swiss procedure. The regional approach allows the boundary conditions for hazard mapping to be set using an entrainment model, and also provides statistical uncertainty of the design release and erosion depth, thus aiding in applying uncertainty analysis to hazard-mapping procedure.

1. Introduction

Mitigation of avalanche risk involves land-planning measures enforcing use restriction for some territories. These require hazard maps, often based on design avalanche with given return period T. In Switzerland, hazard map assessment has to be carried out using the Swiss procedure (Sp), which defines high-hazard red zones (impact pressure p ≥ 30 kPa and T ≤ 300 years, or p ≤ 30 kPa and T ≤ 30 years) and moderate-hazard blue zones (impact pressure p ≤ 30 kPa and 30 ≤ T ≤ 300 years). The Sp couples statistical analysis of snow depth in the release zone, equated with the 3 day snowfall before the event, H ₇₂, with numerical calculation to predict runout distance and impact pressure of the avalanche (BFF/SLF, 1984; Salm and others, 1990). Therefore the Sp requires estimation of H ₇₂ for T = 30 and T = 300, for which local distribution fitting may not be accurate. Unfortunately, the Sp neither provides a method to evaluate uncertainties therein, nor indicates how to deal with it in hazard mapping. To draw accurate 300 year H ₇₂ estimates for avalanche-hazard mapping, one needs about 150 years of data (e.g. Katz and others, 2002). This problem can be overcome using a regional approach based on index value (e.g. Bocchiola and others, 2006; Bocchiola and Rosso, 2008), which evaluates the frequency distribution across all sites within a given homogeneous area or region (Hosking and Wallis, 1993). Key features of this approach are the gain in accuracy due to increased sample dimensionality, together with well-assessed estimation uncertainty, according to a widely known literature concerning regional evaluation of extremes (e.g. De Michele and Rosso, 2001, and references therein). In this paper, an index value approach is applied to the territory of Switzerland. The avalanche-dynamics models proposed by the Sp (Salm and others, 1990; Bartelt and others, 1999) sketch the avalanche mass only considering the release zone. Recent works demonstrate that including mass entrainment along the avalanche path may improve our knowledge of the avalanche process as well as its simulation (Sovilla and others, 2006, 2007). In this paper, we perform avalanche simulations with the two-dimensional (2-D) RAMMS:: avalanche model (e.g. Christen and others, 2007), with mass entrainment. Previous studies on snow entrainment show that it is governed by along-track snow mass availability as well as by snow-cover structure, while topographic features and flow variables are of less importance (e.g. Sovilla and others, 2006). Therefore, evaluation of erodible snow cover is of primary importance in modelling mass entrainment. The results of the regional analysis are here extended to the erodible snow cover, providing the boundary conditions for mass-uptake calculation. We use as a case study the Ariefa/ Samedan avalanche, taken from the Sp (Salm and others, 1990). Simulations with and without entrainment are carried out and compared. An analysis of the sensitivity of the mapped runout zone to the uncertainty in evaluation of H ₇₂ is then carried out.

2. Regional Approach for Frequency Estimation of H₇₂

2.1. Study area and available database

The case study is the Swiss mountain region, where over 200 snow-gauging stations have been installed by SLF and by the Swiss Meteorological Institute. For the present study we used only data collected by the manual stations network (Fig. 1), featuring more than 10 years of measurements. We identified 116 stations, with altitude ranging from 660m at Bellinzona station (6BE) to 2910m at Felskinn station (4FK) (e.g. Bianchi Janetti and Gorni, 2007).

Fig. 1. Case-study area and homogeneous regions.

2.2. Regional estimation of H₇₂

The regional approach is widely adopted in the field of statistical hydrology (e.g. Burn, 1997; Katz and others, 2002). The basic assumption is that values of a hydrological variable observed at a specific site i, scaled by an index value, have identical frequency distributions across all the sites within a given homogeneous region (e.g. Hosking and Wallis, 1993). The first step is the identification of homogeneous regions, here in terms of the variable H ₇₂. Traditionally, Switzerland is divided into seven climatological regions, sketched considering morphological features and political boundaries (see, e.g., Laternser, 2002; Laternser and Schneebeli, 2003). Several authors investigated different regionalization against morphoclimatical features. Laternser (2002) grouped Swiss snow stations applying a cluster analysis to daily snow-depth data and proposed a hierarchical division up to 14 regions. Building on that study and using an iterative approach based on statistical tests (as in, e.g., Bocchiola and others, 2006) for regional homogeneity and homogeneous dependence of H ₇₂ upon altitude, seven homogeneous regions were identified (see Bianchi Janetti and Gorni, 2007). Their boundaries were drawn by including the related gauges, following as much as possible the boundaries of the watersheds and the river valleys, which provide natural divides. The northern slope of the Swiss Alps is divided into four independent areas:

region 1, the northeast belt crossing the country from west to east;

region 2W, the western part of Valais;

region 2E, the eastern part of Valais; and

region 5, the northern part of Grison.

The southern part is divided into three independent areas:

region 3, the southern part of Valais;

region 4W covering Ticino; and

region 4E covering the eastern part of Grison.

The index value approach can be adopted for evaluation of H ₇₂ in each region. The index value is estimated by the single-site sample average (for a discussion on the use of different index values (e.g. the median), see Bocchiola and others, 2003), at each specific site i:

(1)

where Yi is the number of years of observation and the suffix y indicates the yth year. The related standard error of estimation is:

(2)

with σ_{H 72} i the sample standard deviation of H ₇₂ at the specific site i. The scaled value of H ₇₂ at each specific site i is therefore defined as:

(3)

The symbol F indicates the regional cumulated probability distribution of H ₇₂ ^*, valid at each site inside the homogeneous region. The average of H _72i ^* is obviously 1, and the remaining moments need to be estimated from data and are the same at each site of the same region. For the seven considered regions, the sample frequency of H ₇₂ ^* is well accommodated using a generalized extreme value (GEV) (including GEV1, or Gumbel) distribution (e.g. Kottegoda and Rosso, 1997), which is often used to model H ₇₂ (e.g. Barbolini and others, 2002). The GEV quantiles featuring T-years return period are evaluated as:

(5)

with ε, α and κ the location, scale and shape parameter, respectively, and y_T the Gumbel variable: y_T = –lnf–ln[(T–1)/T]g. For κ = 0 (i.e. Gumbel distribution) the quantiles featuring T-years return period are evaluated as:

(5)

The parameters of the considered distributions are shown in Table 1 for the seven regions and for Switzerland as a whole.

Table 1. GEV parameters for Switzerland and for regions 1–5. n _tot is number of equivalent years

These are estimated using L-moments, i.e. linear combinations of the probability weighted moments (PWM) commonly adopted to evaluate distribution parameters (e.g. Hosking, 1990). L-moments are more convenient than PWM because they can be directly interpreted as measures of the distribution shape for regional analysis.

The confidence limits of a given quantile can be evaluated as: where Фis the quantile of standard normal distribution, N.ST(0;1) and σ_{H 72} ^* (T) is the standard deviation of the estimated value of H ₇₂ ^*(T), calculated as (De Michele and Rosso, 2001):

(6)

The value of σ_{H 72} ^*(T ) depends, for the regional case, on the whole sample size n _tot (Table 1). When single-site distribution fitting is carried out, Equation (6) still holds and n _tot is changed with the number of years of observations at the single site, Y_i . The uncertainty of the T-years single-site quantile estimation is calculated as follows (e.g. Bocchiola and others, 2006):

(7)

Notice that σ_μH ₇₂ i (Equation (2)) always depends on the number of observed years, Y_i , while the value of σ_{H 72} ^*( ^T ) depends for regional distribution fitting on the whole sample size n _tot and, for single-site distribution fitting, on the number of sampled years, Y_i . Bocchiola and others (2006) showed that the value of σ_Ti is considerably smaller for the regional case than for the single-site case, particularly for the 300 year return period, so the regional estimates are more accurate than those at single sites. Figure 2 shows, for illustrative purposes, the case of Samedan (7SD) station, the closest to the Ariefa/Samedan avalanche site, used here as a case study. Notice the very high uncertainty in quantile estimation for considerable return periods when single-site distribution fitting is carried out, as compared to the regional approach.

Fig. 2. Regional and local estimation of T-years quantiles of H ₇₂ and related uncertainty (confidence level α = 5%) for the Samedan (7SD) station within region 4E. Notice the plotting position of the local observed values of H72^* in the regional sample, showing good agreement with the regional distribution.

2.3. Evaluation of μH₇₂ for ungauged sites

For hazard mapping, an evaluation of H ₇₂ is carried out in the release zone of the avalanche, where often no snow gauges are available. Evaluation of μ_{H 72} at an ungauged site can be made based on altitude, A, as it is possibly the factor that most influences the distribution of snowfall in space (e.g. Barbolini and others, 2002; Bocchiola and Rosso, 2007). The Sp approach suggests a mean increase of H ₇₂ (for T = 300 years) of 5 cm with 100m in altitude (Salm and others, 1990). A more detailed estimate of the index value is obtained here, using in each region a proper linear regression against altitude, i.e.

(8)

with A _un the altitude of the ungauged site and μ ₀ the intercept for A = 0. Its standard deviation is

(9)

with σ̂ _E [ _{H 72]}the observed (sample) standard deviation of the index value for the sites used to build the A to μ_{H 72} equation, and R² the determination coefficient, or percentage explained variance. The details of the regression analysis are reported in Table 2.

Table 2. Relationship between μH ₇₂ and A in the proposed regions. N _s is the number of stations. Italic indicates non-significant dependence upon altitude

We found significant increase of μ_{H 72} against A in all regions but 2W and 5 (Fig. 3). For region 5, these results are consistent, for example, with the analysis by Laternser (2002). Therein, the author reports that the gradient of the new snow (strictly correlated to H ₇₂) against altitude in the Davos region (region 5 here) is not significant. Concerning region 2W, although no comparison studies were found, the very weak dependence on altitude might be explained by the observation that most of the stations in these regions are above 1300 ma.s.l. in the inner-alpine dry region, and so are not affected by strong orographic precipitation, as confirmed by SLF Davos personnel in charge of snow station management and maintenance (personal communication from C. Marty, 2007). As a comparison, by taking here the average rate of increase of μ_{H 72} for Switzerland, namely c = 2.01 cm(100m)^–1, multiplied by the 300 year quantile of the GEV distribution for Switzerland, namely H ₇₂*(300) = 2.63, one obtains c(300) = 5.28 cm(100m)^–1, consistent with the Sp. The proposed approach therefore shows that the rate of increase of H ₇₂(300) proposed by the Sp is valid on average on the Swiss territory, but a different rate of increase of snow depth H ₇₂(300) with altitude should be considered for design in different areas.

Fig. 3. Altitude to average depth μH ₇₂ relationships: (a) regions 1–3–5; (b) regions 2E–2W; and (c) regions 4E–4W.

3. A Case Study: the Ariefa/Samedan Avalanche

3.1. Available database

As a case study of hazard mapping using the regional approach, we use the Ariefa/Samedan avalanche in canton Grison, Switzerland. The avalanche is located inside region 4E. The greatest avalanche observed at this site occurred on 21 January 1951. The Ariefa avalanche is one of the four calculation examples of the Sp (Salm and others, 1990; Bartelt and others, 1999) and is labeled with an approximate return period of 300 years. It is documented by an indicative release volume and observed runout distance and an approximate avalanche width. The release width was estimated to be 100–180 m, and the release area (Fig. 4) was located between 2000 and 2340 ma.s.l. The release volume was estimated to be V ₀ ≈ 8104^m3. The avalanche flowed on an open slope for about 2000 m. The average slope was 25^˚ in the upper part and 8^˚ in the runout zone.

Fig. 4. Ariefa/Samedan avalanche. Calibration of RAMMS 2-D with entrainment. Release area (light polygon) and simulated maximum flow heights (scale of gray). The dark line shows the main profile.

3.2. RAMMS 2-D model calibration

We used the 2-D numerical model RAMMS::avalanche (Christen and others, 2007), including mass entrainment. The model requires the definition of the release volume, i.e. of the release area and the fracture depth. The release area A ₀ ≈ 8104^m2 is defined with the aid of a semi-automatic procedure. A Geographic Information System program individuates all the possible release areas in a certain domain, i.e. all the cells characterized by a slope of 30–50^˚ (Maggioni and Gruber, 2003). According to the Sp, statistical definition of the 300 year release fracture depth is necessary to evaluate the runout distance, velocity and impact pressure of the avalanche (e.g. Bartelt and others, 1999). Here we obtained a regional growth factor H ₇₂ ^* = 2.85, as from calculation with the Gumbel distribution of region 4E for a return period of 300 years, and an index value μ̂_{H 72} = 0:58 m, calculated at the release-zone altitude by Equation (8). These values result in H ₇₂ = 1.64 m. To calculate the 300 year release fracture depth, we corrected this value for the slope in the release zone, obtaining h _r = 1.04 m. This provides a value very similar to that indicated by the Sp, i.e. 1.00m (e.g. Bartelt and others, 1999). Since the model includes entrainment, it is also necessary to define the snow-cover depth H _e along the avalanche track (Sovilla and others, 2006). The erodible snow-cover H _e is here evaluated according to the regional analysis, i.e. coinciding in the release zone with h _r and decreasing by 7 cm every 100m (rate of increase with altitude for region 4E for a return period of 300 years). The RAMMS model has four calibration parameters: the Coulombian friction μ, the turbulent friction ξ, the active/passive pressure coefficient and the threshold pressure value for erodible snow-cover erosion, τ (see Sovilla and others, 2006, for full explanation of a one-dimensional (1-D) entrainment model based on use of τ, also implemented in RAMMS). While the μ, ξ and parameters and their influence on avalanche simulation results have been thoroughly assessed (e.g. Buser and Frutiger, 1980; Salm, 1993; Bartelt and others, 1999), less is known about the τ parameter, critical for modelling snow entrainment (Sovilla and others, 2006, 2007). Here we introduce a simplified approach (explained extensively in Bianchi Janetti and Gorni, 2007) to evaluate τ against the size of the avalanche (i.e. h _r, A ₀ and H _e) and its dynamic properties (i.e. μ and ξ).

A preliminary τ calibration was carried out by back calculation of two case-study avalanches in Vallée de la Sionne, Canton Valais, Switzerland, (Sovilla, 2004; Sovilla and others, 2006) where avalanche depth, volume and entrainment data were available, as well as runout and along-track snow depth and velocity data. We also sketched a simplified topography of the flowing zone of the Vallée de la Sionne avalanche (i.e. an inclined plane with a slope of 30˚; for a similar approach see, e.g., Sovilla and others, 2007), and simulated several synthetic avalanche events with changing values of μ, ξ, h _r, A ₀ and H _e. Optimal values of τ were then evaluated by eyeball assessment of the likelihood of the eroded area morphology, according to the observed events. In so doing, we obtained a rule to evaluate τ as a regular function of the mentioned variables (Bianchi Janetti and Gorni, 2007). Although this approach is approximated for a number of reasons (e.g. change in avalanche track morphology, snow properties), it provides a rule to sketch a reasonable value of τ when no entrainment data are at hand.

This approach was then used to calibrate the model for the Ariefa/Samedan avalanche. The parameters μ, ξ and τ were evaluated so as to match the observed runout, both without (no τ) and with entrainment. In the simulation with entrainment, μ and ξ were tentatively set and τ calculated accordingly, until the historical end-mark was reached. For all simulations we set snow density _r = 300 kgm^–3 and erodible snow-cover density _e = 200kgm^–3. These are standard values adopted in Switzerland when applying the Sp procedure (Salm and others, 1990); they have been observed in field tests (Sovilla and others, 2006) and are also suggested for avalanche-hazard mapping in Italy (Barbolini and others, 2004).

A summary of the main results is given in Table 3. The values of μ and ξ without entrainment are slightly different from those used by Bartelt and others (1999) for the same avalanche. However, they used a 1-D model, requiring, in practice, different coefficients with respect to the present model RAMMS 2-D, even without entrainment. Figure 5 depicts the maximum flow heights along the avalanche main profile. With entrainment, the avalanche increases its mass by 2.1 times. The mass is entrained at the avalanche front, leading to a maximum flow height of 8.0 m, almost twice as much as in the no-entrainment simulation (Fig. 5a). Figure 6 depicts the greatest velocities along the main profile for all simulations. Those calculated with entrainment in the first part of the track are lower than those calculated without entrainment (Fig. 6a). Entrainment initially leads to decreased velocity, due to acceleration of the entrained mass up to the avalanche velocity. Then the entrained mass increases, thus resulting in higher velocity in the second part of the track (cf., e.g., Sovilla and Bartelt, 2002). After 1400m from the release, the velocity of the entrainment simulation exceeds that of the no-entrainment simulation. The erosion occurs in the upper and steeper part of the track, where the snow cover is deeper. Bartelt and others (1999) predicted a greatest flow depth of about 4 m, similar to that obtained here without entrainment. The greatest velocity predicted by Bartelt and others (1999) is U _max = 37ms^–1, reached in the vicinity of the P point marking the beginning of the runout zone. In the no-entrainment mode, RAMMS here predicts 38.6ms^–1, dropping to 37.2ms^–1 when entrainment is considered. Also reported in Table 3 are the results of a simulation with entrainment carried out with RAMMS using the same values of μ and ξ as in the no-entrainment simulation. This is to evaluate the influence of entrainment on the simulation. Notice the increased runout distance and the greatest flow depth and velocity, due to mass uptake.

Fig. 5. Comparison of the predicted maximum flow heights along the main profile. (a) Entrainment/no-entrainment calibration against the observed runout. (b) Uncertainty level for simulations with entrainment.

Fig. 6. Same as Figure 5, but for predicted maximum velocities.

Table 3. Model parameters and simulation inputs and results. Runout is calculated after point P, which marks the runout zone (Fig. 5). Italic indicates tuning of the model parameters to match the historical end-mark

3.3. Hazard-mapping procedure

Hazard-mapping procedure is carried out here with and without entrainment, to evaluate the influence of entrainment on the resulting maps. Uncertainty analysis is also shown, based on H ₇₂. Full uncertainty analysis of hazard maps might involve uncertainty of the dynamic model parameters, displacement area and H ₇₂, among others (e.g. Barbolini and others, 2002; Sovilla and Bartelt, 2002; Ancey and others, 2003, 2004). Here we focus only on the uncertainty involved in the definition of H ₇₂, as compared to its definition using the Sp, to illustrate the effectiveness of the regional approach. For hazard mapping, it is necessary to calculate the avalanche impact pressure, p = 0:5Cv ², where C is a pressure coefficient, set here to C = 2, as recommended in the field of avalanche hazard mapping, according to the Sp (Salm and others, 1990), and = 300 kgm^–3 is the snow density inside the avalanche body. The results of the uncertainty analysis are reported in Figures 5–7. Considering T=300 years, the standard deviation of H ₇₂ is σˆ _{H 72} = 0.15 m, which seems relatively small because the analysis is based on a large number of data. To evaluate the influence of uncertainties in the model output, two more simulations are performed, with the following snow heights:

Fig. 7. Variability of the red and blue zones for hazard mapping of Samedan/Ariefa avalanche.

Entrainment –σ. H ₇₂ – σˆ _{H 72} , and the erodible snow cover according to the regional distribution, that is, coinciding in release zone with H ₇₂ – σ̂ _{H 72} and decreasing 7 cm every 100 m (rate of increase with altitude for region 4E for a return period of 300 years). Friction parameters are ξ = 2500ms^–2, μ = 0.15, = 2.5 and τ = 195 kPa, calculated according to the new snow depth.

Entrainment +σ. H ₇₂ + σ̂ _{H 72} , and the erodible snow cover according to the regional distribution, coinciding in release zone with H ₇₂ + σ̂ _{H 72} and decreasing 7 cm every 100 m. Friction parameters are ξ = 2500 ms^–2, μ = 0.15, = 2.5 and τ = 220 kPa calculated according to the new snow depth.

Uncertainties in estimation of H ₇₂ only slightly influence maximum flow heights (Fig. 5b), while providing runout distances varying between 837 and 915 m. The velocity distribution also shows changes due to the statistical uncertainties (Fig. 6b). Once the avalanche pressure is estimated, it is possible to evidence the change of the red and blue zones for the simulations (Fig. 7). Notice that in all cases with entrainment the red zone is moved forward with respect to the no-entrainment case (from 25m to about 100m for H ⁷² – σ̂ _{H 72} and H ₇₂ + σ̂ _{H 72} , respectively).

4. Conclusions

We sketched a division of Switzerland into seven homogeneous regions for the purpose of regional evaluation of the 3 day snowfall depth for avalanche hazard mapping, which has not previously been carried out in Switzerland. The method accounts for regional features of the GEV distributions and considers explicitly the rate of increase of H ₇₂ with altitude, both changing from region to region, contrary to what is suggested by the Sp. The method provides well-assessed uncertainty bounds, which are not available through the Sp. More accurate design of H ₇₂ is attained, because the method is based on a larger number of data than traditional single-site analysis. The regional approach provides design values for both depth at release and erodible snow cover for high-return periods, that can be used to feed a 2-D avalanche dynamic model with entrainment, thus achieving a more accurate description of avalanche dynamics. A case study of avalanche-hazard mapping is shown for one of the example sites provided by the Sp. The regional evaluation of snow-cover height allows simulation of avalanche dynamics with entrainment and definition of uncertainty bounds of the snow-height input and is therefore a valuable tool for hazard-maps design based on uncertainty analysis.