# Theoretical Analysis of Rockfall Impacts on the Soil Cushion Layer of Protective Structures

Authors: Xing Wang; Yongxu Xia; Tianyue Zhou
Journal: Advances in Civil Engineering
 (2018)

---

## Abstract

During the Wenchuan Earthquake, with a magnitude of 5.12, collapses and rockfall hazards persisted for a long time after the initial investigations carried out by research fellow S. M. He and his team at the scene of the disaster in October 2008. It is possible that additional incidents of rockfalls in large quantities may continue in the same areas over the next ten to fifteen years. Furthermore, in the vast mountainous region of western China, the topographic relief is evident, and earthquakes occur frequently. Therefore, it is difficult to effectively defend against rockfall hazards. When designing protective structures, the key issue is the analysis of the mechanical response mechanism of the soil cushion layer of the upper cushion when subjected to the impact of rockfall. As such, a theoretical method was used to perform such an analysis. The cavity expansion and energy conservation model were adopted. Analytical solutions for the impact force and penetration depth were then derived. Furthermore, the impact force and penetration depth of rockfall were studied with the LS-DYNA software to obtain values for the impact forces and the penetration depth. Finally, the reliability of the theoretical method was evaluated using the cavity expansion, energy conservation, numerical simulation, Hertz, Japanese, Swiss, Australian, B. S. Guan, tunnel manual, and subgrade methods based on an engineering model. The results show that the cavity expansion and the energy conservation methods yielded consistent results. Meanwhile, the cavity expansion and the energy conservation methods also yielded consistent results with the numerical simulation, Japanese (obtained by laboratory experiment), Swiss (obtained by laboratory experiment), and Australian (obtained by field experiment) methods. The relevant methods and conclusions shall therefore be applied to the design of rockfall protection structure in future investigations.


---

## Body

## 1. Introduction

As a key aspect in the engineering of highway routes, the technology for building tunnels has gradually improved over time [1–3]. Shed tunnel protection structures are special tunnel structures. These structures and others such as open cut tunnels and passive protective structures in hazardous areas are designed primarily to protect against the impact forces of rockfall. In order to avoid immediate impact load on protected structures, a soil layer of a certain thickness is generally used [4]. During the rockfall impact on the soil cushion layer, the time interval is short and the energy conversion is complicated. In addition, elastoplastic deformation of the soil cushion layer occurs. Therefore, the impact process of the rockfall is a complex impact dynamics problem. In the design of protective structures, the key issue is the analysis of the mechanical response mechanism of the soil cushion layer under the impact of rockfall.Many scholars and research institutes are directing increasing attention to rockfall disaster research. The classic Hertz contact theory and Thornton elastoplastic impact theory have been applied widely [5–8]. Johnson described in his book “Contact Mechanics” the details of the mechanical properties of an object under the effects of static contact, sliding, rolling, and impact [9]. Labiouse et al. presented a semiempirical and semitheoretical method to analyze rockfall impact forces based on impact experiments (hereinafter referred to as the Swiss method) [10]. Kishi and Ikeda studied the mechanical response characteristics of a retaining wall under rockfall impact based on field experiments and numerical simulations [11]. Heidenreich and Labiouse studied the mechanical characteristics of a slope with a loose structure under the impact of rockfall using small-scale model experiments [12]. Calvetti and Prisco analyzed the mechanical response characteristics of gravel soil under rockfall impact based on experiments and numerical simulations [13]. Based on the field test method and the theorem of impulse, Pichler and Hillmich developed a calculation method for impact force and penetration depth associated with rockfall (hereinafter referred to as the Australian method) [7]. Mougin et al. examined the mechanical response of rockfall on a concrete plate under direct impact conditions based on model experiments [14]. Delhomme et al. adopted model experiments and numerical simulations to analyze the energy dissipation and seismic mitigation effects of the new type of shed tunnel structures [15]. The Japanese Road Association developed a semiempirical and semitheoretical calculation method for rockfall impact force using model experiments (hereinafter referred to as the Japanese method) [16]. Calvetti and Prisco analyzed the design method for shed tunnel structures based on the uncoupled calculation method [17]. Ye and Chen performed a comparative study on rockfall impact force values by employing the Japanese, Swiss, and Hertz methods in conjunction with an engineering model [18].Ye et al. analyzed the rockfall impact forces under oblique impact on the soil cushion layer [19]. Qi et al. set up an analysis model for rockfall impact forces based on a theoretical approach and studied the amplification coefficient of rockfall impact forces [20]. Yuan et al. analyzed the influence of the weight, shape, falling height, and impact angle of rockfall on the impact force based on experiments [21]. Sun et al. attempted to add discarded tires to the soil cushion layer and investigated the impact resistance of the new cushion structure via an experimental approach [22]. Zhang et al. investigated the mechanical response mechanism of oil/gas pipe under rockfall impact using a numerical simulation method [23]. Hu et al. discovered the law of jumping range of rockfall based on three influential factors: shape, weight, and the release altitude of rockfall using an experimental method [24]. Lam et al. studied the overturning stability of L-shaped rigid barriers subjected to rockfall impacts [25]. Basharat et al. analyzed the effects of volume and topographic parameters on rockfall travel distance based on a case study from the NW Himalayas, Pakistan [26].Zhu et al. analyzed particle size and the thickness of the gravel cushion layer on the restitution coefficient of rockfall in rockfall impact procession [27]. Luo et al. experimentally studied the rockfall impact force of the frame shed tunnel structure [28]. Asteriou and Tsiambaos analyzed the effect of rockfall impact velocity, block mass, and hardness on rockfall restitution coefficients [29]. Toe et al. adopted a metamodel hazard mitigation strategy to improve the effectiveness of rockfall protection barriers [30]. Wang et al. adopted an experimental approach to analyzing the mechanical response law of an EPS cushion layer under rockfall impact [31]. Bhatti combined experimental and numerical simulation methods to study the energy dissipation property of the EPS cushion layer under the impact of rockfall [32]. Bi et al. adopted a numerical simulation method to analyze the mechanical response mechanism of a bumper plate under the rockfall impact [33]. Yu et al. studied the maximum impact force of rockfall using an experimental approach [34]. Zhang et al. devised an estimation model for calculating the maximum impact force of rockfall based on contact theory and verified this model experimentally [35]. Yan et al. utilized numerical simulations to analyze the mechanical response mechanism of reinforced concrete slabs under the impact of rockfall [36].Prof. B. S. Guan proposed a calculation method for impact forces and penetration depth based on experiments. Further, he considered the effect on the thickness of the soil cushion layer (hereinafter referred to as the B. S. Guan method) [37]. Using the theory of momentum, the impact force and penetration depth are given in the Technical Manual for Railway Engineering Design Tunnel (Revised Edition) (hereinafter referred to as the tunnel manual method). In addition, in Code for Design of Highway Subgrade of China, the rockfall impact force and penetration depth are given based on the work-energy principle (hereinafter referred to as the subgrade method).The Japanese, Swiss (based on indoor experiments), and Australian methods (based on field tests) are three classical methods for the design of protection structures. However, these three methods are all semitheoretical and semiempirical. Moreover, few scholars have focused on theoretical derivations to obtain the penetration depth and impact force. Meanwhile, some researchers have indicated that the value of the impact force is small in China [18].Based on the dynamic cavity expansion model and the energy conservation model, a theoretical mechanical system of rockfall impacting the soil cushion layer is established and analytic solutions of rockfall impact force and penetration depth are obtained. In addition, the LS-DYNA software is used to analyze the change laws of impact force and penetration depth during rockfall impact processes. Finally, this study will compare and analyze the calculation results for penetration depth and impact forces calculated based on the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, Australian, Hertz, B. S. Guan, tunnel manual, and subgrade methods.


## 2. Conventional Classical Calculation Theories

### 2.1. Hertz Method

In 1985, Hertz derived expressions for the maximum deformation and maximum impact based on elastic theory.(1)Impact force:F=32327ER453Mv02128R1/2E3/5,Penetration depth:L=453Mv02128R1/2E2/5,where F is the rockfall impact force, E is the elasticity modulus of soil cushion layer, R is the radius of the rockfall, M is the mass of the rockfall, v0 is the initial velocity of rockfall, and L is the penetration depth of the rockfall impact.


### 2.2. Swiss Method

In 1996, a Swiss scholar, Labiouse, developed the relevant empirical method for rockfall impact forces based on rockfall experiments.(2)Impact force:F=1.765⋅ME2/5⋅R1/5⋅QH3/5,where ME is the deformation modulus of soil cushion layer, Q is the quality of rockfall, and H is the falling height of the rockfall.


### 2.3. Japanese Method

In 2000, based on rockfall experiments and combined with Hertz elastic theory, the Japanese Road Association presented the relevant semiempirical and semitheoretical calculation method of rockfall impact forces.(3)Impactforce:F=2.108⋅Mg2/3⋅λ2/5⋅H3/5,lameconstant:λ=υE1+υ1−2υ,where g is the gravitational acceleration, λ is the lame constants, and υ is the Poisson’s ratio of soil cushion layer.


### 2.4. Australian Method

In 2005, Australian scholar Pichler proposed the semiempirical and semitheoretical calculation method of rockfall impact force and penetration depth combined with the rockfall field test.(4)Impactforce:F=2Mv0tw=Mv02L,Impacttime:tw=2Lv0,Penetrationdepth:L=D103500HRc+19180H,LD≤1.257,L=D1.518ln1+19182.39H/Rc1.414+1.257,LD≥1.257,where tw is the impact time of the rockfall, D is the diameter of the rockfall, and Rc is the indentation resistance of target materials.


### 2.5. B. S. Guan Method

In 1996, the scholar B. S. Guan established the empirical method of impact force and penetration depth based on laboratory experiments. This method considers the influence of the thickness of the soil cushion layer on impact force.(5)Impactforce:Fmax=ζMa,Acceleration:a=2gHtw,Impacttime:tw=11000.097Mg+2.21h+0.045H+1.2,Penetrationdepth:L=κQH10,where ζ is the correction coefficient of the rockfall impact force, κ is the coefficient of rockfall penetration depth, a is the acceleration of the rockfall impact, and h is the thickness of soil cushion layer.


### 2.6. Method of Tunnel Manual

