Reliability-Based Design of Protection Net Fences: Influence of Rockfall Uncertainties through a Statistical Analysis

Abstract

Net fences are among the most widespread passive protective measures to mitigate the risk from rockfall events. Despite the current design approach being based on partial safety factors, a more efficient time-dependent reliability approach has been recently introduced by the authors. This method is influenced by various parameters related to the geometry and the kinematics of the block, i.e., the uncertainty related to the distribution of the size of the impacting block, its occurrence probability, and the shape of the right-tail of the distributions of its velocity and trajectory height at the location of the net fence. Furthermore, the block size distribution of the deposit greatly affects the results. The present work focuses on the possible range of such parameters to encompass the great majority of real events. The obtained results are compared with the current design approaches based on fixed partial safety factors. It emerges that the choice of the characteristic mass of the block and the failure probability greatly influence the results. Moreover, if a set of partial safety factors is assigned to different sites, an intrinsic variability in the failure probability has to be accepted. Suggestions for an accurate procedure and future developments are provided.

1. Introduction

Rockfalls are among the most hazardous landslide phenomena, due to their abrupt occurrence and the very high involved energies [1,2,3,4]. Consequently, effective mitigation measures are required for reducing the hazard and, thus, the risk. This need especially holds for high valuable and vulnerable elements at risk as transportation routes [5,6,7], railways [8] or even urban settlements [9]. In recent decades, the knowledge and the technology on passive mitigation measures have been significantly improved, with particular reference to rockfall flexible barrier, also called net fences [10,11]. A considerable number of new devices and assembly methods have been designed, differing not only for the resistant static scheme and the way of dissipating and transferring impact loads [12,13,14], but aiming also at increasing the effectiveness of the protection and the efficiency during the expected working life [15]. Since the eighties, net fences have been considered as protection kits, whose performance has to be assessed in relation to their essential characteristics, e.g., energy absorption capacity, height, and maximum elongation [16,17]. For this purpose, a European codified method for the assessment of net fences as rockfall protection kits was developed, first in ETAG 027 [18] (”Guideline for European Technical Approval of Falling Rock Protection Kits”), and now in EAD 340059-00-0106 [19] (“Falling Rock Protection Kits”).

The common design relies on a performance based procedure. In other words, from the results of the rockfall propagation analyses, an appropriate commercial product is selected [17] among the ones whose energy absorption capacity and the nominal height, evaluated with standard procedure [19], are greater than the correspondent impact energy and passing height of the simulated blocks. This approach considers two possible failure modes: the falling block trajectory is higher than barrier height, or because its energy is larger than the barrier capacity. Nevertheless, the procedure does not directly consider the impacting conditions, e.g., the position and the angle of impact with respect of the net or the posts, as well as the shape of the impacting block [20,21]. To encompass all these uncertainties, each protection kit is classified according to energy absorption categories. More complex analyses and studies, performed through FEM (Finite Element Method) or DEM (Discrete Element Method) modelling approaches [22,23,24,25,26,27,28], provide an enhancement in understanding the dynamic interaction between the block and the structure, and the transfer load process. Nonetheless, these computationally demanding approaches offer profitable solutions for the design of more efficient devices or the assembly of the single components, resulting more useful for the producer rather than being affordable for the designer. Similarly, meta-models techniques [21,29] requiring outputs derived from either accurate numerical models or experiments are tailored on specific products and, hence, can be inconvenient for general design purposes.

At the base of any modern design code, such as the European standards EN 1990:2002 [30], the design value of the effects of the action $E_d$ has to be lower or equal to the design value of the effects of the resistance $R_d$. The standardised performance-based design procedure defined within the national recommendations UNI 11211-4 [31] in Italy and ONR 24810 [32] in Austria suggest to evaluate the design parameters of the actions from a single characteristic value of the distribution of the impacting energy (or velocity) and height, resulting from the trajectory analyses [10,33,34,35]. In the framework of the previously mentioned national recommendations, these design values are the product between characteristic values and partial safety factors. In particular, the UNI 11211-4 [31] considers the 95th percentile of the distribution for both velocity and height, and the related partial safety factors depend on the quality of the adopted trajectory model type and input parameters. The characteristic value of the impacting block mass is defined as a percentile (95th or larger) of the distribution of the masses. The ONR 24810 [32] adopts the 99th percentile for the energy and the 95th for the height, and the partial safety factors are related with the occurrence probability of a rockfall event and the consequence class of the elements at risk. A distribution of possible falling blocks has to be defined and a percentile of this distribution has to be considered as the characteristic value.

Some critical aspects in the application of the European standards EN 1990:2002 [30] for the design of rockfall net fences have been recently highlighted and fixed through a time-integrated reliability based approach [36]. Choosing a suitable failure probability, the design kinetic energy and passing height of the impacting block are determined from the results of the performed trajectory analyses. The reliability approach allows to consider the random nature of the parameters affecting the performance and the safety of a system [37] providing, in this specific case, probability distributions of the parameters of the impact and, theoretically, of the resistance. The method also accounts for the variability in the magnitude and the temporal occurrence probability of an impact, which are thoroughly site specific [1,38,39].