The method of tunnel manual refers to an approach recommended inTechnical Manual for Railway Engineering Design Tunnel (Revised Edition), and in essence, it is also an approximation method of the theorem of momentum.(6)Impactforce:F=Qv0gtw,Impacttime:tw=2hc,Compressionwavevelocity:c=1−υ1+υ1−2υ⋅Eρ0,Penetrationdepth:L=v0Q2gγψ×12tg445°+φ∗/2−1,where c is the reciprocating velocity of compression waves in the soil cushion layer, γ is the specific gravity of the soil cushion layer, ψ is the cross-sectional area of an equivalent sphere for rockfall, φ∗ is the inner friction angle of soil cushion layer, and ρ0 is the density before deformation of spherical cavity microbody.


### 2.7. Subgrade Method

Code for Design of Highway Subgrade provides a computational formula for rockfall impact force and penetration depth based on the work-energy principle.(7)Impactforce:F=2γL2tg445°+φ∗/2−1ψ,Penetrationdepth:L=v0Q2gγψ×12tg445°+φ∗/2−1.


## 2.1. Hertz Method

In 1985, Hertz derived expressions for the maximum deformation and maximum impact based on elastic theory.(1)Impact force:F=32327ER453Mv02128R1/2E3/5,Penetration depth:L=453Mv02128R1/2E2/5,where F is the rockfall impact force, E is the elasticity modulus of soil cushion layer, R is the radius of the rockfall, M is the mass of the rockfall, v0 is the initial velocity of rockfall, and L is the penetration depth of the rockfall impact.


## 2.2. Swiss Method

In 1996, a Swiss scholar, Labiouse, developed the relevant empirical method for rockfall impact forces based on rockfall experiments.(2)Impact force:F=1.765⋅ME2/5⋅R1/5⋅QH3/5,where ME is the deformation modulus of soil cushion layer, Q is the quality of rockfall, and H is the falling height of the rockfall.


## 2.3. Japanese Method

In 2000, based on rockfall experiments and combined with Hertz elastic theory, the Japanese Road Association presented the relevant semiempirical and semitheoretical calculation method of rockfall impact forces.(3)Impactforce:F=2.108⋅Mg2/3⋅λ2/5⋅H3/5,lameconstant:λ=υE1+υ1−2υ,where g is the gravitational acceleration, λ is the lame constants, and υ is the Poisson’s ratio of soil cushion layer.


## 2.4. Australian Method

In 2005, Australian scholar Pichler proposed the semiempirical and semitheoretical calculation method of rockfall impact force and penetration depth combined with the rockfall field test.(4)Impactforce:F=2Mv0tw=Mv02L,Impacttime:tw=2Lv0,Penetrationdepth:L=D103500HRc+19180H,LD≤1.257,L=D1.518ln1+19182.39H/Rc1.414+1.257,LD≥1.257,where tw is the impact time of the rockfall, D is the diameter of the rockfall, and Rc is the indentation resistance of target materials.


## 2.5. B. S. Guan Method

In 1996, the scholar B. S. Guan established the empirical method of impact force and penetration depth based on laboratory experiments. This method considers the influence of the thickness of the soil cushion layer on impact force.(5)Impactforce:Fmax=ζMa,Acceleration:a=2gHtw,Impacttime:tw=11000.097Mg+2.21h+0.045H+1.2,Penetrationdepth:L=κQH10,where ζ is the correction coefficient of the rockfall impact force, κ is the coefficient of rockfall penetration depth, a is the acceleration of the rockfall impact, and h is the thickness of soil cushion layer.


## 2.6. Method of Tunnel Manual

The method of tunnel manual refers to an approach recommended inTechnical Manual for Railway Engineering Design Tunnel (Revised Edition), and in essence, it is also an approximation method of the theorem of momentum.(6)Impactforce:F=Qv0gtw,Impacttime:tw=2hc,Compressionwavevelocity:c=1−υ1+υ1−2υ⋅Eρ0,Penetrationdepth:L=v0Q2gγψ×12tg445°+φ∗/2−1,where c is the reciprocating velocity of compression waves in the soil cushion layer, γ is the specific gravity of the soil cushion layer, ψ is the cross-sectional area of an equivalent sphere for rockfall, φ∗ is the inner friction angle of soil cushion layer, and ρ0 is the density before deformation of spherical cavity microbody.


## 2.7. Subgrade Method

Code for Design of Highway Subgrade provides a computational formula for rockfall impact force and penetration depth based on the work-energy principle.(7)Impactforce:F=2γL2tg445°+φ∗/2−1ψ,Penetrationdepth:L=v0Q2gγψ×12tg445°+φ∗/2−1.


## 3. Theoretical Analysis Based on Dynamic Cavity Expansion Model

The spherical expansion model is widely applied to penetration mechanical model analysis. In this section, the impact of rockfall on a soil cushion layer is simplified by considering the normal impact of a rigid ball on a semi-infinite soil cushion layer.By referring to Figure1, the normal stress of point A on the soil cushion layer σA could be calculated using σrqA of the cavity expansion sphere at the same position A, so a series of spherical cavities originating from the bottom of the sphere could be used to simulate the process of crater formation due to rockfall impact on the soil cushion layer. The radius of these cavity spheres gradually increases from 0 to rq, thus forming the so-called dynamic spherical cavity expansion model.Figure 1
Force and displacement of spherical cavity microelement.The radius of cavity expansion sphererqA is(8)rqA=xAsinα,where rqA is the radius of the spherical cavity, xA is the x coordinate value of point A, and α is the acute angle between the normal line and the vertical line at point A.During the normal impact of rockfall on the soil cushion layer at a velocity ofvh, the velocity of the dynamic spherical cavity expansion is(9)vq=r˙q=vhcosα,where vq is the velocity of spherical cavity expansion, r˙q is the expansion velocity of spherical cavity, and vh is the velocity of rockfall impact.And the acceleration of the dynamic cavity expansion is(10)v˙q=r¨q=v˙hcosα,where v˙q is the acceleration of spherical cavity expansion, r¨q is the expansion acceleration of spherical cavity, and v˙h is the acceleration of a rockfall impact.


### 3.1. Mechanical Analysis of Microunit with Spherical Coordinates

As shown in Figure2, the Euler coordinate of spherical coordinates is employed to perform mechanical analysis on the selected microbody.Figure 2
Spherical cavity.By referring to Figure2, a mechanical analysis of the microbody selected yields equation (7) according to the theorem of momentum:(11)σr+dσrr+drdθ2−σrrdθ2−4σθdrrdθsindθ2dt=−ρ0dvrdθ2dr,where σr is the radial stress of the spherical cavity microbody, σθ is the hoop stress of the spherical cavity microbody, r is the distance from the coordinate origin to the mass point of the spherical cavity microbody, v is the moving velocity of the mass point of the spherical cavity microbody, and θ is the opening angle of spherical cavity microbody.The velocity of the mass point of the microbody in this coordinate is(12)v=∂s∂t+v∂s∂r,where s is the displacement of the mass point of the spherical cavity microbody and t is the moving time of the mass point of the spherical cavity microbody.And the acceleration of the mass point could be expressed by formula (8):(13)dvdt=∂v∂t+v∂v∂r.For the smallerdθ, sindθ/2 could be expressed as approximately dθ/2. Substituting this simplified relation into formula (11) and dividing both sides of the simplified equation by r2dr dθ2, if the high order trace is omitted, considering equations (11) and (13), we obtain(14)∂σr∂r+2σr−σθr=−ρ0∂v∂t+v∂v∂r.Figure3 represents the force and displacement of a spherical cavity microelement.Figure 3
Force and displacement of spherical cavity microelement.Combined with the law of conservation of mass, we obtain(15)rdθ2dr⋅ρ0=ρr−sdθ2dr−s,where ρ is the density after deformation of spherical cavity microbody.After the partial derivative of equation (15), we obtain(16)13∂∂rr−s3=ρ0ρr2.Taking the integral of equation (16) and assuming the medium of the plastic zone is incompressible, we obtain(17)r−s3=r3+ft,rq≤r.Further, combined with the boundary conditions ofr=rq, s=rq, we obtain the displacement of the microbody, as follows:(18)s=r1−1−rq3r31/3,where rq is the radius of spherical cavity.Combined with equations (18) and (12), we obtain an expression for the velocity of the microbody, as follows:(19)v=rq2r2drqdt.According to the continuous medium assumption, the expression of for the velocity shall be applied to the medium of both elastic and plastic zones.


### 3.2. Mechanical Analysis of the Elastic Zone

After rockfall impact on the soil cushion layer, the elastic zone formed within this layer has the following relationships:(1)
Physical equation(20)εr=σr−2υσθE,εθ=1−υσθ−υσrE,(2)
Geometric equation(21)εr=−∂s∂r,εθ=−sr,where εr is the radial strain of spherical cavity microbody and εθ is the hoop strain of spherical cavity microbody.Substituting the expression for the displacement of the microbodys into equation (21) and using a Taylor series to expand, only the first two terms of the expansion are substituted into equation (20) on the premise of ensuring the effectiveness of the calculation results. We then obtain(22)σr−σθ=E1+υrqr3.Combining equations (22), (19), and (14), we obtain(23)∂σr∂r+2Erq31+υr4=−ρ02rqr˙q+rq2r¨qr2−2rq4r˙q2r5.Considering the boundary conditionsr ⟶ ∞ and σr = 0 and taking the integral of equation (23), we obtain(24)σr=2Erq231+υr3+2ρ0rqr˙q2+ρ0rq2r¨qr−ρ0r˙q2rq42r4.


### 3.3. Mechanical Analysis of Plastic Zone

For the ideal elastic-plastic soil cushion layer in the case of shear failure of the soil cushion layer, combined with Thornton’s theory, equation (25) is established:(25)σr−σθ=2τ0.Combined with the dynamic strength theory proposed by Richart,τ0 can be determined using equation (26):(26)τ0=Kτ,where τ0 is the shear strength of soil cushion layer, τ is the undrained shear strength of soil cushion layer, and K is the shear strength coefficient of soil cushion layer.Considering equations (22), (19), and (14), we obtain(27)∂σr∂r+4τ0r=−ρ02rqr˙q+rq2r¨qr2−2rq4r˙q2r5.Taking the boundary conditionsr=rq and σr=σrq into account, we obtain(28)σr=σrq+4τ0lnrqr−3ρ0r˙q22+ρ0rqr¨q−2ρ0rqr˙q2+ρ0rq2r¨qr+ρ0rq4r˙q22r4,where σrq is the radial stress on the spherical cavity surface.Using the boundary conditions ofr=rp and σrq−=σrq+ according to the assumption of the soil cushion layer being in the elastic-plastic zone and substituting the boundary conditions into equation (24) and (28), respectively, we can obtain an expression for the normal stress on the surface of cavity σrq.(29)σrq=2E31+υrqrp3−4τ0lnrqrp+3ρ0r˙q22+ρ0rqr¨q,where rp is the plastic region radius of spherical cavity.Combining equations (22) and (25), we obtain(30)rqrp=2τ01+υE3.Substituting equation (30) into (29), we obtain(31)σrq=4τ03lneE2τ01+υ+3ρ0vq22+rqρ0dvqdt.The process of rockfall impact on cratering can be imitated through a series of infinitesimal cavity expansion originating from the bottom of the sphere. As shown in Figure4, for the horizontal infinitesimal cavity that lies at point G, the cavity radius gradually increases from 0 to rqG, and then forms the final infinitesimal cavity QG. The formation process for the infinitesimal cavity QN is also the same. Further, we perform an integration to analyze the entire process for the impact of rockfall on Earth mass.Figure 4
Dynamic cavity of analysis of rockfall.From Figure4, the velocity and acceleration of the dynamic cavity expansion are(32)vq=r˙q=vhcosφ,v˙q=r¨q=v˙hcosφ,where φ is the acute angle between the tangent at the spherical cavity surface through point G with x axis.The value of cosφ is determined using relevant geometrical relations. A tangent line of the cavity is made through point G to obtain the geometrical relation expression as follows:(33)y=R−R2−x2,where x is the x coordinate value of point G on spherical cavity and y is the y coordinate value of point G on spherical cavity.SinceL1 is the tangent line through point G, we obtain(34)y1′=tanφ,where y1′ is the slope of the tangent line L1.According to the trigonometric function relation,(35)sec2φ−1=tan2φ,secφ=cos−1φ,and we can combine equations (33)–(35) to obtain(36)cosφ=11+y1′2.


### 3.4. Impact Force

Equation (31) is the expression of the normal stress on the surface of cavity σrq, over which the integral is performed so as to obtain the impact force of the rockfall.(37)F=2π∫0R1xσrqdxdl⋅dl=2π∫0R1xσrqdx,where R1 is the radius of the horizontal infinitesimal cavity at the upper surface of soil cushion layer, and l is the width value of infinitesimal cavity.R
1 is calculated according to the equation below:(38)R1=R2−R−L2=2RL−L2.Combined equations (31), (32), (37), and (38), and given that when rq = R and Lmax ≤ R, the impact force of rockfall on the soil cushion layer shall be calculated according to equation (39):(39)F1=4πτ032RL−L2lneE2τ01+υ+3πρ0R4−R−L4vh24R2+2πρ0R3−R−L33dvhdt,where F1 is the rockfall impact force.Generally, the velocity of rockfall for the calculation is 4–24 m/s. For the convenience of calculation, we omit3ρ0vq2/2 and rqρ0dvq/dt, in equation (31), so equation (31) can be rewritten as(40)σrq=4τ03lneE2τ01+υ,and equation (39) for the calculation of the impact force can be simplified into(41)F2=π2RL−L24τ03lneE2τ01+υ,where F2 is the rockfall impact force.During the rockfall impact process of the rolling stone, whenLmax ≥ R, the value of R can be substituted into that of L in equation (37) and (41).


### 3.5. Penetration Depth

Combined with Newton’s second law, we obtain the relation between the impact of velocity of rockfall and penetration depth:(42)dvhdt=vhdvhdL.During the impact process, if the velocity and acceleration of the rockfall are not considered for the influence on penetration depth, and whenLmax ≤ R, we multiply the rockfall mass M on both ends of equation (42) and obtain the equation for the impact force as follows:(43)F2=−MvhdvhdL,where Lmax is the maximum penetration depth of the rockfall by the cavity method.That is,(44)4πτ032RL−L2lneE2τ01+υ=−MvhdvhdL.Combined with the initial conditions ofvh=v0 and L = 0 and taking the integral of equation (44), we obtain(45)4πτ03MRL2−L33lneE2τ01+υ=v02−vh22.From equation (45), we can obtain the penetration depth of rockfall at any time when the velocity is vh, during the period in which impact velocity decreases from v0 to 0. When vh = 0, we obtain the maximum penetration depth of the rockfall Lmax.(46)4πτ03MRLmax2−Lmax33lneE2τ01+υ=v022.If the calculation result satisfiesLmax ≥ R, in equation (44), we set L = R.If the influence of rockfall velocity and acceleration during the impact process is considered and whenLmax ≤ R, we obtain(47)F1=−MvhdvhdL.Substituting equations (39) into (47), we obtain(48)dvhdL=−4τ0/3lneE/τ01+υ2πRL−πL2+3πρ0R4−R−L4/4R2vh2M−2πρ0/3R3−R−L3vh.Considering the initial conditions ofvh=v0 and L=0, when the calculation result for the penetration depth does not satisfy Lmax ≥ R, in equation (45), we set L = R. Using an analytical method, we can obtain the relation between L and vh.Using equation (48), we obtain(49)F=M4τ0/3lneE/2τ01+υ2πRL−πL2+3πρ0R4−R−L4/4R2vh2M−2πρ0/3R3−R−L3.We can substitute the values ofL and vh obtained using the analytical method into equation (49), to calculate the value of the impact force F.


## 3.1. Mechanical Analysis of Microunit with Spherical Coordinates

As shown in Figure2, the Euler coordinate of spherical coordinates is employed to perform mechanical analysis on the selected microbody.Figure 2
Spherical cavity.By referring to Figure2, a mechanical analysis of the microbody selected yields equation (7) according to the theorem of momentum:(11)σr+dσrr+drdθ2−σrrdθ2−4σθdrrdθsindθ2dt=−ρ0dvrdθ2dr,where σr is the radial stress of the spherical cavity microbody, σθ is the hoop stress of the spherical cavity microbody, r is the distance from the coordinate origin to the mass point of the spherical cavity microbody, v is the moving velocity of the mass point of the spherical cavity microbody, and θ is the opening angle of spherical cavity microbody.The velocity of the mass point of the microbody in this coordinate is(12)v=∂s∂t+v∂s∂r,where s is the displacement of the mass point of the spherical cavity microbody and t is the moving time of the mass point of the spherical cavity microbody.And the acceleration of the mass point could be expressed by formula (8):(13)dvdt=∂v∂t+v∂v∂r.For the smallerdθ, sindθ/2 could be expressed as approximately dθ/2. Substituting this simplified relation into formula (11) and dividing both sides of the simplified equation by r2dr dθ2, if the high order trace is omitted, considering equations (11) and (13), we obtain(14)∂σr∂r+2σr−σθr=−ρ0∂v∂t+v∂v∂r.Figure3 represents the force and displacement of a spherical cavity microelement.Figure 3
Force and displacement of spherical cavity microelement.Combined with the law of conservation of mass, we obtain(15)rdθ2dr⋅ρ0=ρr−sdθ2dr−s,where ρ is the density after deformation of spherical cavity microbody.After the partial derivative of equation (15), we obtain(16)13∂∂rr−s3=ρ0ρr2.Taking the integral of equation (16) and assuming the medium of the plastic zone is incompressible, we obtain(17)r−s3=r3+ft,rq≤r.Further, combined with the boundary conditions ofr=rq, s=rq, we obtain the displacement of the microbody, as follows:(18)s=r1−1−rq3r31/3,where rq is the radius of spherical cavity.Combined with equations (18) and (12), we obtain an expression for the velocity of the microbody, as follows:(19)v=rq2r2drqdt.According to the continuous medium assumption, the expression of for the velocity shall be applied to the medium of both elastic and plastic zones.


## 3.2. Mechanical Analysis of the Elastic Zone

After rockfall impact on the soil cushion layer, the elastic zone formed within this layer has the following relationships:(1)
Physical equation(20)εr=σr−2υσθE,εθ=1−υσθ−υσrE,(2)
Geometric equation(21)εr=−∂s∂r,εθ=−sr,where εr is the radial strain of spherical cavity microbody and εθ is the hoop strain of spherical cavity microbody.Substituting the expression for the displacement of the microbodys into equation (21) and using a Taylor series to expand, only the first two terms of the expansion are substituted into equation (20) on the premise of ensuring the effectiveness of the calculation results. We then obtain(22)σr−σθ=E1+υrqr3.Combining equations (22), (19), and (14), we obtain(23)∂σr∂r+2Erq31+υr4=−ρ02rqr˙q+rq2r¨qr2−2rq4r˙q2r5.Considering the boundary conditionsr ⟶ ∞ and σr = 0 and taking the integral of equation (23), we obtain(24)σr=2Erq231+υr3+2ρ0rqr˙q2+ρ0rq2r¨qr−ρ0r˙q2rq42r4.


## 3.3. Mechanical Analysis of Plastic Zone