In a tentative to merge the time-integrated reliability based approach with the current practice, the calculation can be simplified by considering the values of the CE certified maximum energy level absorption capacity (MEL) $E_B$, the nominal height $h_B$, and the tolerance t as Dirac-$\delta$ distributions [36]. If more information derived from field tests or numerical modelling would be available, different probability distributions could be assigned to the variables. With the aim of providing a compelling solution for the designer, De Biagi et al. [36] implemented their approach into the current practice [30]: equivalent partial safety factors for the impacting block energy, mass, velocity, and height were derived, and properties and the affecting variables were analysed.

Starting from the evidence that rockfall phenomena are extremely variable in terms of magnitude and kinematics, i.e., type and direction of motion, velocities, and height, in the present paper the Authors investigate the possible scenarios which can occur in a real case. The results lead to the identification of a range of values of partial safety factors for each of the two failure modes.

The paper is structured as follow. Section 2 summarises the reliability-based design approach, considering the possible failure mechanisms. Section 3 explains the adopted procedure to investigate all the possible rockfall scenarios, while the results of the parametric analyses are introduced in Section 4. Section 5 provides an attempt to find a range of values of partial safety factors to adopt. Finally, the suggestions and future possible development are discussed in Section 6.

2. Basics of Reliability of Rockfall Protection Structures

As mentioned in Section 1, the possible failure modes can be related to the height $h_b$ and the energy $E_b$ of the impacting block, which can be higher than the barrier intercepting height $h_B$ and absorption capacity $E_B$. Lower case ${}_b$ stands for block, while capital ${}_B$ for the barrier. The failure probability of the system $p_f$ can be considered as the sum of the failure probabilities related to each of the two scenarios, i.e., $F_h$ and $F_k$, for excessive height and for excessive kinetic energy, respectively. The failure probability $F_h$ is assumed to be equal to the probability that a block has a trajectory higher than barrier height, disregarding its kinetic energy. Similarly, the failure probability $F_k$ is assumed to be equal to the probability that a block impacts with a kinetic energy greater than the barrier capacity, disregarding its height. This assumption allows to obtain a conservative value of $p_f$ and, thus, a precautionary design procedure [36].

Each failure mode accounts for a specific failure scenario, which can be mathematically described through a limit function representing the boundary between failure and safety, through which the reliability can be assessed. The two limit functions are derived as solutions, i.e., imposing to zero, of a state functions, which describe if the system is in a safe or unsafe condition. For a specific value of the reliability, i.e., a failure probability, the design values of the variables (denoted with subscript ${}_d$) composing the state function are obtained.

The proposed approach is time-integrated, considering that all the statistical properties of the variables have to be related to a time period. This is particularly important since the possible rockfall impacting impacting mass is associated to a temporal occurrence. According to De Biagi et al. [36], for a given return period T, it is assumed that the mass is normally distributed with a mean value $m_{50}\left(T\right)$ computed through a power-law rule, as suggested in [40]:

$$m_{50}\left(T\right)=M_{th}{\left(\lambda T\right)}^{1/\alpha },$$

(1)

where $\alpha$ represents the shape coefficient of a Pareto Type I distribution and accounts for the heterogeneity of the size of the blocks, $M_{th}$ is the threshold mass whose occurrence frequency is $\lambda$. The spread of the distribution is determined through the coefficient of variation, $CO{V}_m$, i.e., the ratio between the standard deviation and the mean value.

2.1. Failure Due to Excessive Height

This section deals with the failure due to excessive height. In the proposed approach, the state function accounting for the exceeding of the barrier height, $h_B$, is [36]:

$$H\left(\mathbf{h}\right)=H\left(\begin{array}{c}{h}_b\\ h_B\\ t\end{array}\right)=\left(h_B-t\right)-h_b,$$

(2)

where $h_b$ is the height of the trajectory in the position in which a net fence is supposed to be installed, and t is the tolerance due to the block size [31], considering that $h_b$ is measured in the centre of mass of the impacting block. The failure probability $p_{fa,h}$, related to height state function, is computed as:

$$p_{fa,h}=P\left(H\left(\mathbf{h}\right)\le 0\right)={\int \int \int }_{H\left(\mathbf{h}\right)\le 0}{f}_H\left(\mathbf{h}\right)d{h}_{b}d{h}_{B}dt,$$

(3)

where $f_H$ is the joint probability density function of the trajectory, of the barrier height, and of the tolerance. This expression is consistent with a lumped-mass assumption for the trajectory analyses, i.e., absence correlation between the size of the falling block and its trajectory. The failure probability due to excessive height during the period of analysis $\tau$, e.g., the lifetime of the system or one year, is obtained as [36]:

$$F_h\left(\tau \right)=p_{f,h}\left(\tau \right)\approx 1-exp\left(-\nu \tau p_{fa,h}\right),$$

(4)

where $\nu$ is the mean rate of occurrence of an event, which can be assumed equal to $\lambda$. As already mentioned, a Dirac-$\delta$ distribution at the value of the CE (marking) certified nominal height $h_B$ is assumed. For sake of simplicity, the same assumption holds for the tolerance t.

2.2. Failure due to Excessive Kinetic Energy