For the ideal elastic-plastic soil cushion layer in the case of shear failure of the soil cushion layer, combined with Thornton’s theory, equation (25) is established:(25)σr−σθ=2τ0.Combined with the dynamic strength theory proposed by Richart,τ0 can be determined using equation (26):(26)τ0=Kτ,where τ0 is the shear strength of soil cushion layer, τ is the undrained shear strength of soil cushion layer, and K is the shear strength coefficient of soil cushion layer.Considering equations (22), (19), and (14), we obtain(27)∂σr∂r+4τ0r=−ρ02rqr˙q+rq2r¨qr2−2rq4r˙q2r5.Taking the boundary conditionsr=rq and σr=σrq into account, we obtain(28)σr=σrq+4τ0lnrqr−3ρ0r˙q22+ρ0rqr¨q−2ρ0rqr˙q2+ρ0rq2r¨qr+ρ0rq4r˙q22r4,where σrq is the radial stress on the spherical cavity surface.Using the boundary conditions ofr=rp and σrq−=σrq+ according to the assumption of the soil cushion layer being in the elastic-plastic zone and substituting the boundary conditions into equation (24) and (28), respectively, we can obtain an expression for the normal stress on the surface of cavity σrq.(29)σrq=2E31+υrqrp3−4τ0lnrqrp+3ρ0r˙q22+ρ0rqr¨q,where rp is the plastic region radius of spherical cavity.Combining equations (22) and (25), we obtain(30)rqrp=2τ01+υE3.Substituting equation (30) into (29), we obtain(31)σrq=4τ03lneE2τ01+υ+3ρ0vq22+rqρ0dvqdt.The process of rockfall impact on cratering can be imitated through a series of infinitesimal cavity expansion originating from the bottom of the sphere. As shown in Figure4, for the horizontal infinitesimal cavity that lies at point G, the cavity radius gradually increases from 0 to rqG, and then forms the final infinitesimal cavity QG. The formation process for the infinitesimal cavity QN is also the same. Further, we perform an integration to analyze the entire process for the impact of rockfall on Earth mass.Figure 4
Dynamic cavity of analysis of rockfall.From Figure4, the velocity and acceleration of the dynamic cavity expansion are(32)vq=r˙q=vhcosφ,v˙q=r¨q=v˙hcosφ,where φ is the acute angle between the tangent at the spherical cavity surface through point G with x axis.The value of cosφ is determined using relevant geometrical relations. A tangent line of the cavity is made through point G to obtain the geometrical relation expression as follows:(33)y=R−R2−x2,where x is the x coordinate value of point G on spherical cavity and y is the y coordinate value of point G on spherical cavity.SinceL1 is the tangent line through point G, we obtain(34)y1′=tanφ,where y1′ is the slope of the tangent line L1.According to the trigonometric function relation,(35)sec2φ−1=tan2φ,secφ=cos−1φ,and we can combine equations (33)–(35) to obtain(36)cosφ=11+y1′2.


## 3.4. Impact Force

Equation (31) is the expression of the normal stress on the surface of cavity σrq, over which the integral is performed so as to obtain the impact force of the rockfall.(37)F=2π∫0R1xσrqdxdl⋅dl=2π∫0R1xσrqdx,where R1 is the radius of the horizontal infinitesimal cavity at the upper surface of soil cushion layer, and l is the width value of infinitesimal cavity.R
1 is calculated according to the equation below:(38)R1=R2−R−L2=2RL−L2.Combined equations (31), (32), (37), and (38), and given that when rq = R and Lmax ≤ R, the impact force of rockfall on the soil cushion layer shall be calculated according to equation (39):(39)F1=4πτ032RL−L2lneE2τ01+υ+3πρ0R4−R−L4vh24R2+2πρ0R3−R−L33dvhdt,where F1 is the rockfall impact force.Generally, the velocity of rockfall for the calculation is 4–24 m/s. For the convenience of calculation, we omit3ρ0vq2/2 and rqρ0dvq/dt, in equation (31), so equation (31) can be rewritten as(40)σrq=4τ03lneE2τ01+υ,and equation (39) for the calculation of the impact force can be simplified into(41)F2=π2RL−L24τ03lneE2τ01+υ,where F2 is the rockfall impact force.During the rockfall impact process of the rolling stone, whenLmax ≥ R, the value of R can be substituted into that of L in equation (37) and (41).


## 3.5. Penetration Depth

Combined with Newton’s second law, we obtain the relation between the impact of velocity of rockfall and penetration depth:(42)dvhdt=vhdvhdL.During the impact process, if the velocity and acceleration of the rockfall are not considered for the influence on penetration depth, and whenLmax ≤ R, we multiply the rockfall mass M on both ends of equation (42) and obtain the equation for the impact force as follows:(43)F2=−MvhdvhdL,where Lmax is the maximum penetration depth of the rockfall by the cavity method.That is,(44)4πτ032RL−L2lneE2τ01+υ=−MvhdvhdL.Combined with the initial conditions ofvh=v0 and L = 0 and taking the integral of equation (44), we obtain(45)4πτ03MRL2−L33lneE2τ01+υ=v02−vh22.From equation (45), we can obtain the penetration depth of rockfall at any time when the velocity is vh, during the period in which impact velocity decreases from v0 to 0. When vh = 0, we obtain the maximum penetration depth of the rockfall Lmax.(46)4πτ03MRLmax2−Lmax33lneE2τ01+υ=v022.If the calculation result satisfiesLmax ≥ R, in equation (44), we set L = R.If the influence of rockfall velocity and acceleration during the impact process is considered and whenLmax ≤ R, we obtain(47)F1=−MvhdvhdL.Substituting equations (39) into (47), we obtain(48)dvhdL=−4τ0/3lneE/τ01+υ2πRL−πL2+3πρ0R4−R−L4/4R2vh2M−2πρ0/3R3−R−L3vh.Considering the initial conditions ofvh=v0 and L=0, when the calculation result for the penetration depth does not satisfy Lmax ≥ R, in equation (45), we set L = R. Using an analytical method, we can obtain the relation between L and vh.Using equation (48), we obtain(49)F=M4τ0/3lneE/2τ01+υ2πRL−πL2+3πρ0R4−R−L4/4R2vh2M−2πρ0/3R3−R−L3.We can substitute the values ofL and vh obtained using the analytical method into equation (49), to calculate the value of the impact force F.


## 4. Theoretical Derivation Based on Energy Conservation Model

In this section, the impact of rockfall on the soil cushion layer is simplified into the normal impact of a rigid ball on a semi-infinite soil cushion layer, and further derivation is based on the principle of conservation of energy.During the process of rockfall impact on the soil cushion layer, there are four main forms of energy:(1)
Initial kinetic energy of rockfall at the moment of contact with the cushion layer(2)
Plastic deformation and friction generate forms of energy dissipation(3)
Microgravitational potential energy of rockfall during the impact process


### 4.1. Surface Stress Analysis of Cylindrical Cavity

After rockfall impact on the cushion layer, we take the cross section of a cylindrical cavity with a radius ofrq, and the plastic zone and elastic zone are shown in Figure 5.Figure 5
Cylindrical cavity.We then consider the microbody as shown in Figure6 during mechanical analysis. The value of displacement of the microbody is positive in the outward direction, and the values of the stress σr∗ and σθ∗ and strain εr∗ and εθ∗ are positive with pressure.Figure 6
Force and displacement of cylindrical cavity microelement.So, we obtain the equation below:(50)∂σr∗∂r∗+σr∗−σθ∗r∗=0,where σr∗ is the radial stress of the cylindrical cavity microbody, σθ∗ is the hoop stress of the cylindrical cavity microbody, and θ∗ is the opening angle of cylindrical cavity microbody.Further, in the elastic zone, the geometric relation is given by equation (51):(51)εr∗=−∂s∗∂r∗,εθ∗=−s∗r∗,where εr∗ is the radial strain of cylindrical cavity microbody and εθ∗ is the hoop strain of cylindrical cavity microbody.For the medium of the elastic zone, according to Hooke’s law and the referenceCavity Expansion Methods in Geomechanics [38], the elastic relation expression is given by equation (52):(52)σr∗−σθ∗=Eεr∗−εθ∗1+υ.Further, according to the law of the conservation of mass, we obtain(53)12∂∂r∗r∗−s∗2=ρ0∗ρ∗r∗,where r∗ is the distance from the coordinate origin to the mass point of the cylindrical cavity microbody, s∗ is the displacement of the mass point of the cylindrical cavity microbody, ρ0∗ is the density before deformation of cylindrical cavity microbody, and ρ∗ is the density after deformation of cylindrical cavity microbody.Here, we consider that only the medium of the elastic zone is compressed, so we obtain(54)∂∂r∗r∗−s∗2=2r∗.Considering the boundary conditionsr∗ = rq∗, s∗ = rq∗ and taking the integral of equation (65), we obtain the expression for the displacement of a microbody:(55)s∗=r∗1−1−rq∗2r∗21/2,where rq∗ is the radius of the cylindrical cavity.Substituting equation (55) in equation (51), and given that the radius of the plastic zone of the soil cushion layer rq∗ is far less than the radius of elastic zone r∗, we obtain(56)εr∗−εθ∗=rq∗r∗21−rq∗r∗2−1/2≈rq∗r∗2.Combining equations (56) and (52), we obtain(57)σr∗−σθ∗=E1+υrq∗r∗2.Substituting equation (57) in (50) and considering the boundary conditions r∗ ⟶ ∞ and σr∗ = 0, we obtain(58)σr∗=E1+υrq∗r∗2.Assuming the soil cushion layer is an ideal elastic-plastic body, we obtain the failure criterion according to equation (59):(59)σr∗−σθ∗=2τ0=σy.Substituting equation (59) into (50), we obtain(60)∂σr∗∂r∗+2τ0r∗=0.Further, consider the boundary conditionsr∗ = rp∗, σr∗ = σrq∗, we obtain(61)σr∗=σrq∗+2τ0lnrq∗r∗,where σrq∗ is the radial stress on the cylindrical cavity surface.According to the boundary conditionsr∗ = rp∗ and σrq∗− = σrq∗+ for the assumption of a continuous medium in the elastic-plastic zone and combined with the yield criterion of equation (59), we obtain t equation (62) for the surface stress of a cylindrical cavity:(62)σrq∗=τ0lneE2τ01+υ,where σp∗ is the plastic region radius of cylindrical cavity.


### 4.2. Energy Dissipation due to Crater Formation

Using the slices method, the energy dissipation due to crater formation is equal to the sum of energy works of the surface stress of each cylindrical cavity in the stripe microbody, as shown in Figure7.(63)We=∫0rqf∗2πrq∗σrq∗drq∗=πσrq∗rqf∗2,where x∗ and y∗ are expressed by the equations below:(64)rqf∗=x∗,y∗=R−R2−x∗2,where x∗ is the x coordinate value of point G∗ on the rockfall surface and y∗ is the y coordinate value of point G∗ on the rockfall surface.Figure 7
Relative displacement diagram of thin layer.Substitute equation (64) into (63), we obtain(65)We=πτ02Ry∗−y∗2lneE2τ01+υ.Using the above equations, we can obtain the final energy dissipation due to crater formation of the rockfall:(66)Wte=∫0Lmax∗Wedz=∫0Lmax∗πτ02Ry∗−y∗2lneE2τ01+υdy∗,where We is the energy dissipated due to the deformation of thin layer of soil cushion layer, Wte is the total energy dissipated due to the deformation of soil cushion layer, and Lmax∗ is the maximum penetration depth of the rockfall by energy method.


### 4.3. Energy Dissipation due to Friction

Combined with Figures7 and 8, the slices method is adopted to analyze the internal energy dissipation due to the friction of rockfall and soil cushion layer.Figure 8
Stress analysis of thin layer.In Figure9, the friction of thin layer of microsoil cushion layer and rockfall Ff is(67)Ff=2πrq∗μσrq∗cosβ,where Ff is the friction between the thin layer of soil cushion layer and the rockfall, μ is the friction coefficient between soil cushion layer and the rockfall, and β is the acute angle between the inclined plane of the thin soil layer with the vertical line.Figure 9
“Impact velocity-contact radius” curve.Combined with Figure7, when the thin layer of soil cushion layer has the displacement rqf∗ relative to rockfall, the energy work of friction Ff is(68)Whf=2πμ∫0rqf∗rq∗σrq∗dzcos2β,where rqf∗ is the horizontal movement distance of the thin soil layer, dz is the thickness of the soil layer and Whf is the energy dissipated due to friction between the rockfall and the soil cushion layer.Further, the expression for the side profile of rockfall is(69)z=frq∗drq∗.Combining equations (62), (68), and (69), we obtain(70)Whf=2πμτ0lneE2τ01+υ∫0rqf∗rq∗f′rq∗drq∗cos2β.In equation (70), drq∗ is the distance of movement of a thin layer of soil cushion relative to the rockfall. This can be determined using equation (71):(71)dy∗=x∗R2−x∗2dx∗.To determinecosβ as represented in Figure 7, a tangent line L2 is made through the contact point of the thin layer of soil cushion and rockfall and the included angle of L2. The horizontal axis is β, so we obtain(72)y2′=tanχ,where y2′ is the slope of the tangent line L2 and χ is the acute angle between the inclined plane of the thin soil layer with the horizontal line.Further, based on the trigonometric function relation,(73)sec2χ−1=tan2χ,secχ=cos−1χ,and combining the equations (72) and (73), we obtain(74)cosχ=11+y2′2.Combining equations (71) and (72), we obtain(75)y2′=dy∗dx∗=x∗R2−x∗2.Combining equations (74) and (75), we obtain(76)sin2χ=x∗2R2.According to Figure8, we can write the equation (77):(77)χ+β=π2,and we can state that(78)cos2β=x∗2R2.According to the above equations, we obtain the final expression of equation (79):(79)Whf=2πR2μτ0arcsin2Ry∗−y∗2RlneE2τ01+υ.Therefore, when rockfall has the maximum penetration depth ofLmax∗, the energy dissipation due to friction is(80)Wth=∫0Lmax∗Whfdy∗,where Wth is the total energy dissipated due to friction between the rockfall and the soil cushion layer.


### 4.4. Energy Dissipation of Shed Tunnel Structure

Under the impact effect of rockfall load, the roof of a shed tunnel structure will have slight deflection deformation due to rockfall impact, and therefore, energy is dissipated. The deflection of the roof of the shed tunnel under an impact forcePf is [39](81)δf=Pfl∗348E1I,where δf is the roof deflection of the shed tunnel structure, Pf is the Given rockfall impact load, E1 is the elasticity modulus of the roof of the shed tunnel structure, I is the inertia moment of the roof of the shed tunnel structure, and l∗ is the span of the shed tunnel structure.According to structural mechanics and equation (81), we obtain the relevant energy of deformation as(82)Wf=Pf2l∗396E1I.Further,Wf can be determined using equation (83):(83)Wf=π2σy22RLmax−Lmax∗22l∗396E1I,where σy is the yield strength of soil cushion layer and Wf is the total energy dissipated due to the deflection deformation of the roof of the shed tunnel.


### 4.5. Penetration Depth

The initial kinetic energy of rockfall during the impact process isMv02/2; the energy of deformation due to crater Wte shall be calculated using equation (66), the energy of friction Wth using equation (80), and the energy of deflection deformation of the roof of the shed tunnel Wf using equation (83). The work done by gravity during rockfall impact process is MgLmax∗; therefore, according to the law of energy conservation, we obtain:(84)Mv022=Wte+Wth+Wf−MgLmax∗.There is only one unknown parameterLmax∗ in equation (76), so we obtain the solution using equation (84). The value of the impact force F is obtained based on Lmax∗.


### 4.6. Impact Force

According to the principle of work and energy, the relation expression during the contact of the rockfall and the soil cushion layer is(85)Fa∗⋅Lmax∗=Mv022+MgLmax∗,where Fa∗ is the average impact force of rockfall.TheFa∗ value determined using equation (85) is the average of the impact force. However, the rockfall impact process is an impulse type, Yang et al. [37], and with the vibration curve of the impact force over time during the impact process that simulates the process of changes of impact force approximately in the form of a sine function, we obtain(86)F=Fmax∗sinω0t1,where Fmax∗ is the maximum impact force of the rockfall, ω0 is the angular frequency, and t1 is the time.Since the impulse due to actual impact force is equal to that of the average impact force of rockfall, we obtain equation (87):(87)Fa∗⋅tw=∫0twFt1dt1.Combining the equations (86) and (87), we obtain(88)Fa∗=∫0twFmax∗sinπt1/twdt1tw.Further, we have(89)Fmax∗=πFa∗2,and the time duration of rockfall impact is(90)tw=2Mv0Fa∗.


## 4.1. Surface Stress Analysis of Cylindrical Cavity

After rockfall impact on the cushion layer, we take the cross section of a cylindrical cavity with a radius ofrq, and the plastic zone and elastic zone are shown in Figure 5.Figure 5
Cylindrical cavity.We then consider the microbody as shown in Figure6 during mechanical analysis. The value of displacement of the microbody is positive in the outward direction, and the values of the stress σr∗ and σθ∗ and strain εr∗ and εθ∗ are positive with pressure.Figure 6
Force and displacement of cylindrical cavity microelement.So, we obtain the equation below:(50)∂σr∗∂r∗+σr∗−σθ∗r∗=0,where σr∗ is the radial stress of the cylindrical cavity microbody, σθ∗ is the hoop stress of the cylindrical cavity microbody, and θ∗ is the opening angle of cylindrical cavity microbody.Further, in the elastic zone, the geometric relation is given by equation (51):(51)εr∗=−∂s∗∂r∗,εθ∗=−s∗r∗,where εr∗ is the radial strain of cylindrical cavity microbody and εθ∗ is the hoop strain of cylindrical cavity microbody.For the medium of the elastic zone, according to Hooke’s law and the referenceCavity Expansion Methods in Geomechanics [38], the elastic relation expression is given by equation (52):(52)σr∗−σθ∗=Eεr∗−εθ∗1+υ.Further, according to the law of the conservation of mass, we obtain(53)12∂∂r∗r∗−s∗2=ρ0∗ρ∗r∗,where r∗ is the distance from the coordinate origin to the mass point of the cylindrical cavity microbody, s∗ is the displacement of the mass point of the cylindrical cavity microbody, ρ0∗ is the density before deformation of cylindrical cavity microbody, and ρ∗ is the density after deformation of cylindrical cavity microbody.Here, we consider that only the medium of the elastic zone is compressed, so we obtain(54)∂∂r∗r∗−s∗2=2r∗.Considering the boundary conditionsr∗ = rq∗, s∗ = rq∗ and taking the integral of equation (65), we obtain the expression for the displacement of a microbody:(55)s∗=r∗1−1−rq∗2r∗21/2,where rq∗ is the radius of the cylindrical cavity.Substituting equation (55) in equation (51), and given that the radius of the plastic zone of the soil cushion layer rq∗ is far less than the radius of elastic zone r∗, we obtain(56)εr∗−εθ∗=rq∗r∗21−rq∗r∗2−1/2≈rq∗r∗2.Combining equations (56) and (52), we obtain(57)σr∗−σθ∗=E1+υrq∗r∗2.Substituting equation (57) in (50) and considering the boundary conditions r∗ ⟶ ∞ and σr∗ = 0, we obtain(58)σr∗=E1+υrq∗r∗2.Assuming the soil cushion layer is an ideal elastic-plastic body, we obtain the failure criterion according to equation (59):(59)σr∗−σθ∗=2τ0=σy.Substituting equation (59) into (50), we obtain(60)∂σr∗∂r∗+2τ0r∗=0.Further, consider the boundary conditionsr∗ = rp∗, σr∗ = σrq∗, we obtain(61)σr∗=σrq∗+2τ0lnrq∗r∗,where σrq∗ is the radial stress on the cylindrical cavity surface.According to the boundary conditionsr∗ = rp∗ and σrq∗− = σrq∗+ for the assumption of a continuous medium in the elastic-plastic zone and combined with the yield criterion of equation (59), we obtain t equation (62) for the surface stress of a cylindrical cavity:(62)σrq∗=τ0lneE2τ01+υ,where σp∗ is the plastic region radius of cylindrical cavity.


## 4.2. Energy Dissipation due to Crater Formation

Using the slices method, the energy dissipation due to crater formation is equal to the sum of energy works of the surface stress of each cylindrical cavity in the stripe microbody, as shown in Figure7.(63)We=∫0rqf∗2πrq∗σrq∗drq∗=πσrq∗rqf∗2,where x∗ and y∗ are expressed by the equations below:(64)rqf∗=x∗,y∗=R−R2−x∗2,where x∗ is the x coordinate value of point G∗ on the rockfall surface and y∗ is the y coordinate value of point G∗ on the rockfall surface.Figure 7
Relative displacement diagram of thin layer.Substitute equation (64) into (63), we obtain(65)We=πτ02Ry∗−y∗2lneE2τ01+υ.Using the above equations, we can obtain the final energy dissipation due to crater formation of the rockfall:(66)Wte=∫0Lmax∗Wedz=∫0Lmax∗πτ02Ry∗−y∗2lneE2τ01+υdy∗,where We is the energy dissipated due to the deformation of thin layer of soil cushion layer, Wte is the total energy dissipated due to the deformation of soil cushion layer, and Lmax∗ is the maximum penetration depth of the rockfall by energy method.


## 4.3. Energy Dissipation due to Friction