The failure due to excessive kinetic energy has to account for the variability of the impacting block size, which differs from event to event. The relationship between the rockfall block mass, to which the kinetic energy is proportional, and its temporal frequency, i.e., its return period T, is expressed by Equation (1).

Given the occurrence of an event, the probability of failure due to excessive kinetic energy $p_{fa,k}$ considers that the maximum energy depends on the occurrence of certain values of the mass and of the velocity, chosen inside their probability distributions as representative (or characteristic) as [36]:

$$p_{fa,k}={\int }_0^{\infty }{\int }_0^{\infty }{p}_{fe}|\left(m_k=\mu \mathrm{and}{v}_k=\omega \right)f_{m_k,v_k}\left(\mu ,\omega \right)d\mu d\omega ,$$

(5)

where $f_{m_k,v_k}\left(\mu ,\omega \right)$ is the joint probability density function of the characteristic values of the mass and the velocity, defined through the random variables $\mu$ and $\omega$, respectively. Under the lumped-mass assumption, the velocity and the mass are independent and, thus, the joint probability density function is the product of the single probability density functions, i.e., $f_{m_k}$ and $f_{v_k}$. Since the lumped-mass propagation analysis outputs a single distribution of the velocities, the characteristic value of the velocity is unique: $f_{v_k}$ is a Dirac-$\delta$ distribution at $v_k$. Hence, Equation (5) turns into:

$$p_{fa,k}={\int }_0^{\infty }{p}_{fe}|\left(m_k=\mu \right)f_{m_k}\left(\mu \right)d\mu .$$

(6)

Conversely, for each return period T a normal distribution centred in the mean value $m_{50}\left(T\right)$ computed through Equation (1) can be evaluated, and the probability density function $f_{m_k}\left(\mu \right)$ related to the characteristic value of the mass, variable with T, can be described through a Pareto Type I distribution according to De Biagi et al. [39] as:

$$f_{m_k}\left(\mu \right)=\left\{\begin{array}{cc}0\hfill & \mu <M_{th}\hfill \\ \frac{\alpha }{M_{th}}{\left(\frac{\mu }{M_{th}}\right)}^{-\alpha -1}\hfill & \mu \ge M_{th}\hfill \end{array}\right..$$

(7)

The conditional failure probability $p_{fe}|\left(m_k=\mu \right)$ is studied introducing the state function accounting for the exceeding of barrier energy capacity $E_B$, expressed as:

$$K\left(\mathbf{e}\right)=K\left(\begin{array}{c}{E}_B\\ m\\ v\end{array}\right)=E_B-\frac{1}{2}m{v}^2,$$

(8)

where m and v are the impacting block mass and velocity, respectively. The conditional failure probability $p_{fe}|\left(m_k=\mu \right)$ related to the state function is computed as:

$$p_{fe}|\left(m_k=\mu \right)=P\left(K\left(\mathbf{e}\right)|\left(m_k=\mu \right)\le 0\right)={\int \int \int }_{K\left(\mathbf{e}\right)\le 0}{f}_E|\left(m_k=\mu \right)\left(\mathbf{e}\right)d{E}_{B}dmdv,$$

(9)

where $f_E|\left(m_k=\mu \right)$ is the joint probability density function of the mass (with a characteristic value equal $m_k=\mu$) the velocity of the impacting block, and the capacity of the barrier. As reported in [36], the failure probability in the chosen time period $\tau$ due to excessive kinetic energy is obtained as:

$${\mathrm{F}}_k\left(\tau \right)=p_{f,k}\left(\tau \right)\approx 1-exp\left(-\nu \tau p_{fa,k}\right).$$

(10)

As in Equation (4), $\nu$ can be assumed equal to $\lambda$.

2.3. From the Reliability Design Approach to the Semi-Probabilistic Method: Equivalent Partial Safety Coefficients

Similarly for other structural engineering design problems, the reliability approach can be alternatively considered both for the design, and, for a given reliability, for finding the partial safety factors to be adopted in the current semi-probabilistic methods (code calibration procedure [37]). In the present problem, three partial safety factors, namely, ${\gamma }_h$, ${\gamma }_m$ and ${\gamma }_v$, can be defined. The safety factors are the ratios between the design values of the parameter and their characteristic values referred with the subscript ${}_k$, in accordance with EN 1990:2002 [30]. Referring to the Italian recommendation UNI 11211-4 [31], the 95th percentiles of velocity and passing height are considered as characteristic values, i.e., $v_k=v_{95}$ and $h_k=h_{95}$. To evaluate the partial safety factor ${\gamma }_m$ it is necessary to consider a reference return period $T_k$, according to which a characteristic value of the mass $m_{kT}$ is computed as the mean value $m_{kT}=m_{50}\left(T_k\right)$. Hence, it results:

$${\gamma }_h=\frac{h_{bd}}{h_k},$$

(11)

$${\gamma }_m=\frac{m_d}{m_{kT}},$$

(12)

$${\gamma }_v=\frac{v_d}{v_k}.$$

(13)

The previous can be combined into an energy partial safety factor:

$${\gamma }_E={\gamma }_m{\gamma }_v^2.$$

(14)

The design values of height, mass and velocity are obtained applying the proposed time dependent reliability approach, once defined a failure probability.

3. Methodology

Rockfall events are extremely variable, and their random nature is affected by the geometrical and mechanical characteristics of both the falling block and the slope [35,41]. Among them, the topography and the size of the block are the most important [42,43,44]. In the lumped-mass approach assumption, the variability of the mass of the block does not affect its kinematics but it directly influences its kinetic energy.

The current semi-probabilistic design method accounts only for a single characteristic value and thus, the correspondent $\gamma$s (Equations (11)–(13)) have to encompass this variability. In a previous paper, the Authors [36] highlighted the extreme variability in the shapes of the cumulative probability functions of velocity and height. Consequently, two values were adopted to characterise the trend of the distributions, namely the 95th and the 99th percentiles and a normal-tail approximation was performed. The Authors proved that the resulting $\gamma$s are affected by the ratios between these two percentiles, i.e., $h_{99}/h_{95}$ and $v_{99}/v_{95}$, rather than by each single value (95th or 99th percentile).

Furthermore, (12) and (13) are affected by multiple variables: (i) the dispersion around $m_{kT}$, i.e., $CO{V}_m$, (ii) the probability density function $f_{m_k}\left(\mu \right)$ of the possible impacting masses, and (iii) the mass-return period relationship. The parameters are highly site specific through the parameter $\alpha$. In detail, the parameter $\alpha$ refers to the heterogeneity of the sizes of the blocks masses: small values of $\alpha$ relate to a large ratio between large and small blocks (large heterogeneity). For a given design energy, the value of $m_{kT}$ is inversely proportional to ${\gamma }_E$.

The following approaches were adopted to investigate the range encompassing for the great majority of real cases of the above mentioned parameters:

• The possible ranges of $h_{99}/h_{95}$ and $v_{99}/v_{95}$ were investigated with a procedure widely adopted in the literature [42,45,46]: a synthetic slope was considered and parametric simulations performed;
• The possible values of $\alpha$ were derived from a set of surveys performed by the Authors in the Northwestern Italian Alps;
• $CO{V}_m$ was investigated through an uncertainty-based analytical approach related to the rockfall volume frequency now proposed by De Biagi et al. [40].

This activity represents a first step of a complete code calibration for a semi-probabilistic design approach [37].

3.1. Height and Velocity

A series of 2D numerical simulations with a lumped mass software were performed. RocFall software [47] was used for the analyses. A double slope profile was adopted to mimic a typical topography of a rock face, as depicted in Figure 1, with a global height H of 210 m.

The leading idea is that the topography is the parameter that mostly affects the results [42,48,49,50]. This is why, aiming at representing all the most frequent rockfall hazardous situations, the simulations were performed varying only the slope profile, keeping constant the material parameters of both slope and rock blocks, and the model settings. A linear profile, i.e., without any topographic roughness, representing the worst scenario in terms of runout, was considered. The slope angles ${\vartheta }_1$ and ${\vartheta }_2$ of the lower and the upper sections of the profile, respectively, were varied, one by one, between 45${}^{\circ }$ and 90${}^{\circ }$, with an interval of 15${}^{\circ }$, in order to consider typical rockfall prone configurations. A vertical line seeder above the upper section was used to specify a set of source locations for all the simulated rock falls. The source locations were generated with a uniform distribution along the length of the seeder, equal to 10 m. The initial velocity of the falling rocks was set as a deterministic quantity equal to 1 m/s in the vertical downward direction. The assumption of an initial velocity is assumed as precautionary, to consider possible different triggering situations, i.e., earthquake-triggered rockfall, or a complex phenomenon of sliding evolving in free fall of one or more blocks, or detachment due to the impact of blocks from upslope. Normal distributions were adopted for the normal ($R_n$) and the tangential ($R_t$) coefficients of restitution, with mean values of 0.35 and 0.85, respectively, and standard deviation of 0.04 both. On the contrary, the friction coefficient was set as a deterministic quantity equal to $tan{30}^{\circ }$. All the adopted parameters were chosen to represent a typical rock slope, according to Bar et al. [51], following the suggestion provided by Pfeiffer and Bowen [52] and Hoek and Bray [53] for bedrock or boulders with little soil or vegetation or just for bedrock in case of lumped mass models. As specifically tailored, the adopted coefficients take into account the limitation in the model kinematics evaluation, e.g., the rotational velocity of the impacting block. These assumptions entail the worst case in terms of propagation, i.e., larger runout distance and lower energy absorption during impact. To increase the accuracy of the results, the number of rocks thrown was set equal to 10,000 with a Monte Carlo sampling technique around the distributions of the restitution coefficients. A collector (C) was positioned at the slope toe, assumed as a flat deposit area, normal to ${\vartheta }_1$, as to take representative and comparable values of height and velocity.

3.2. Alpha

Several surveys were performed for searching the possible range of the parameter $\alpha$. The parameter is an index of the spread of the block sizes distribution in the location of the net fence: the more heterogeneous the deposit, the lower $\alpha$. The nature of the deposit largely depends on the lithology, the degree of fracturing of the rock mass, and the topography, i.e., slope angle and presence of obstacles on the path [54]. All these aspects influence the kinematics and the possible fragmentation of the blocks. Eight test sites located in the in the Northwestern Italian Alps and characterised by different lithology, exposure, and mean altitude of the source zone were investigated, as reported in

Table 1. Main characteristics of the test sites and correspondent $\alpha$ values.

	Lithology	Exposure	Altitude (m a.s.l.)	$\mathit{\alpha }$
Site 1 [57]	Metagabbros	SW	800	0.995
Site 2	Metabasalts	SW	1750	1.412
Site 3	Prasinite	S	620	0.746
Site 4	Limestones	S	700	0.691
Site 5	Gneiss	SW	1350	1.242
Site 6	Gneiss	W	1100	1.533
Site 7	Gneiss	W	1150	1.252
Site 8 [55]	Polimetamorphic complex	S	1530	1.019

. The surveys were performed manually by counting the size of the blocks, according to the procedure reported in [55,56].

Figure 2 shows the views of the deposit area of test sites No.3 and No.6 as representative of two opposite situations in terms of heterogeneity in block sizes of the deposit. The figure highlights a wide size variability in test site No.3, while test site No.6 is characterised by a deposit quite homogeneous. The value of $\alpha$ was determined through the maximum likelihood estimator [58]:

$$\alpha =N{\left[\sum _{i=1}^{N}ln\left(\frac{m_i}{M_{th}}\right)\right]}^{-1},$$

(15)

where N is the number of surveyed blocks and $m_i$ is the mass of each block.

3.3. Coefficient of Variation of the Mass

The coefficient of variation $CO{V}_m$ is the ratio between the standard deviation and the mean value of a distribution. In the present discussion, this value was considered representative of the spread of the distribution of the block mass at a given return period T, whose relationship follows Equation (1). The parameters of Equation (1), namely $\alpha$ and $M_{th}$ are determined through the observation of past events and by means of a detailed survey of the masses of the blocks where the rockfall net fence is expected to be installed. An uncertainty can be associated to the value of $\alpha$ depending on the number of surveyed blocks N. Following Malik [59], the real value of the parameter follows a Chi-square distribution with $\left(2N-2\right)$ degrees of freedom [60]. In other words, given the number of surveyed blocks, the maximum likelihood estimator (m.l.e.) provides a “best” value of the parameter, see Equation (15), that, in reality, is distributed as:

$$\check{\alpha }\sim \frac{\alpha }{2N}{\chi }_{2N-2}^2,$$

(16)

where $\check{\alpha }$ is the variate, i.e., the value of the shape parameter, and $\alpha$ is the estimate through the m.l.e. Asymptotic normality is proven for $\check{\alpha }$, see [60].

A distribution of masses can be associated to each return period, determined from $\check{\alpha }$ and Equation (1). The left plot in Figure 3a shows, in black, the distribution of the masses obtained with the procedure previously described for a return period equal to $T=20$ years ($\alpha =1.2$, $\lambda =0.1$ events/yr, $N=1000$, $M_{th}=1350$ kg). The red curve is the normal distribution that best fits the obtained values. A good agreement between the real distribution of the masses at a given return period and a normal distribution is proven for other return periods, as shown in Figure 3b. This similarity is the hypothesis at the base of the proposed methodology for computing the coefficient of variation of the mass $CO{V}_m$ at a given return period. Furthermore, assuming a normal distribution centred in the mean value, instead of a right-skewed Chi-squared distribution, overestimating the impacting block mass, is considered by the authors as precautionary.

To this purpose, (i) the cumulative distribution of the masses was determined, as shown in Figure 4. (ii) Two percentiles were identified, namely 15.87th and 84.13th, which represent the bounds of the mean ± one standard deviation of an equivalent normal distribution. (iii) The masses corresponding to those percentiles are identified, say $m_{0.1587}$ and $m_{0.8413}$. (iv) The standard deviation of the equivalent normal distribution was computed as:

$${\sigma }_{eq}=\frac{m_{0.8413}-m_{0.1587}}{2}.$$

(17)

(v) Supposing that the mean value of the equivalent normal distribution is the mass computed through Equation (1), $CO{V}_m$ was computed as

$$CO{V}_m=\frac{{\sigma }_{eq}}{M_{th}{\left(\lambda T\right)}^{1/\alpha }}.$$

(18)

4. Results and Discussion

In this section, the results of the performed analyses are shown and the identified range of parameters is discussed. At the end, the obtained results are adopted to evaluate a range of $\gamma$ values associated to different failure probabilities.

4.1. Height and Velocity

Figure 5a,b plot the empirical cumulative distribution of the values of the passing heights and the velocities normalised to the 95th percentile of each parameter. In analogy with what highlighted in Section 3, the results, in terms of both velocity and height, display a huge variety of probability distributions shapes, generally far from a normal one. A large sort of trends is displayed. As a limit situation, in the case in which no rebound occurs (the parabolic solving equation relates to deterministic parameters, only), the cumulative probability density functions of both height and velocity tend to a step curve at $h/h_{95}=1$ or $v/v_{95}=1$ (a Dirac-$\delta$ distribution), respectively.

Furthermore, to deeply investigate the statistics of $h_{99}/h_{95}$ and $v_{99}/v_{95}$, Figure 5a,b display the boxplot of both variables. These representations graphically reveal the variation in samples of the statistical populations neglecting their distributions, highlighting the degree of dispersion, the skewness in the data, and showing the outliers. The central red mark indicates the mean, and the bottom and the top edges indicate the 25th ($Q_1$) and 75th ($Q_3$) percentiles, respectively. The top of the upper whisker is located at $Q_3+1.5\left(Q_3-Q_1\right)$. Larger values are considered as outliers (red cross). A different trend can be observed between the considered variables: the height displays right-skewed boxplots, i.e., a positive skew, while the velocity boxplot a slightly left-skewed, i.e., a negative skew. This aspect highlights that both $h_{99}/h_{95}$ and $v_{99}/v_{95}$ cannot be related to normal distributions. Each boxplot exhibits an outlier: 2.5 for $h_{99}/h_{95}$ and 1.27 for $v_{99}/v_{95}$. Ignoring the outliers, the 75th percentile of $h_{99}/h_{95}$ and $v_{99}/v_{95}$, are 1.364 and 1.025, respectively. The upper whiskers are located at 1.870 and 1.06 for $h_{99}/h_{95}$ and $v_{99}/v_{95}$, respectively. Considering all the conservative assumption made in modelling the synthetic profiles, the Authors suggest to consider the upper quartile (75th percentile) of each set of ratios as a conservative reference value of $h_{99}/h_{95}$ and $v_{99}/v_{95}$.

4.2. Alpha

Table 1 reports the values of the parameter $\alpha$ for each test site. The results span from about 0.7 to 1.6, revealing a great site specificity of this parameter. Nevertheless, comparing the results with those found in [61,62,63], this range encloses the great majority of the cases. Figure 6 plots the complementary cumulative distributions of block volumes related to the eight study cases (left) and the best fitting power law (right).

4.3. Coefficient of Variation of the Mass

The procedure proposed in Section 3.3 was applied to define the value of $CO{V}_m$ for different T, $\alpha$, $\lambda$ and N. Besides, since the coefficient of variation is a ratio, it can be proved that it does not depend on $M_{th}$. With reference to the case reported in Figure 4, it results that $m_{0.1587}=2364$ kg and $m_{0.8413}=2453$ kg, from which it results ${\sigma }_{eq}=44.5$ kg. Substituting into Equation (1), one gets $m_{50}\left(20\mathrm{years}\right)=2405$ kg. Thus, the coefficient of variation of the mass is:

$$CO{V}_m=\frac{44.5\mathrm{kg}}{2405\mathrm{kg}}=0.0185.$$

(19)

The same procedure was repeated for various values of T, $\alpha$, $\lambda$ and N.

Table 2. Values of $CO{V}_m$ (Coefficient of variation of the mass) for various T, $\alpha$, $\lambda$ and N.

$\mathit{\alpha }$ = 0.7	N				$\mathit{\alpha }$ = 1.0	N
$\mathit{\lambda }\mathit{T}$	100	1000	10,000	100,000	$\mathit{\lambda }\mathit{T}$	100	1000	10,000	100,000
10	0.367	0.105	0.033	0.010	10	0.249	0.073	0.023	0.007
100	0.830	0.213	0.066	0.021	100	0.537	0.148	0.046	0.015
1000	1.456	0.324	0.099	0.031	1000	0.884	0.224	0.069	0.022
$\mathit{\alpha }$ = 1.25	$\mathit{N}$				$\mathit{\alpha }$ = 1.5	$\mathit{N}$
$\mathit{\lambda }\mathit{T}$	100	1000	10,000	100,000	$\mathit{\lambda }\mathit{T}$	100	1000	10,000	100,000
10	0.197	0.059	0.018	0.006	10	0.163	0.049	0.015	0.005
100	0.416	0.118	0.037	0.012	100	0.340	0.098	0.031	0.010
1000	0.667	0.178	0.055	0.017	1000	0.537	0.148	0.046	0.015

reports the values of $CO{V}_m$ corresponding to different combinations of the aforementioned parameters. It is shown that the variability of $CO{V}_m$ is large and the coefficient is affected by all the four parameters. By interpolation, an approximate expression for $CO{V}_m$ was found:

$$CO{V}_m\approx 1.3606\frac{{\left(\lambda T\right)}^{0.3}}{N^{0.525}\alpha }.$$

(20)

5. Towards a Range of Gamma Values

The values of the equivalent partial safety factors are affected by the ratios $v_{99}/v_{95}$ and $h_{99}/h_{95}$, which are related to the slope and rockfall propagation, the parameter $\alpha$, which is related to the impacting block mass distribution, and the coefficient of variation of the mass $CO{V}_m$, which directly depends on the number of surveyed blocks.

The range of such parameters was investigated in the previous section (Section 4). The upper quartiles of the obtained distributions of $h_{99}/h_{95}$ and $v_{99}/v_{95}$, i.e., 1.364 and 1.025, respectively, were considered as representative of the great majority of the propagation kinematics variability. The analysed block mass distributions lead to assume that $\alpha$ spans between 0.7 and 1.6. A variable $CO{V}_m$ computed according to Equation (20), given the number of surveyed blocks N and $\lambda$, was implemented in the reliability calculations, in particular in computing the failure probability $p_{fe}|\left(m_k=\mu \right)$.

Within the above mentioned ranges, the partial safety factors corresponding to various failure probabilities were computed accounting for ${\mathrm{F}}_h={\mathrm{F}}_k=0.5p_f$, i.e., a failure probability equally distributed between a failure related with the excessive height and a failure related to excessive kinetic energy. In the present comparative study, $p_f$ extends from $1.3\times {10}^{-6}$ per year to ${10}^{-2}$ per year. The lower probability corresponds to the annual reliability index (4.7) adopted for civil ordinary structures with a consequence class CC2, according to [30]. To calculate an annual failure probability, $\tau$ was set equal to 1 year.

The simulations were performed considering a mean frequency of events $\lambda$ in the range 0.1 and 1 events per year, while the number of surveyed blocks N spans between 200 and 1000. It is assumed that a number of surveyed blocks smaller than 200 does not constitute a statistically representative sample, while larger than 1000 is difficult to achieve. The threshold mass was set to $M_{th}=2700$ kg. Although the common definition of threshold mass, i.e., the minimum value of the fallen block mass that has always been observed and recorded in the location in which the net fence is expected to be installed [60], this quantity can also be associated to the minimum mass of the impacting block which can be relevant for the design. Two different reference block masses were adopted to compute ${\gamma }_m$, i.e., the expected mean masses at a return period of $T_k=50$ and $T_k=200$ years, respectively.

Referring to the failure due to excessive kinetic energy, Figure 7 and Figure 8 depict the value of ${\gamma }_E$ obtained through Equation (14). The values of ${\gamma }_m$ and ${\gamma }_v$ are reported in the Appendix A. It clearly emerges that the value of ${\gamma }_E$ increases as soon as a more safe condition is considered, i.e., the annual failure probability decreases. The heterogeneity of the size of the blocks that can impact the barrier largely affects the value of the coefficients. For $\alpha$ in the range 0.7 to 1.0, large values of ${\gamma }_E$ are observed. This is due to heterogeneity of the impacting blocks sizes resulting in stricter conditions in evaluating the failure probability. On the contrary, a less abrupt, but similar, trend is shown for $\alpha$ in the range 1.0 to 1.6. The trend is marginal for $T_k=200$ years. For a given failure probability, the larger the $\alpha$ the smaller the ${\gamma }_E$. The influence of the average rockfall frequency $\lambda$ is highlighted and it is similar to the effect of the number of surveyed blocks N. The increase of the safety coefficient when the database of a surveyed blocks has low cardinality (i.e., reduced number) can be attributed to the inverse proportionality between N and the coefficient of variation of the mass, $CO{V}_m$. Thus, increasing N, $CO{V}_m$ reduces. Similar trends can be observed on both ${\gamma }_m$ and ${\gamma }_v$. The effects of N and $\lambda$ are more emphasised on the values of the coefficient of the mass rather than of the velocity. As already mentioned, such trend is a consequence of the influence of N and $\lambda$ on $CO{V}_m$. This variation presupposes that a unique value of each partial safety factor cannot be identified, but the problem is strictly site dependent.

Comparing the results obtained with different return periods $T_k$ through which the characteristic mass $m_{kT}$ is defined, it clearly emerges that higher partial safety factors are obtained for $T_k=50$ years. A rough comparison between the ${\gamma }_E$ associated to $T_k=50$ years and the ones related to $T_k=200$ years highlights that the latter is more or less 5 or 6 times smaller than the former. The higher rate of decrease is observed for smaller $p_f$ and $\alpha$. The same trend emerges for ${\gamma }_v$, while ${\gamma }_m$ halves its value at $T_k=200$ years, almost independently from $\alpha$ and $p_f$.

Focusing on the value, in the case of $T_k=50$ years, the plots show values of ${\gamma }_m$, ${\gamma }_v$, and, consequently of ${\gamma }_E$, larger if compared to those proposed by the Italian and the Austrian national standards. On the contrary, for $T_k=200$ years, $p_f$ equal to 0.01, the correspondent coefficient ${\gamma }_E$ is close to one, independently from $\alpha$, N, and $\lambda$.

It is worth mentioning that the suggested values of the above mentioned national standards do not explicitly refer to a failure probability, while they account for a specific time period in which the net fence has to be considered safe, i.e., service life time equal to 25 years [19]. On the contrary, the suggested approach requires to define an annual failure probability not to be exceeded. Furthermore, the standards do not provide any suggestion related to the reference period for the characteristic mass ($T_k$), computed as a conservative value in the UNI 11211-4 [31] or as a percentile of the cumulative block size distribution in ONR 24810 [32], independently from the occurrence frequency. It reveals that increasing the return period for the characteristic mass the coefficients ${\gamma }_m$, ${\gamma }_v$, and, consequently, ${\gamma }_E$ decrease.

Referring to the failure due to excessive height, Figure 9 depicts the values of ${\gamma }_h$ obtained through Equation (11) as a function of the failure probability $p_f$. Accounting for the uncertainty on the height of the trajectory of the impacting block, the value of ${\gamma }_h$ is independent from $\alpha$, as, for the present approach, the tolerance t was assumed as a separate variable, with a Dirac-$\delta$ distribution. The obtained range spans from 1.50 to 2.73, increasing $p_f$.

From all these results, it emerges that a unique value of each partial safety factor cannot be identified. Nevertheless, as a first approximation, to partially avoid some uncertainties, the use of a block mass referred to a return period of almost 200 years can be suggested, as well as a minimum number of surveyed blocks N of about 500.

6. Conclusions

A novel time-dependent reliability approach was introduced for the design of rockfall net fences [36]. This method takes into account the variability in time of the the potential impacting mass, as well as the uncertainties related to the size and the kinematics of the impacting block. These parameters were treated as statistical quantities to which a frequency distribution was associated. This design approach was compared to the current practice based on partial safety factors applied to the characteristic values of the impacting block mass, velocity, and height. In De Biagi et al. [36], four variables related to the impacting block were identified as the most influencing parameters, i.e., the ratios $h_{99}/h_{95}$ and $v_{99}/v_{95}$, $CO{V}_m$, and the parameter $\alpha$. Different types of analyses were performed in order to evaluate the possible ranges of the above mentioned parameters representative of the great majority of situations that can be observed in a real rockfall prone area. Rockfall trajectory analyses on synthetic profiles were conducted to investigate $h_{99}/h_{95}$ and $v_{99}/v_{95}$, while onsite surveys served to evaluate the $\alpha$ parameter. The coefficient of variation of the mass $CO{V}_m$ was investigated assuming a normal distribution of the mass at a given return period. The range of $CO{V}_m$ was studied as a function of the number of surveyed blocks N in the area on which the net fence has to be installed, the frequency of the events $\lambda$, and $\alpha$, resulting in an approximate expression, Equation (20).

The obtained results were adopted to evaluate the correspondent $\gamma$ values for different failure probabilities. The Authors suggest to adopt $h_{99}/h_{95}=1.364$ and $v_{99}/v_{95}=1.025$, i.e., the 75th percentile of their distributions, while a range spanning from 0.7 to 1.6 was considered for $\alpha$. The calculations have been performed for annual failure probability $p_f$ from $1.3\times {10}^{-6}$ to ${10}^{-2}$.

It emerges that the partial safety factor ${\gamma }_h$ spans from 1.50 to 2.73 according to the failure probability, while, ${\gamma }_E$ stringently depends on $\alpha$, N, $\lambda$, and on the return period $T_k$ on which the mass $m_{kT}$ is evaluated. Decreasing $\alpha$ from 1.0 to 0.7 the rate of change of ${\gamma }_E$ increases significantly. For a given failure probability, the larger the $\alpha$, the smaller the ${\gamma }_E$. A less abrupt, but similar, trend is shown for $\alpha$ in the range 1.0 to 1.6. The trend is marginal for $T_k=200$ years. It can be observed that for $\alpha$ greater than 1.1, in case of $T_k=200$ years, the influence of the failure probability is quite negligible. The influence of $\lambda$ is present and it is similar to the effect of the number of surveyed blocks N: the lower the $\lambda$ and the higher the N, the lower the ${\gamma }_E$. Comparing the results obtained for different return periods $T_k$, higher partial safety factors are obtained for $T_k=50$ years: the ${\gamma }_E$ associated to $T_k=50$ years are more or less 5 or 6 times higher than the ones related to $T_k=200$ years. At $T_k=50$ years, the plot shows values of ${\gamma }_m$, ${\gamma }_v$, and, consequently of ${\gamma }_E$, very high compared to those proposed by the Italian and the Austrian national standards, which do not account for an annual failure probability or for another target value.

The obtained results highlights that the design value for the height of the block can be several times the characteristic value. The same can be observed for the energy in case of short reference return period $T_k$, and for low value of $\alpha$. This is justified with the fact that the value of the mass changes with the return period.In addition, the number of surveyed blocks N affects the uncertainties related to the characteristic mass and, consequently, the $\gamma$ values.

All these considerations result in the fact that, for a given failure probability, quantifying a set of partial safety factors valid for different sites is not possible since they are strictly site dependent. On the contrary, if a set of partial safety factors is assigned to different sites, an intrinsic variability in the failure probability is accepted. As a first suggestion, to partially avoid some uncertainties, the use of a block mass referred to a return period of almost 200 years can be suggested, as well as a surveyed block number N of about 500. Referring to the ratios between velocities and heights, namely $v_{99}/v_{95}$ and $h_{99}/h_{95}$, the Authors suggest to adopt the values herein proposed, which are representative of the 75% of the simulated synthetic profiles. The provided results should be accurately used, accordingly with the required degree of safety. Further developments should be done in order to consider the influence of the mass in the evaluation of the design height, i.e., introducing the block size in the applied tolerance value. Additional scenarios for the trajectory analyses related to different slope profile or input parameters, e.g., restitution coefficients, would be analysed in future researches.