Combined with Figures7 and 8, the slices method is adopted to analyze the internal energy dissipation due to the friction of rockfall and soil cushion layer.Figure 8
Stress analysis of thin layer.In Figure9, the friction of thin layer of microsoil cushion layer and rockfall Ff is(67)Ff=2πrq∗μσrq∗cosβ,where Ff is the friction between the thin layer of soil cushion layer and the rockfall, μ is the friction coefficient between soil cushion layer and the rockfall, and β is the acute angle between the inclined plane of the thin soil layer with the vertical line.Figure 9
“Impact velocity-contact radius” curve.Combined with Figure7, when the thin layer of soil cushion layer has the displacement rqf∗ relative to rockfall, the energy work of friction Ff is(68)Whf=2πμ∫0rqf∗rq∗σrq∗dzcos2β,where rqf∗ is the horizontal movement distance of the thin soil layer, dz is the thickness of the soil layer and Whf is the energy dissipated due to friction between the rockfall and the soil cushion layer.Further, the expression for the side profile of rockfall is(69)z=frq∗drq∗.Combining equations (62), (68), and (69), we obtain(70)Whf=2πμτ0lneE2τ01+υ∫0rqf∗rq∗f′rq∗drq∗cos2β.In equation (70), drq∗ is the distance of movement of a thin layer of soil cushion relative to the rockfall. This can be determined using equation (71):(71)dy∗=x∗R2−x∗2dx∗.To determinecosβ as represented in Figure 7, a tangent line L2 is made through the contact point of the thin layer of soil cushion and rockfall and the included angle of L2. The horizontal axis is β, so we obtain(72)y2′=tanχ,where y2′ is the slope of the tangent line L2 and χ is the acute angle between the inclined plane of the thin soil layer with the horizontal line.Further, based on the trigonometric function relation,(73)sec2χ−1=tan2χ,secχ=cos−1χ,and combining the equations (72) and (73), we obtain(74)cosχ=11+y2′2.Combining equations (71) and (72), we obtain(75)y2′=dy∗dx∗=x∗R2−x∗2.Combining equations (74) and (75), we obtain(76)sin2χ=x∗2R2.According to Figure8, we can write the equation (77):(77)χ+β=π2,and we can state that(78)cos2β=x∗2R2.According to the above equations, we obtain the final expression of equation (79):(79)Whf=2πR2μτ0arcsin2Ry∗−y∗2RlneE2τ01+υ.Therefore, when rockfall has the maximum penetration depth ofLmax∗, the energy dissipation due to friction is(80)Wth=∫0Lmax∗Whfdy∗,where Wth is the total energy dissipated due to friction between the rockfall and the soil cushion layer.


## 4.4. Energy Dissipation of Shed Tunnel Structure

Under the impact effect of rockfall load, the roof of a shed tunnel structure will have slight deflection deformation due to rockfall impact, and therefore, energy is dissipated. The deflection of the roof of the shed tunnel under an impact forcePf is [39](81)δf=Pfl∗348E1I,where δf is the roof deflection of the shed tunnel structure, Pf is the Given rockfall impact load, E1 is the elasticity modulus of the roof of the shed tunnel structure, I is the inertia moment of the roof of the shed tunnel structure, and l∗ is the span of the shed tunnel structure.According to structural mechanics and equation (81), we obtain the relevant energy of deformation as(82)Wf=Pf2l∗396E1I.Further,Wf can be determined using equation (83):(83)Wf=π2σy22RLmax−Lmax∗22l∗396E1I,where σy is the yield strength of soil cushion layer and Wf is the total energy dissipated due to the deflection deformation of the roof of the shed tunnel.


## 4.5. Penetration Depth

The initial kinetic energy of rockfall during the impact process isMv02/2; the energy of deformation due to crater Wte shall be calculated using equation (66), the energy of friction Wth using equation (80), and the energy of deflection deformation of the roof of the shed tunnel Wf using equation (83). The work done by gravity during rockfall impact process is MgLmax∗; therefore, according to the law of energy conservation, we obtain:(84)Mv022=Wte+Wth+Wf−MgLmax∗.There is only one unknown parameterLmax∗ in equation (76), so we obtain the solution using equation (84). The value of the impact force F is obtained based on Lmax∗.


## 4.6. Impact Force

According to the principle of work and energy, the relation expression during the contact of the rockfall and the soil cushion layer is(85)Fa∗⋅Lmax∗=Mv022+MgLmax∗,where Fa∗ is the average impact force of rockfall.TheFa∗ value determined using equation (85) is the average of the impact force. However, the rockfall impact process is an impulse type, Yang et al. [37], and with the vibration curve of the impact force over time during the impact process that simulates the process of changes of impact force approximately in the form of a sine function, we obtain(86)F=Fmax∗sinω0t1,where Fmax∗ is the maximum impact force of the rockfall, ω0 is the angular frequency, and t1 is the time.Since the impulse due to actual impact force is equal to that of the average impact force of rockfall, we obtain equation (87):(87)Fa∗⋅tw=∫0twFt1dt1.Combining the equations (86) and (87), we obtain(88)Fa∗=∫0twFmax∗sinπt1/twdt1tw.Further, we have(89)Fmax∗=πFa∗2,and the time duration of rockfall impact is(90)tw=2Mv0Fa∗.


## 5. LS-DYA Numerical Simulation Analysis

In order to evaluate the theoretical calculation result and explore the detailed process of rockfall impacting on the soil cushion layer, LS-DYA software was used to analyze the impact process. A typical working condition is taken as an example: size of the soil cushion layer is 6 m (length) × 6 m (width) × 3 m (height); the radius of the rockfall is 0.5 m; and impact velocity of rockfall is 20 m/s. Face-to-face contact is set between the rockfall and the soil cushion layer, and grid refinement is conducted within the range of 1.5 m × 1.5 m in the center of the soil cushion layer. All the nodes in the bottom of the cushion are restricted in three directions: 0.06 s is set as the total time of the analysis step and 200 steps are chosen as the total calculation steps. The calculation parameters are shown in Tables1 and 2.Table 1
Slope soil parameters of the calculation model.
Volume-weight (kN/m3)
Elasticity modulus (MPa)
Poisson’s ratio
Yield strength (kPa)
Friction coefficient of rock and soil
20.5
40
0.3
400
0.4Table 2
Rockfall and shed tunnel beam parameters of the calculation model.
Rockfall
Shed tunnel beam
Volume-weight (kN/m3)
Elasticity modulus (MPa)
Radius (m)
Falling height (m)
Elasticity modulus (MPa)
Poisson’s ratio
Profile moments of inertia (m4)
Span (m)
25
28
0.5
20
28.5
0.19
0.057
10Figure10 shows the results of the LS-DYNA numerical model and the numerical simulation results at 0.02 s, 0.04 s, and 0.06 s. Figure 11 shows the variation curves of impact force and penetration depth of the rockfall with time. After rockfall impact on soil cushion layer, the impact force of the rockfall reaches a maximum value of 1559801 N at 0.004 s (14/200 steps), then reduces gradually, and finally reaches a value of 0 at 0.44 s (147/200 steps). However, the curve shows a stronger vibration between 0.018 s and 0.026 s. The penetration depth of the rockfall gradually increases in the range 0∼0.035 s and reaches a maximum value 0.35 m at 0.035 s (117/200 step), but approximately 9% spring back quantity is generated between 0.035 s and 0.05 s, forming the ultimate penetration depth of 0.312 m.Figure 10
Numerical model and results.
(a)
(b)
(c)
(d)Figure 11
Impact force and penetration depth of simulation.


## 6. Engineering Example

### 6.1. Influence of Rockfall Velocity

We have investigated the effects of the collisions of a rockfall block, 0.5 m in radius, with the soil cushion layer when the impact speed of the block is 4, 6, 8, 10, 12, 16, 20, and 24 m/s.


#### 6.1.1. “Impact Velocity-Penetration Depth” Relation Analysis

Figure12 shows the plot of “impact velocity vs penetration depth.” With an increase in the rockfall impact velocity, the results of these eight methods all reveal a trend of linear growth. The calculation results for the rockfall penetration depth obtained using the cavity expansion, energy conservation, numerical simulation, Australian, and B. S. Guan methods all demonstrated good consistency. Comparatively, the Hertz method has smaller results, while the results from the tunnel manual method and subgrade method are obviously larger.Figure 12
“Impact velocity-penetration depth” curve.


#### 6.1.2. “Impact Velocity-Contact Radius” Relation Analysis

Figure9 shows a plot of “impact velocity vs contact radius.” The “contact radius” in this study is defined as follows: when a rock falls on the earth, there is a certain contact area between them. The contact area is calculated and converted into an equivalent circle, and the radius of the circle is the contact radius. From Figure 9, the variation trend for the results of the eight methods is basically consistent with that in Figure 12. In addition, the calculation results for the contact radius based on the cavity expansion, energy conservation, numerical simulation, Australia, and B. S. Guan methods reveal better consistency.


#### 6.1.3. “Impact Velocity vs Impact Force” Relation Analysis

Figure13 represents the plot of “Impact velocity vs impact force.” With the increase in rockfall impact velocity, the calculation results for the 10 methods all reveal a linear growth trend. The calculation results for the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods demonstrate better consistency. In contrast, the calculation results for the B. S. Guan, tunnel manual, and subgrade methods are smaller overall, while the calculation result for the Hertz method is slightly larger.Figure 13
“Impact velocity vs impact force” curve.


#### 6.1.4. “Penetration Depth-Impact Force” Relation Analysis

Figure14 shows the curve of “penetration depth vs impact force.” The difference between the curves for the eight methods is evident. The calculation results based on the cavity expansion, energy conservation, numerical simulation, and Australian methods demonstrate better consistency. In comparison, the B. S. Guan method exhibited a slightly smaller impact force with an appropriate penetration depth, the tunnel manual method, and subgrade method exhibited weaker impact forces with larger penetration depth, while the Hertz Method resulted in a greater force with a smaller penetration depth.Figure 14
“Penetration depth vs impact force” curve.


#### 6.1.5. “Contact Radius vs Impact Force” Relation Analysis

Figure15 shows the curve for the “contact radius vs impact force.” According to the figure, the difference between the plots for the eight methods is basically the same as that of Figure 14. The calculation results for the cavity expansion, energy conservation, numerical simulation, and Australian methods demonstrated better consistency and moderate positions.Figure 15
“Contact radius-impact force” curve.


### 6.2. Influence of Rockfall Radius

The impact velocity of 20 m/s and rockfall radii of 0.2 m, 0.3 m, 0.4 m, 0.5 m, 0.6 m, 0.7 m, and 0.8 m were considered to analyze the impact effect of rockfall on the soil cushion layer.


#### 6.2.1. “Rockfall Radius vs Penetration Depth” Relation Analysis

Figure16 shows the curve of “rockfall radius vs penetration depth.” The calculation results for the eight methods all reveal linear growth trends. The calculation results for the penetration depth based on the cavity expansion, energy conservation, numerical simulation, Australian, and B. S. Guan methods demonstrated better consistency. In contrast, the results for the Hertz method are smaller, while that of the tunnel manual method and subgrade method are obviously larger.Figure 16
“Rockfall radius-penetration depth” curve.


#### 6.2.2. “Rockfall Radius vs Impact Force” Relation Analysis

Figure17 represents the curve of “rockfall radius vs impact force.” According to Figure 17, with the increase in the rockfall radius, the calculation results for the 10 methods all reveal a linear growth trend. The calculation results for the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods demonstrated better consistency, while the calculation results for the B. S. Guan, tunnel manual, and subgrade methods are all smaller. In addition, the calculation result for the Hertz method is slightly larger.Figure 17
“Rockfall radius vs impact force” curve.


#### 6.2.3. “Rockfall Radius vs Contact Radius” Relation Analysis

Figure18 shows the curve of “rolling stone radius vs contact radius.” As shown in this figure, the variation trend of the results for the 8 methods is basically consistent with that shown in Figure 7. Moreover, the calculation results for the contact radius obtained using the cavity expansion, energy conservation, Australian, numerical simulation, and B. S. Guan methods are relatively consistent.Figure 18
“penetration depth vs impact force” curve.As observed from the preceding analysis, the calculation results for the rockfall impact force and penetration depth based on the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods reveal better consistency, which explains why the cavity expansion method and energy conservation method have better accuracy. In comparison, the penetration depth obtained for the B. S. Guan method is accurate, but the impact force is slightly smaller. The calculation results for the impact force obtained using the tunnel manual method and subgrade method are smaller, while the calculation result for the penetration depth is larger. In addition, the calculation results for the penetration depth based on the Hertz method are smaller, while that of the impact force is larger.


### 6.3. Energy Dissipation Analysis

#### 6.3.1. Analysis on Impact Velocity Influence on Energy Dissipation of Rockfall

We analyzed the energy dissipation of the rockfall impact with the rockfall radius of 0.5 m and impact velocity of 4, 6, 8, 10, 12, 16, 20, and 24 m/s.Figure19 represents the “influence of impact velocity on energy dissipation”, from which it can be seen that with an increase in the impact velocity, the energy dissipation due to crater formation and due to friction increase rapidly. The increase rate and tendency of the two curves are basically consistent. Furthermore, the gross energy dissipation due to friction is always higher than that due to crater formation, and the energy dissipation of the shed tunnel structure is small overall. Figure 19 also shows the proportion of three kinds of energy dissipations under impact velocity, from which it can be seen that the proportion of energy dissipation due to crater formation gradually rises from approximately 25% to approximately 40% and then becomes stable. Energy dissipation due to friction decreases from approximately 70% to 60% and becomes stable. The energy dissipation of the shed tunnel structure also decreased from approximately 5% to approximately 0%. Therefore, the energy dissipation of the shed tunnel structure could be ignored.Figure 19
Influence of impact velocity on energy dissipation.


#### 6.3.2. Analysis on the Influence of the Rockfall Radius on Energy Dissipation

Futher, the impact velocity of 20 m/s and rockfall radii of 0.2 m, 0.3 m, 0.4 m, 0.5 m, 0.6 m, 0.7 m, and 0.8 m were considered to analyze the energy dissipation of the rockfall impact.Figure20 represents the “Influence of rockfall radius on energy dissipation.” Based on this figure, it can be seen that the energy dissipation due to crater formation and friction increases rapidly with an increase in the rockfall radius. In addition, the increasing tendency of the two curves is consistent. Furthermore, energy dissipation due to friction is always higher than that due to crater formation. Figure 20 also shows the proportion of energy dissipation for different rockfall radii. It can be seen that energy dissipation due to crater formation is maintained at approximately 40% and that friction is approximately 60%, which further indicates that the energy dissipation of a conventional shed tunnel structure can be ignored.Figure 20
Influence of rockfall radius on energy dissipation.


## 6.1. Influence of Rockfall Velocity

We have investigated the effects of the collisions of a rockfall block, 0.5 m in radius, with the soil cushion layer when the impact speed of the block is 4, 6, 8, 10, 12, 16, 20, and 24 m/s.


### 6.1.1. “Impact Velocity-Penetration Depth” Relation Analysis

Figure12 shows the plot of “impact velocity vs penetration depth.” With an increase in the rockfall impact velocity, the results of these eight methods all reveal a trend of linear growth. The calculation results for the rockfall penetration depth obtained using the cavity expansion, energy conservation, numerical simulation, Australian, and B. S. Guan methods all demonstrated good consistency. Comparatively, the Hertz method has smaller results, while the results from the tunnel manual method and subgrade method are obviously larger.Figure 12
“Impact velocity-penetration depth” curve.


### 6.1.2. “Impact Velocity-Contact Radius” Relation Analysis

Figure9 shows a plot of “impact velocity vs contact radius.” The “contact radius” in this study is defined as follows: when a rock falls on the earth, there is a certain contact area between them. The contact area is calculated and converted into an equivalent circle, and the radius of the circle is the contact radius. From Figure 9, the variation trend for the results of the eight methods is basically consistent with that in Figure 12. In addition, the calculation results for the contact radius based on the cavity expansion, energy conservation, numerical simulation, Australia, and B. S. Guan methods reveal better consistency.


### 6.1.3. “Impact Velocity vs Impact Force” Relation Analysis

Figure13 represents the plot of “Impact velocity vs impact force.” With the increase in rockfall impact velocity, the calculation results for the 10 methods all reveal a linear growth trend. The calculation results for the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods demonstrate better consistency. In contrast, the calculation results for the B. S. Guan, tunnel manual, and subgrade methods are smaller overall, while the calculation result for the Hertz method is slightly larger.Figure 13
“Impact velocity vs impact force” curve.


### 6.1.4. “Penetration Depth-Impact Force” Relation Analysis

Figure14 shows the curve of “penetration depth vs impact force.” The difference between the curves for the eight methods is evident. The calculation results based on the cavity expansion, energy conservation, numerical simulation, and Australian methods demonstrate better consistency. In comparison, the B. S. Guan method exhibited a slightly smaller impact force with an appropriate penetration depth, the tunnel manual method, and subgrade method exhibited weaker impact forces with larger penetration depth, while the Hertz Method resulted in a greater force with a smaller penetration depth.Figure 14
“Penetration depth vs impact force” curve.


### 6.1.5. “Contact Radius vs Impact Force” Relation Analysis

Figure15 shows the curve for the “contact radius vs impact force.” According to the figure, the difference between the plots for the eight methods is basically the same as that of Figure 14. The calculation results for the cavity expansion, energy conservation, numerical simulation, and Australian methods demonstrated better consistency and moderate positions.Figure 15
“Contact radius-impact force” curve.


## 6.1.1. “Impact Velocity-Penetration Depth” Relation Analysis

Figure12 shows the plot of “impact velocity vs penetration depth.” With an increase in the rockfall impact velocity, the results of these eight methods all reveal a trend of linear growth. The calculation results for the rockfall penetration depth obtained using the cavity expansion, energy conservation, numerical simulation, Australian, and B. S. Guan methods all demonstrated good consistency. Comparatively, the Hertz method has smaller results, while the results from the tunnel manual method and subgrade method are obviously larger.Figure 12
“Impact velocity-penetration depth” curve.


## 6.1.2. “Impact Velocity-Contact Radius” Relation Analysis

Figure9 shows a plot of “impact velocity vs contact radius.” The “contact radius” in this study is defined as follows: when a rock falls on the earth, there is a certain contact area between them. The contact area is calculated and converted into an equivalent circle, and the radius of the circle is the contact radius. From Figure 9, the variation trend for the results of the eight methods is basically consistent with that in Figure 12. In addition, the calculation results for the contact radius based on the cavity expansion, energy conservation, numerical simulation, Australia, and B. S. Guan methods reveal better consistency.


## 6.1.3. “Impact Velocity vs Impact Force” Relation Analysis

Figure13 represents the plot of “Impact velocity vs impact force.” With the increase in rockfall impact velocity, the calculation results for the 10 methods all reveal a linear growth trend. The calculation results for the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods demonstrate better consistency. In contrast, the calculation results for the B. S. Guan, tunnel manual, and subgrade methods are smaller overall, while the calculation result for the Hertz method is slightly larger.Figure 13
“Impact velocity vs impact force” curve.


## 6.1.4. “Penetration Depth-Impact Force” Relation Analysis

Figure14 shows the curve of “penetration depth vs impact force.” The difference between the curves for the eight methods is evident. The calculation results based on the cavity expansion, energy conservation, numerical simulation, and Australian methods demonstrate better consistency. In comparison, the B. S. Guan method exhibited a slightly smaller impact force with an appropriate penetration depth, the tunnel manual method, and subgrade method exhibited weaker impact forces with larger penetration depth, while the Hertz Method resulted in a greater force with a smaller penetration depth.Figure 14
“Penetration depth vs impact force” curve.


## 6.1.5. “Contact Radius vs Impact Force” Relation Analysis

Figure15 shows the curve for the “contact radius vs impact force.” According to the figure, the difference between the plots for the eight methods is basically the same as that of Figure 14. The calculation results for the cavity expansion, energy conservation, numerical simulation, and Australian methods demonstrated better consistency and moderate positions.Figure 15
“Contact radius-impact force” curve.


## 6.2. Influence of Rockfall Radius

The impact velocity of 20 m/s and rockfall radii of 0.2 m, 0.3 m, 0.4 m, 0.5 m, 0.6 m, 0.7 m, and 0.8 m were considered to analyze the impact effect of rockfall on the soil cushion layer.


### 6.2.1. “Rockfall Radius vs Penetration Depth” Relation Analysis

Figure16 shows the curve of “rockfall radius vs penetration depth.” The calculation results for the eight methods all reveal linear growth trends. The calculation results for the penetration depth based on the cavity expansion, energy conservation, numerical simulation, Australian, and B. S. Guan methods demonstrated better consistency. In contrast, the results for the Hertz method are smaller, while that of the tunnel manual method and subgrade method are obviously larger.Figure 16
“Rockfall radius-penetration depth” curve.


### 6.2.2. “Rockfall Radius vs Impact Force” Relation Analysis

Figure17 represents the curve of “rockfall radius vs impact force.” According to Figure 17, with the increase in the rockfall radius, the calculation results for the 10 methods all reveal a linear growth trend. The calculation results for the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods demonstrated better consistency, while the calculation results for the B. S. Guan, tunnel manual, and subgrade methods are all smaller. In addition, the calculation result for the Hertz method is slightly larger.Figure 17
“Rockfall radius vs impact force” curve.


### 6.2.3. “Rockfall Radius vs Contact Radius” Relation Analysis

Figure18 shows the curve of “rolling stone radius vs contact radius.” As shown in this figure, the variation trend of the results for the 8 methods is basically consistent with that shown in Figure 7. Moreover, the calculation results for the contact radius obtained using the cavity expansion, energy conservation, Australian, numerical simulation, and B. S. Guan methods are relatively consistent.Figure 18
“penetration depth vs impact force” curve.As observed from the preceding analysis, the calculation results for the rockfall impact force and penetration depth based on the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods reveal better consistency, which explains why the cavity expansion method and energy conservation method have better accuracy. In comparison, the penetration depth obtained for the B. S. Guan method is accurate, but the impact force is slightly smaller. The calculation results for the impact force obtained using the tunnel manual method and subgrade method are smaller, while the calculation result for the penetration depth is larger. In addition, the calculation results for the penetration depth based on the Hertz method are smaller, while that of the impact force is larger.


## 6.2.1. “Rockfall Radius vs Penetration Depth” Relation Analysis

Figure16 shows the curve of “rockfall radius vs penetration depth.” The calculation results for the eight methods all reveal linear growth trends. The calculation results for the penetration depth based on the cavity expansion, energy conservation, numerical simulation, Australian, and B. S. Guan methods demonstrated better consistency. In contrast, the results for the Hertz method are smaller, while that of the tunnel manual method and subgrade method are obviously larger.Figure 16
“Rockfall radius-penetration depth” curve.


## 6.2.2. “Rockfall Radius vs Impact Force” Relation Analysis

Figure17 represents the curve of “rockfall radius vs impact force.” According to Figure 17, with the increase in the rockfall radius, the calculation results for the 10 methods all reveal a linear growth trend. The calculation results for the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods demonstrated better consistency, while the calculation results for the B. S. Guan, tunnel manual, and subgrade methods are all smaller. In addition, the calculation result for the Hertz method is slightly larger.Figure 17
“Rockfall radius vs impact force” curve.


## 6.2.3. “Rockfall Radius vs Contact Radius” Relation Analysis

Figure18 shows the curve of “rolling stone radius vs contact radius.” As shown in this figure, the variation trend of the results for the 8 methods is basically consistent with that shown in Figure 7. Moreover, the calculation results for the contact radius obtained using the cavity expansion, energy conservation, Australian, numerical simulation, and B. S. Guan methods are relatively consistent.Figure 18
“penetration depth vs impact force” curve.As observed from the preceding analysis, the calculation results for the rockfall impact force and penetration depth based on the cavity expansion, energy conservation, numerical simulation, Japanese, Swiss, and Australian methods reveal better consistency, which explains why the cavity expansion method and energy conservation method have better accuracy. In comparison, the penetration depth obtained for the B. S. Guan method is accurate, but the impact force is slightly smaller. The calculation results for the impact force obtained using the tunnel manual method and subgrade method are smaller, while the calculation result for the penetration depth is larger. In addition, the calculation results for the penetration depth based on the Hertz method are smaller, while that of the impact force is larger.


## 6.3. Energy Dissipation Analysis

### 6.3.1. Analysis on Impact Velocity Influence on Energy Dissipation of Rockfall

We analyzed the energy dissipation of the rockfall impact with the rockfall radius of 0.5 m and impact velocity of 4, 6, 8, 10, 12, 16, 20, and 24 m/s.Figure19 represents the “influence of impact velocity on energy dissipation”, from which it can be seen that with an increase in the impact velocity, the energy dissipation due to crater formation and due to friction increase rapidly. The increase rate and tendency of the two curves are basically consistent. Furthermore, the gross energy dissipation due to friction is always higher than that due to crater formation, and the energy dissipation of the shed tunnel structure is small overall. Figure 19 also shows the proportion of three kinds of energy dissipations under impact velocity, from which it can be seen that the proportion of energy dissipation due to crater formation gradually rises from approximately 25% to approximately 40% and then becomes stable. Energy dissipation due to friction decreases from approximately 70% to 60% and becomes stable. The energy dissipation of the shed tunnel structure also decreased from approximately 5% to approximately 0%. Therefore, the energy dissipation of the shed tunnel structure could be ignored.Figure 19
Influence of impact velocity on energy dissipation.


### 6.3.2. Analysis on the Influence of the Rockfall Radius on Energy Dissipation

Futher, the impact velocity of 20 m/s and rockfall radii of 0.2 m, 0.3 m, 0.4 m, 0.5 m, 0.6 m, 0.7 m, and 0.8 m were considered to analyze the energy dissipation of the rockfall impact.Figure20 represents the “Influence of rockfall radius on energy dissipation.” Based on this figure, it can be seen that the energy dissipation due to crater formation and friction increases rapidly with an increase in the rockfall radius. In addition, the increasing tendency of the two curves is consistent. Furthermore, energy dissipation due to friction is always higher than that due to crater formation. Figure 20 also shows the proportion of energy dissipation for different rockfall radii. It can be seen that energy dissipation due to crater formation is maintained at approximately 40% and that friction is approximately 60%, which further indicates that the energy dissipation of a conventional shed tunnel structure can be ignored.Figure 20
Influence of rockfall radius on energy dissipation.


## 6.3.1. Analysis on Impact Velocity Influence on Energy Dissipation of Rockfall

We analyzed the energy dissipation of the rockfall impact with the rockfall radius of 0.5 m and impact velocity of 4, 6, 8, 10, 12, 16, 20, and 24 m/s.Figure19 represents the “influence of impact velocity on energy dissipation”, from which it can be seen that with an increase in the impact velocity, the energy dissipation due to crater formation and due to friction increase rapidly. The increase rate and tendency of the two curves are basically consistent. Furthermore, the gross energy dissipation due to friction is always higher than that due to crater formation, and the energy dissipation of the shed tunnel structure is small overall. Figure 19 also shows the proportion of three kinds of energy dissipations under impact velocity, from which it can be seen that the proportion of energy dissipation due to crater formation gradually rises from approximately 25% to approximately 40% and then becomes stable. Energy dissipation due to friction decreases from approximately 70% to 60% and becomes stable. The energy dissipation of the shed tunnel structure also decreased from approximately 5% to approximately 0%. Therefore, the energy dissipation of the shed tunnel structure could be ignored.Figure 19
Influence of impact velocity on energy dissipation.


## 6.3.2. Analysis on the Influence of the Rockfall Radius on Energy Dissipation

Futher, the impact velocity of 20 m/s and rockfall radii of 0.2 m, 0.3 m, 0.4 m, 0.5 m, 0.6 m, 0.7 m, and 0.8 m were considered to analyze the energy dissipation of the rockfall impact.Figure20 represents the “Influence of rockfall radius on energy dissipation.” Based on this figure, it can be seen that the energy dissipation due to crater formation and friction increases rapidly with an increase in the rockfall radius. In addition, the increasing tendency of the two curves is consistent. Furthermore, energy dissipation due to friction is always higher than that due to crater formation. Figure 20 also shows the proportion of energy dissipation for different rockfall radii. It can be seen that energy dissipation due to crater formation is maintained at approximately 40% and that friction is approximately 60%, which further indicates that the energy dissipation of a conventional shed tunnel structure can be ignored.Figure 20
Influence of rockfall radius on energy dissipation.


## 7. Conclusion

In this study, the impact force and penetration depth of rockfall were analyzed theoretically. The main conclusions are as follows:(1)
This investigation began with an analysis of the mechanical response mechanism of rockfall impact on the soil cushion layer. Based on the theoretical derivation and the combined cavity expansion model and energy conservation model, respectively, two analytical solutions for the penetration depth and impact forces were obtained. Further, the LS-DYNA software was employed to study the variation law of penetration depth and impact force.(2)
For the impact force of the rockfall, either with the change of impact velocity or its radius, the calculation results for the cavity expansion, energy conservation, and numerical simulation methods demonstrated better consistency and coincided well with the results for the Japanese, Swiss, and Australian methods. The results of the Japanese method and the Swiss method were obtained based on laboratory tests. The results for the Australian method were acquired based on field tests.(3)
For the penetration depth and contact radius of the rockfall either with the change of impact velocity or its radius, the calculation results for the cavity expansion, energy conservation, numerical simulation, and B. S. Guan methods exhibited better consistency and coincided well with the results for the Australia method.(4)
For the relation between the “penetration depth-impact force” and the “contact radius-impact force,” the calculation results for the cavity expansion, energy conservation, numerical simulation, and Australian methods exhibited better consistency. The penetration depth obtained using the Hertz Method was smaller, and the value for the impact force was greater. In comparison, the tunnel manual method and ubgrade method showed inverse results. The penetration depth obtained using the B. S. Guan method is accurate, but the impact force is slightly smaller. Thus, this indicates that the cavity expansion method and energy conservation method have better accuracy.(5)
When the impact velocity or radius of the rockfall increases, both the deformation energy consumption and friction energy consumption increase in the form of a parabola. Further, the friction energy consumption is a bit larger than the deformation energy consumption. The impact velocity of the rockfall will influence the energy dissipation ratio for various energies, while the radius of the rockfall has a slight influence. A normal shed hole structure rarely influences the energy dissipation process.



---

*Source: 9324956-2018-12-26.xml*