Numerical modeling of the rockslide break-runout-impact on the flexible barrier using DEM-FDM coupled method

Abstract

Rockslides can cause severe damage to infrastructure, while flexible barriers offer an effective prevention measure. The protection and control effect of flexible barriers can be estimated from the process of rockslide break-runout-impact on the flexible barrier, which can be simulated using the DEM-FDM coupled method in the finite difference code FLAC3D. However, accurate modeling of this process has been hindered by the inability to adequately describe the mechanical behavior of flexible barriers in FLAC3D. To address this issue, this study presents an equivalent continuum model to describe the mechanical behavior of a flexible barrier under rockslide impact. This equivalent model reduces computational complexity while maintaining accuracy. The proposed theoretical model incorporates tensile and bending components along with yield criteria governing stress state changes. Furthermore, the model was implemented in FLAC3D through secondary development. Three experiments were numerically simulated to verify the effectiveness of the proposed model. Utilizing the proposed flexible barrier model and the DEM-FDM coupled method, the process of rockslide break-runout-impact on flexible barriers was numerically investigated. The proposed method was applied to predict the rockslide break-runout-impact process on a flexible barrier in Sichuan Province, China. The numerical results offer scientific guidance for the selection of appropriate flexible barriers.

Keywords

rockslide flexible barrier numerical modeling DEM-FDM coupled secondary development

Introduction

Rockslides can cause destructive impacts on infrastructure, which require adequate protection and control measures (Chen et al., 2023; Cui et al., 2023; Scuderi et al., 2023; Shao et al., 2023; Tomás et al., 2023; Verrucci et al., 2023). Flexible barriers have been widely adopted for rockslide protection and control due to their cost-effectiveness and environmental friendliness (Coulibaly et al., 2019; Tang et al., 2024; Xu et al., 2018; Yang et al., 2019, 2024), as shown in Figure 1. The efficiency of flexible barriers significantly depends on the rockslide break-runout-impact process. For instance, the breaking behavior relates to the risk level of potential rockslides, providing a basis for implementing protective measures. (Li et al., 2023). The runout trajectory governs the kinetic energy and jump height of the rockslide mass, which are critical parameters for designing hazard-mitigation structures (Asteriou and Tsiambaos, 2018). Predicting barrier effectiveness based on their dynamic response is highly complex and challenging (Yu et al., 2019). Therefore, accurately simulating the rockslide break-runout-impact process on flexible barriers is essential for their effective design.

Figure 1.

Flexible barriers intercepting rocks: (a) field tests, (b) details of flexible barriers.

Even though the process of rockslide break-runout-impact on flexible barriers can provide reliable suggestions during protection design, it has not yet been accurately modeled. The accurate numerical modeling of rockslide break-runout-impact process remains a significant challenge due to several factors, including: the complex mechanical behavior and large deformations of the flexible barrier, the dynamic of interaction between the discrete rock mass and the continuum structure, the need for efficient computation of large domains, the accurate numerical modeling of slope, and the fragmentation of rocks during the fall, etc. Researchers have separately conducted numerical modeling of the rockslide break-runout process and its impact on flexible barriers. For example, the DEM-FDM coupled method effectively simulated the rockslide break-runout process. The finite difference code FLAC3D is a common platform for numerically modeling the rockslide break-runout process using the DEM-FDM coupled method. (di Prisco et al., 2023; Shi et al., 2021; Wang et al., 2022). Similarly, FEM (Escallón et al., 2015; Gentilini et al., 2012; Qi et al., 2018), DEM (Albaba et al., 2017), Coupled Arbitrary Lagrangian Eulerian and Finite Element Method (ALE-FEM) (Kwan et al., 2016), and DEM-FEM (Leonardi et al., 2016), have been widely adopted to numerically model the rockslide impact on flexible barriers. However, few studies have simultaneously modeled both the rockslide break-runout process and its impact on flexible barriers. For example, Zhao et al. (2020) investigated the complex dynamic interaction between a debris avalanche and flexible barrier using a coupled DEM-FEM. This model does not consider the actual 3D terrain, potentially leading to discrepancies in the actual deformation and runout behavior. Each method has its own advantages and limitations. FEM is well-suited for modeling structural responses with high accuracy but often struggles with simulating extreme deformations, fragmentation, and discontinuous processes. In contrast, DEM excels in capturing granular flow, block interactions, large displacements, and fragmentation phenomena, but is computationally expensive for large-domain or long-duration simulations. The coupled DEM-FDM approach, adopted in this study, effectively integrates the strengths of both methods: DEM efficiently simulates the complex behavior of the rockslide mass, including disintegration, runout, and impact, while FDM accurately represents the complex 3D bedrock topography and facilitates computationally efficient analysis of structural responses. However, a major challenge within the DEM-FDM coupled method lies in the representation of the flexible barrier itself. Modeling the barrier with DEM would require an excessive number of particles, resulting in high computational cost. Although representing the barrier within FDM is computationally efficient, existing constitutive models in codes such as FLAC3D are unable to accurately capture the large deformations, membrane behavior, and progressive failure typical of flexible barriers under impact.

This limitation highlights the critical need for a novel constitutive model that is both efficient and mechanically appropriate for flexible barriers. Although the model in this study is implemented in FLAC3D—due to its computational efficiency and precision, the proposed equivalent continuum approach is general and could be adapted to other continuum-based numerical platforms. To this end, this study theoretically developed a transversely isotropic plastic hardening constitutive model for flexible barriers and implemented it in FLAC3D through secondary development. The model was validated using experimental data from the literature and subsequently applied to analyze a case study in Jiuzhaigou County, China, examining the influence of protective energy capacity and barrier placement on protection effectiveness.

This rest of this paper is organized as follows. “Methodology” describes the processing of proposed theoretical model for flexible barriers. “Numerical example for validation” verifies the proposed method. “Application case” presents an application case. “Discussion and conclusion” draws the discussion and conclusion in this study.

Methodology

The DEM-FDM coupled method was employed to simulate the rockslide break-runout-impact process on flexible barriers. As illustrated in Figure 2, the numerical model was divided into two parts: the flexible barrier and bedrock were modeled using FDM continuum zones, while the rockslide was modeled using DEM discrete particles. To accurately model the dynamic response of flexible barriers under impact, a theoretical model was developed. The research steps are as follows: (1) theoretical basis and assumptions; (2) theoretical analysis of the tensile ring; (3) implementation in FLAC3D via secondary development. The process of establishing the rockslide and bedrock models involves building the numerical geometry and determining the material parameters. Ultimately, the validated numerical model was successfully applied to a practical case study of an actual rockslide event.

Figure 2.

Schematic flowchart of this study.

Determination of the equivalent model

Precise and equivalent numerical models were adopted to analyze the mechanical behavior of the flexible barriers. Precise numerical models explicitly represent individual components of flexible barriers, such as energy-dissipating devices (Castanon et al., 2018). However, this explicit representation often incurs significant computational overhead and cost. To address this issue, researchers have adopted equivalent numerical models to represent flexible barriers. Equivalent models simplify the complex flexible barrier structure such that the simplified representation maintains consistency with the mechanical properties, deformation behavior, and dynamic response of actual flexible barriers. For example, Tahmasbi et al. (2019) described a hybrid model of a flexible barrier in which the area impacted by rocks was designated as a beam element and the surrounding region was represented as a shell element. Other studies (e.g., Sautter et al. (2021), Mentani et al. (2018), and Dhakal et al. (2011)) have also presented a complex protection structure with a simpler surface description composed of isotropic membrane elements.

Mentani et al. (2018) confirmed that using mechanical test results of the metal rings in flexible barriers provides reasonable and reliable data for developing equivalent models. As shown in Figure 3, flexible barriers can be conceptualized as a “film” material, with the rings serving as the fundamental units constituting the “skeleton”.

Figure 3.

Equivalent model for flexible barriers.

DEM-FDM coupled method

The DEM-FDM coupled method was used to simulate the actual trajectory of rockslides. Figure 4 shows the principle of the DEM-FDM coupled method (Cai et al., 2007). In FDM, force and velocity are transmitted through nodes, which can result in large deformations that are challenging to simulate accurately. Meanwhile, DEM excels at simulating large displacements, rotations, and complete separations, but typically incurs high computational costs. Therefore, the DEM-FDM coupled method leverages the strengths of each approach: bedrock and flexible barriers are modeled using FLAC3D continuum elements (FDM), while the rockslide is modeled using PFC3D discrete particles (DEM). This hybrid approach offers computational efficiency and numerical convenience.

Figure 4.

Principle of DEM-FDM coupled method.

Development of theoretical and numerical model for flexible barriers

Theoretical basis and assumptions

The following theoretical analysis of the ring tensile behavior is based on the experimental work of Nicot and Cambou (2001) and employs a rigid-plastic analytical framework to derive the transition between deformation stages. The theoretical analysis is based on these assumptions:

(1) The ring material is idealized as rigid-perfectly plastic, neglecting elastic deformations.

(2) The cross-sectional diameter d is significantly smaller than the ring radius R (d<<R), allowing it to be neglected in the geometric analysis of the neutral axis.

(3) The total length of the ring remains constant throughout the deformation process.

(4) The collapse mechanism is characterized by the formation of four discrete plastic hinges.

(5) Segments of the ring between plastic hinges undergo rigid body motion.

(6) In the collapse mechanism (Figure 7(b)), $\overset{⌢}{E F}$ is assumed to remain circular (radius R, center O) based on rigid rotation around O. FH is simplified as a straight line, a common simplification in plastic hinge theory that captures the fundamental mechanics while acknowledging it may introduce minor error in the precise deformation path.

(7) The tensile response of the ring is idealized into two distinct mechanical phases: a pure bending-dominated stage where the influence of axial force is neglected, followed by a pure axial-dominated stage where the influence of bending moment is neglected.

Theoretical analysis of the tensile ring

Descriptions of experimental behavior: The mechanical behavior of the tensile ring, which is fundamental to understanding the response of flexible barriers under impact, was investigated experimentally by Nicot and Cambou (2001). As shown in Figure 5(a), a single ring was subjected to diametral tension, and the resulting force-displacement curve (Figure 5(b)) exhibits two distinct deformation stages. In the first stage (OAB domain), the ring undergoes considerable deformation with only a slow increase in tensile force. This response is primarily governed by bending deformation, as the ring forms plastic hinges and undergoes geometric rearrangement without significant axial stretching (Figure 6). In the second stage (BC domain), the force increases rapidly with limited additional deformation. Here, the ring behaves essentially as a tensile member, with resistance dominated by axial action similar to that of a low-carbon steel bar under tension. Therefore, the ring response can be divided into two mechanistic stages: an initial bending-dominated stage, where the bending moment governs the behavior entirely and the axial force has no effect; and a subsequent hardening stage, where the behavior is entirely dominated by the axial force, and the bending moment has no influence. The objective of the following analysis is to theoretically derive the transition point between these two stages and provide a mechanical interpretation for this critical behavior.

Figure 5.

Tensile ring test: (a) tensile loading test on a single ring, (b) experimental behavior of the ring.

Figure 6.

Deformation process of flexible barrier metal rings.

Theoretical derivation of the transition point: An ideal rigid-plastic ring with a neutral axis radius R and cross-sectional thickness d is subjected to a pair of radial tensile forces, as shown in Figure 7(a). Given that d is significantly smaller than R, the cross-sectional thickness d can be neglected. The strength of the metal ring is proportional to d, allowing different wire rope diameters to be represented by varying d. The radial tensile ring requires four plastic hinges to form a collapse mechanism, where F and H are moving plastic hinges. E is the fixed plastic hinge and FH represents the straightened segment, as shown in Figure 7(b). To render the problem analytically tractable, two key simplifications are introduced: (1) EF is treated as a circular arc of radius R, centered at point O, consistent with the assumption of rigid rotation of segments OE and OF; (2) FH is idealized as a straight line, a conventional simplification in plastic mechanism analysis that assumes deformation localizes at discrete hinges with rigid segments in between. These simplifications are justified as they capture the essential kinematics of collapse observed experimentally (Figure 6). The primary error introduced is a slight overestimation of the displacement at the transition point, as the model neglects distributed plasticity and elastic effects that would soften the initial response. However, as will be shown, the model’s prediction aligns well with experimental results, validating its use for defining the critical transition criterion.

Figure 7.

Rigid plastic ring under tension: (a) geometric model, (b) collapse mechanism.

Let $φ$ denote the deflection angle of the neutral axis. A point O is defined inside the ring such that OF = OE = R, and θ is the angle between OE and the perpendicular direction. The initial stage of deformation is described as follows: Initially, the moving plastic hinge F and the fixed plastic hinge E converge. Meanwhile, segment FH continues to enlarge and segment FG gradually decreases. The moving plastic hinge F and the fixed plastic hinge E overlap entirely. Ultimately, the collapse mechanism evolves into that of two parallel wire ropes. The total length of the ring was assumed to remain constant throughout the deformation process. Therefore, the geometric relationship within the collapsed structure can be expressed by equation (1):

{\begin{cases} \begin{array}{l} O^{'} E = R + \frac{δ}{2} \\ E G = R c o s θ \end{array} \\ F G = R - R s i n θ \\ ∠ E O F = \frac{π}{2} - θ \end{cases}

(1)

where

δ

denotes the elongation of the ring, the length of

\hat{E F}

R (\frac{π}{2} - θ)

; and the total length of the ring is calculated using equation (2):

4 E F + 4 O^{'} G = 4 R (\frac{π}{2} - θ) + 4 O^{'} G = 2 π R

(2)

It can be obtained from equation (2):

O^{'} G = θ R

(3)

Due to $E G + O^{'} G = O^{'} E$ , therefore:

R c o s θ + θ R = R + \frac{δ}{2}

(4)

The relationship between the displacement and radius of the ring is obtained as follows:

δ = 2 R (\cos θ + θ - 1)

(5)

or:

\frac{δ}{D} = \frac{δ}{2 R} = \cos θ + θ - 1

(6)

The transition from the initial bending-dominated stage to the subsequent axial-force-dominated hardening stage is ultimately governed by material failure due to excessive bending strain at the plastic hinges, not merely by geometric evolution. Figure 8 shows a partially enlarged view of the tensile ring. The most significant curvature change was in the position of the ring under tension. According to the large deformation theory, this point is a fixed plastic hinge, and the finite length of the plastic hinge is $λ = (2 \sim 5) d$ , where d denotes the diameter of the circular section of the steel wire and 3d is adopted here. The final curvature was $κ = 2 φ / 3 d$ , and the initial curvature of the outermost fiber of the initial ring was 1/R. Therefore, the maximum bending strain of the outermost fiber is:

ε = (κ - \frac{1}{R}) \frac{d}{2} = (\frac{2 φ}{3 d} - \frac{1}{R}) \frac{d}{2} \approx \frac{φ}{3}

(7)

Figure 8.

Partially enlarged view of the plastic hinge region.

Let $ε_{f}$ denote the maximum tensile strain of the tensile low-carbon steel and $ε_{f} \approx 0.3$ . The outermost fiber must be less than the maximum tensile strain of the material, that is $ε < ε_{f}$ , that is $φ < 52 °$ . In addition, the geometrical positions of $θ$ and $φ$ are similar during the tensile ring process. It is assumed that $θ + φ \to π / 2$ was within this limit. Therefore, the critical $θ = \frac{π}{2} - φ = 38 °$ can be calculated. The limit ratio of the displacement to diameter can be obtained as follows:

\frac{δ}{D} = \frac{δ}{2 R} = \cos θ + θ - 1 = \cos 38 ° + \frac{38 °}{180 °} π - 1 = 0.45

(8)

The tensile ring test of Zhao (2018) was used to verify the theoretical analysis, as shown in Figure 9(a)–(e). The results of the tensile load-displacement test are shown in Figure 10(a), where it can be observed that the critical point of the hardening stage was approximately 130 mm corresponding to $\frac{δ}{D} = \frac{130}{300} = 0.43$ , which is almost entirely consistent with the theoretical value proposed in this section.

Figure 9.

Deformation process of tensile ring: (a) start, (b) bending stage, (c) tensile stage, (d) break stage, (e) location of the breakpoint.

Figure 10.

Results of tensile ring tests: (a) force-displacement curve of tensile ring; (b) stress-strain curve of the equivalent model.

Table 1 lists the tensile breaking forces of the three ring types. The ultimate tensile stress of the ring is

σ = \frac{F}{n A}

, where F is the ring tensile breaking force, A is the ring cross-sectional area, and n is the number of steel rings. The ultimate tensile stress of the rings is equivalent to the nominal strength of the steel wire rope (1770 MPa). The computed ultimate stresses (1542 MPa, 1688 MPa, and 1788 MPa) align closely with the nominal strength of the steel wire rope (1770 MPa), indicating that the mechanical behavior of the rings during the hardening stage is consistent with that of the wire rope. This agreement supports the two-stage idealization adopted in the model: the initial stage is dominated purely by bending, with negligible contribution from axial force, while the hardening stage is governed primarily by axial stretching. The initial and hardening equivalent modulus correspond to E₁ and E₂, respectively, as illustrated in Figure 10(b). The initial equivalent modulus (E₁) is consistent across all ring types due to the bending-dominated response in the first stage, which is mainly influenced by the ring geometry and is independent of the number of rings. This is further corroborated by the nearly identical initial slopes observed in the experimental force–displacement curves. Failure is considered to occur when the maximum principal stress under combined bending and tension exceeds the ultimate tensile strength [σ₀] of the low-carbon steel.

Table 1.

Experimental results of the ring under tension.

Items	Type
Items	R5/3/300	R7/3/300	R9/3/300
Break tensile force (kN)	54.50	83.50	113.70
Ultimate tensile stress (MPa)	1542.00	1688.00	1788.00
Initial equivalent modulus (kPa)	40.00	40.00	40.00
Hardening equivalent modulus (MPa)	3.00	4.00	5.50

Theoretical derivation of constitute model: The constitutive model developed here describes the stress-strain relationship of the equivalent flexible barrier material, which is divided into the initial and hardening stages. The constitutive equation of the equivalent film material in the initial stage is:

{\begin{cases} ε = \frac{σ}{E_{e 1}} σ \leq σ_{e 1} \\ ε = \frac{σ_{e 1}}{E_{e 1}} + \frac{1}{E_{p 1}} (σ - σ_{e 1}) σ_{e 1} < σ \leq σ_{s} \end{cases}

(9)

where E_e1 and E_p1 are the elastic and plastic moduli, respectively, for the equivalent film material in the initial stage. Due to the relatively small value of

E_{e}

E_{e} = E_{p}

. Therefore, the increment of the stress-strain of the equivalent film material in the initial stage is:

d ε_{i j}^{0} = C_{i j k l}^{e} d σ_{k l}

(10)

where

d ε_{i j}^{0}

is the strain increment in the initial stage of the equivalent film;

C_{i j k l}^{e}

is the elastic compliance matrix of transversely isotropic material, which is expressed as follows:

C_{i j k l}^{e} = [\begin{array}{l} \frac{1}{E} & - \frac{v}{E} & - \frac{v^{'}}{E^{'}} & 0 & 0 & 0 \\ - \frac{v}{E} & \frac{1}{E} & - \frac{v^{'}}{E^{'}} & 0 & 0 & 0 \\ - \frac{v^{'}}{E^{'}} & - \frac{v^{'}}{E^{'}} & \frac{1}{E^{'}} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{G} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G^{'}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{G^{'}} \end{array}]

(11)

where E and E′ are the elastic modulus of the homogenous surface and the anisotropic surface, respectively; v and v' are the Poisson’s ratio of the homogenous surface and the anisotropic surface, respectively; G and G′ are the shear modulus of the homogenous surface and the anisotropic surface, respectively.

The constitutive model for the equivalent film material is consistent with the ring material in the hardening stage. Based on the Hill yield criterion for anisotropic materials, the corrected yield function of transversely isotropic materials is:

f = \sqrt{a {(σ_{2} - σ_{3})}^{2} + a {(σ_{3} - σ_{1})}^{2} + b {(σ_{1} - σ_{2})}^{2} + c τ_{23}^{2} + d τ_{31}^{2} + e τ_{12}^{2}} - ε_{p}

(12)

where a, b, c, d, e are material constants, and

ε_{p}

is the yield stress, which evolves with

ε_{p}

. The effective plastic strain

ε_{p}

is defined as:

d ε_{p} = \sqrt{d ε_{i j}^{p} d ε_{k l}^{p}}

(13)

The correlation flow rule is adopted for the metal material constituting the ring. The relationship between the plastic strain increment and the yield function is as follows:

d ε_{i j}^{p} = d λ \frac{\partial f}{\partial σ_{k l}}

(14)

where

d λ

is the plastic multiplier, which is determined by the consistency condition

d f = 0

d f = \frac{\partial f}{\partial σ_{i j}} d σ_{k l} + \frac{\partial f}{\partial k} d k = 0

(15)

Combining equations (13)–(16) yields:

d λ = - \frac{\frac{\partial f}{\partial σ_{k l}} d σ_{k l}}{\frac{\partial f}{\partial k} \sqrt{\frac{\partial f}{\partial σ_{i j}} \frac{\partial f}{\partial σ_{k l}}}}

(16)

Numerical parameters

The transverse isotropic plastic hardening constitutive model has 10 numerical parameters, as listed in Table 2. Three state indicators are set to describe the flexible barrier stages, as listed in Table 3. In the equivalent numerical model, the strain component within an element is utilized to simulate the displacement (δ) of the represented ring, while the characteristic size of the element (e.g., length) corresponds to the ring diameter (d). Thus, the zone enters the hardening state when the strain reaches 0.45 times the size of the zone, and the zone enters the breaking state when its maximum principal stress exceeds the breaking force in the ring test.

Table 2.

Numerical parameters of the constitutive model.

Items	Keywords	Descriptions
Initial stage	E-initial-plane	Equivalent deformation modulus in isotropy plane
	E-initial-normal	Equivalent deformation modulus in heterogeneous plane
	Poisson-plane	Equivalent Poisson’s ratio in the plane of isotropy
	Poisson-normal	Equivalent Poisson’s ratio normal to isotropy plane
	Shear-normal	Shear modulus for any plane normal to isotropy plane
Isotropic plane	Dip	Dip angle of the plane of isotropy
Isotropic plane	Dip-direction	Dip direction of the plane of isotropy
Hardening stage	E-harden-plane	Equivalent hardening deformation modulus in isotropy plane
	E-harden-normal	Equivalent hardening deformation modulus in heterogeneous plane
	Poisson-harden	Equivalent hardening Poisson’s ratio
	Break-force	Break force of the ring

Note. The tensile ring test only occurs on homogeneous surfaces, and the parameters of the heterogeneous surfaces are lacking. The authors selected the following appropriate simulation parameters through repeated experiments: Poisson-normal = 0. E-initial-normal = Shear-normal $\approx$ (10-100) $\times$ E-initial-plane. Poisson-plane $\approx$ 0.3∼0.5. Poisson-harden $\approx$ 1/3∼1/2 $\times$ (Poisson-plane).

Table 3.

State indicators of the constitute model.

State indicators	Descriptions
None	Zones in this state are at the initial stage
Harden	Zones in this state are at the hardening stage
Break	Zones in this state are at the broken stage and the actual ring breaks

Parameter calibration and generalization

The proposed constitutive model provides a framework for simulating flexible barriers with varying ring configurations (number of rings, diameter, etc.) and material properties, which requires calibration of key parameters (E-initial-plane, E-harden-plane, Break-force) for specific ring designs. The methodology is as follows:

Experimental calibration: The most reliable approach involves conducting a tensile test on a representative ring or ring assembly. The force-displacement (F-x) curve obtained is used to directly calibrate the parameters. The initial slope is used to determine the initial modulus E₁ (bending-dominated global stiffness). The slope in the hardening stage gives the hardening modulus E₂, and the peak force gives the Break-force.

Theoretical guidance: In the absence of immediate test data, theoretical derivations offer initial estimates. The transition strain between stages is consistently estimated at δ/D ≈ 0.45, as derived in Section 2.3.2. The initial modulus (E₁) can be approximated based on ring diameter and bending theory, as it is primarily governed by geometry. The hardening modulus (E₂) can be initially estimated from the elastic modulus of the wire material.

Numerical inverse analysis: Parameters can be further refined through iterative numerical simulation of the component test, minimizing the difference between the simulated and experimental F-x curves.

Considerations for different design variables: (1) Number of rings (n): the number of rings primarily influences the Break-force and the hardening modulus (E₂), as they scale with the total load-bearing cross-sectional area. The initial modulus (E₁), being dominated by bending geometry, remains largely independent of the number of rings. (2) Ring diameter (D): The ring diameter is a key factor influencing the initial modulus (E₁), as bending stiffness is highly sensitive to geometric size. (3) Wire material properties: different material properties (e.g., yield strength, hardening modulus) fundamentally alter the force-displacement response. While the theoretical transition strain (δ/D ≈ 0.45) is geometry-based and thus less sensitive to material changes, accurate determination of the Break-force and hardening modulus for a new material requires experimental data or material-specific constitutive properties.

Numerical example for validation

Three experimental tests were simulated to verify the proposed constitutive model for the flexible barrier: a single ring under tension (Zhao 2018), multi-ring under tension (Boulaud and Douthe, 2022) and full-scale impact tests (Gentilini et al., 2012; Gottardi and Govoni, 2010).

Single ring under tension

As shown in Figure 11, the size of the established numerical model was 0.3 × 0.03 × 0.3 m. The dimension in the y direction was 0.1 times that of the x and z directions, as the length in the y direction remains unaffected. The upper and lower boundaries were under tension along the diameter with a loading velocity of 8 × 10^-5 m/s. The numerical parameters are listed in Table 4. A detailed comparison of the numerical and experimental force–displacement curves (Figure 11(b)) reveals differences in the precise loading path, particularly in the initial and post-transition stages. To quantitatively evaluate the accuracy of the simulation, key mechanical indicators are compared in Table 5. The numerical model shows agreement with experimental results in terms of the ultimate breaking force (9.4% error) and the displacement at the transition point (3.8% error), demonstrating its capability to accurately predict critical mechanical properties such as load-bearing capacity and the onset of the hardening stage. However, a noticeable discrepancy is observed in the initial stiffness, with a relative error of 44.5%. This deviation can be primarily attributed to the idealized boundary conditions (e.g., frictionless pinned connections) used in the simulation, which differ from the actual physical constraints (e.g., friction and minor gaps in bolted connections) in the experiment. Additionally, the neglect of elastic deformation in the adopted rigid-plastic material model further contributes to this discrepancy in the initial elastic response. Despite these minor discrepancies in the detailed loading path, the overall trend of the simulated force–displacement curve agrees well with experimental measurements. The model reproduces the characteristic two-stage mechanical behavior of the ring—namely, the initial bending-dominated stage with low stiffness, followed by the transition to a high-stiffness, axial-dominated hardening stage. The close agreement in ultimate breaking force and transition displacement, as quantitatively summarized in Table 5, confirms the validity of the proposed model for assessing the overall performance of the ring system.

Figure 11.

The single ring under tension: (a) experiment test, (b) comparison of the numerical modeling and experimental test (the black and red curve correspond to the black and red axis respectively), (c) numerical modeling, (d) deformation states of the ring.

Table 4.

Numerical parameters of the single ring under tension.

Items	Values
Dip (°)	90
Dip-direction (°)	0
E-initial-plane (Pa)	4 × 10⁴
E-initial-normal (Pa)	3 × 10⁶
Poisson-plane	0.45
Poisson-normal	0
Shear-normal (Pa)	3 × 10⁶
E-harden-plane (Pa)	3 × 10⁶
Poisson-harden	0.20
Break-force (N)	5.5 × 10⁴

Table 5.

Comparison of key mechanical indicators between simulation and experimental results for the single-ring test.

Mechanical indicator	Simulation	Experiment	Relative error
Initial stiffness (N/mm)	47.1	32.6	44.5 %
Transition point (mm)	135.0	130.0	3.8 %
Ultimate breaking force (kN)	59.6	54.5	9.4 %

Multi-ring under tension

As shown in Figure 12, the multi-ring under tension consisted of a square rigid support frame with 3 × 3 rings, each with a radius of 13.5 cm. During loading, one edge of the support frame was displaced while constraining normal displacements on the opposite edge; tangential displacements were free on all boundaries. Table 6 lists the numerical parameters used in this study. A quantitative comparison of key mechanical indicators between the simulation and experimental results is summarized in Table 7. The relative errors for the ultimate breaking force and transition displacement are 19.2% and 17.3%, respectively. Although these values are higher than those observed in the single-ring test, they are considered acceptable given the greater complexity of the multi-ring system, which involves intricate ring-to-ring and ring-to-frame interactions. The most significant discrepancy is observed in the initial stiffness (37.9% error), which is consistent with the single-ring findings and can be similarly attributed to the idealized boundary conditions and material model. As shown in Figure 12(b), the transition to the hardening stage occurred at a displacement of approximately 0.15 m, which aligns broadly with the value of 0.12 m predicted by the single-ring theoretical model (Equation (8)). The minor deviation can be attributed to the idealized assumptions of the analytical model, particularly the neglect of elastic deformation and detailed contact interactions. Nevertheless, the overall mechanical behavior shows good agreement with the theoretical prediction. Despite these discrepancies in quantitative values, the numerical results reproduce the key features of the experimental response, including the characteristic two-stage behavior and the ultimate failure mode. As shown in Figure 12(d), the simulation accurately captures the progressive failure sequence, with edge rings entering the hardening state first, followed by the central ring—a pattern also observed experimentally. This consistency between simulation and experiment, particularly in capturing the essential failure mechanisms, demonstrates the ability of the proposed constitutive model to represent the essential mechanical behavior of multi-ring barrier systems.

Figure 12.

Multi-ring under tension: (a) experiment test, (b) comparison of the numerical modeling and experimental tests (the black and red curve correspond to the black and red axis respectively), (c) numerical modeling, (d) deformation states of multi-rings.

Table 6.

Numerical parameters of the multi-ring under tension.

Items	Values
Dip (°)	90
Dip-direction (°)	0
E-initial-plane (Pa)	4 × 10⁴
E-initial-normal (Pa)	3 × 10⁶
Poisson-plane	0.45
Poisson-normal	0
Shear-normal (Pa)	3 × 10⁶
E-harden-plane (Pa)	1 × 10⁶
Poisson-harden	0.20
Break-force (N)	1 × 10⁵

Table 7.

Comparison of key mechanical indicators between simulation and experimental results for the multi-ring test.

Mechanical indicator	Simulation	Experiment	Relative error
Initial stiffness (N/mm)	50.2	80.9	37.9 %
Transition point (mm)	121.5	147	17.3 %
Ultimate breaking force (kN)	84	104	19.2 %

Full-scale impact tests

To further validate the proposed constitutive model and the coupled DEM-FDM approach, a simulation of a full-scale impact test was conducted and compared with field test data reported in the literature (Gentilini et al., 2012; Gottardi and Govoni, 2010) in Fonzaso, Belluno, Italy. The numerical model was configured to replicate the conditions of the referenced field test. The flexible barrier parameters were calibrated based on tensile test data of the OM CTR 50/07/A 5000 kJ flexible barrier (Gottardi and Govoni, 2010), as shown in Figure 13. The mesh size in the numerical model was set to 0.35 m to match the actual ring size of the physical barrier. The dotted line represents the experimental results and the solid line represents the fitted constitutive model curve of the flexible barrier. Therefore, the initial equivalent modulus (E₁) and hardening modulus (E₂) were set to 3.5 MPa and 35 MPa, respectively. Table 8 summarizes the complete set of numerical parameters used for the flexible barrier model. The isotropic and heterogeneous planes of the flexible barrier deflected during rockslide impact; therefore, the numerical parameters were simplified. Owing to the lack of tensile tests for other types of flexible barriers, these numerical parameters were estimated based on OM CTR 50/07/A. The initial modulus (E₁) is designed to be 0.1 times the hardening modulus (E₂). The modulus of the flexible barriers and the variation in local damping are generally positively correlated with the impact energy.

Figure 13.

Flexible barrier (OM CTR 50/07/A, 5000 kJ) tensile test: (a) experiment test; (b) comparison of the numerical modeling and experimental tests.

Table 8.

Numerical parameters of the flexible barrier.

Energy level (kJ)	500	3000	5000
E-initial (MPa)	1.50	4.00	3.50
Poisson-plane	0.20	0.20	0.20
E-harden (MPa)	15.00	40.00	35.00
Poisson-harden	0	0	0
Tension (N)	8 × 10⁵	3 × 10⁶	4 × 10⁶
Local damping	0.68	0.75	0.78

Figure 14 presents a comparison between the field test results and the numerical simulation. The model reproduces the key characteristics of the impact process, including the dynamic response of the barrier and the resulting deformation pattern. The agreement is observed between the simulated and experimental deformation values across all energy levels. For the 500 kJ impact case, both the maximum displacement (2.9 m) and final displacement (2.6 m) match exactly between simulation and experiment. Similarly, for the 3000 kJ case, both maximum (5.5 m) and final (4.9 m) displacements show perfect agreement. This exact match is also maintained in the 5000 kJ case, with both maximum (5.1 m) and final (5.0 m) displacements identical between simulation and experiment. Furthermore, the simulation accurately captures the slight rebound phenomenon of the flexible barrier observed in physical tests. The failure progression sequence observed in the simulation matches exactly with the actual field test recordings. The structural behavior reveals that the area of rings entering the hardening stage expands progressively with increasing impact energy. Under low-energy conditions, energy dissipation occurs primarily through large deformation and bending moment resistance, whereas high-energy impacts are mainly absorbed through axial force-bearing during the hardening stage. The simulation reproduces the actual propagation sequence of hardening: beginning at the impact point, extending to both ends of the flexible barrier, and finally connecting the central and end regions. This progression, consistently observed in both field tests and simulations, confirms that the most vulnerable regions in practical applications are, in order, the impact location, the fixed ends, and the diagonal areas of the flexible barrier. The concurrence between simulated and experimental failure sequences further validates the selected parameters and overall modeling methodology.

Figure 14.

Field test compared with numerical simulation.

Application case

Building the numerical model

The “8.8″ Jiuzhaigou earthquake induced the formation of several perilous rock masses (Hu et al., 2019; Ling et al., 2021). The Huohua Lake perilous rock is located in the Jiuzhaigou Valley Scenic, Sichuan Province (33°13′3.20″N, 103°54′16.99″E) (Figure 15). Currently, the perilous rock is approximately 5.86 × 10⁴ m³ and can easily induce a rockslide under the action of gravity, rainfall, etc. Figure 16 shows the rockslide DEM-FDM coupled numerical model.

Figure 15.

Location of the study area.

Figure 16.

Numerical modeling of rockslides break-runout-impact on flexible barriers: (a) bird view, (b) details of flexible barriers, (c) vertical view.

Determining the numerical parameters

The fragmentation of the rockslide mass during the falling and impact process is a critical mechanism influencing its runout behavior and impact force. In this study, the fragmentation is simulated within the DEM framework (PFC3D) using the Linear Parallel Bond Model (LPBM). This model represents the intact rock mass as an assembly of discrete particles bonded together at their contacts. The breakage of the rock is governed by the following criteria at each bond:

(1) Tensile Failure: A bond breaks if the tensile stress acting on it exceeds its prescribed normal strength (σ_max).

(2) Shear Failure: A bond breaks if the shear stress acting on it exceeds its prescribed shear strength (τ_max).

Once a bond breaks, it is permanently removed, allowing new contacts to form between the now-separate fragments. This process enables the simulation of progressive rock fragmentation upon collision with the terrain or other rocks. The micro-parameters controlling this behavior, including the bond strengths, are summarized in Table 9. These values were selected based on calibration against typical rock properties and values from the literature (Li et al., 2021; Liu and Liu, 2022) to realistically represent the mechanical behavior of the jointed rock mass in the study area.

Table 9.

The numerical modeling of rockslide (DEM) and the bedrock (FDM).

Rockslide		Bedrock
Items	Values	Items	Values
Range of particle radius, R_min (m)	1.00∼2.00	Density, $ρ$ (kg/m³)	2650.00
Particle density, $ρ_{b}$ (kg/m³)	2500.00	Young’s modulus, E (GPa)	28.00
Effective modulus, emod (MPa)	25.00	Poisson’s ratio, v	0.25
Normal-to-shear stiffness ratio, kratio	3.10	Cohesion, c (MPa)	2.50
Normal strength of parallel bonds, (MPa)	20.00	Internal angle of friction, $φ$ (°)	45.00
Shear strength of parallel bonds, (MPa)	40.00
Friction coefficient, fric	0.60
Viscous damping ratio of normal direction $β_{n}$	0.50
Viscous damping ratio of shear direction $β_{s}$	0.50

Table 7 lists the complete set of numerical parameters for the rockslide (DEM) and bedrock (FDM) used for the method verification. The parameters for the flexible barrier were adopted from the validated model presented in Table 8, which was calibrated and verified against full-scale impact test data. This ensures that the barrier component of the coupled model accurately represents the mechanical behavior of real-world protective structures.

Rockslide break-runout-impact process

Figure 17 shows snapshots of the rockslide break-runout-impact process on a 500 kJ flexible barrier. (1) At 15.5 s, the particles first impacted the flexible barrier, and the impact position quickly entered the hardening state. (2) At 16.5 s, the flexible barrier began to enter the broken state as the number of particles reaching the flexible barrier increased. However, the flexible barrier still exhibited a robust load-bearing capacity owing to the small number of broken rings. (3) At 17.0 s, the number of broken flexible barrier rings gradually increased. Numerical modeling in this study did not explicitly simulate penetration mechanisms. Therefore, observed particle passage through the barrier mesh in the simulation may indicate potential penetration in reality. (4) At 18.5 s, a significant deformation of the flexible barrier occurred and many rings were in the broken state.

Figure 17.

Rockslides break-runout-impact flexible barriers process: (a) rockslides break, (b) rockslides runout, (c) rockslides impact on flexible barriers, (d) flexible barriers break.

Influence of the protective level

Figure 18 shows the numerical model of the rockslide break-runout-impact on flexible barriers of 500 kJ and 5000 kJ. Compared to the flexible barriers of 500 kJ, the zones of the flexible barriers of 3000 kJ and 5000 kJ are not in the broken stage, and both can protect and control rockslides. For safety and cost efficiency, the adoption of a flexible barrier of 3000 kJ is recommended.

Figure 18.

Influence of the protective level: (a) barrier with 3000 kJ, (b) barrier with 5000 kJ.

Influence of the location

The location of the flexible barrier is one of the main factors to be considered during the protection design. To this end, numerical models of three locations of the flexible barriers were conducted to investigate the influence of the location: flexible barriers near the rockslide (barrier Ⅰ), at the foot of the slope (barrier Ⅱ), and near the road (barrier Ⅲ). Barrier Ⅲ, the farthest distance allowed for installation, was constructed as described in the previous section. The selected flexible barrier energy was 3000 kJ, and all flexible barriers were vertical to the ground. As shown in Figure 19(a), the zones in the hardening state of barrier Ⅰ cover a much larger area than those of barrier Ⅲ, indicating that the barrier near the source requires lower energy capacity. As shown in Figure 19(b), the zones of barrier Ⅱ were in the broken stage. In addition, some rocks also bounce at the foot of the slope, exceeding the height of barrier Ⅱ. The numerical results indicate that barrier Ⅱ cannot effectively prevent and control the rockslides.

Figure 19.

Influence of the location: (a) near rockslides, (b) at the bottom of the slope.

Discussion and conclusion

This study has developed and validated an equivalent continuum model for simulating the entire process of rockslide break-runout-impact on flexible barriers using a DEM-FDM coupled approach. The discussion herein addresses the key findings, methodological contributions, and practical implications of this research.

Methodological advancements and validation

The proposed model bridges the gap between computationally expensive detailed simulations and oversimplified analytical models. Through theoretical derivation and experimental validation, we have established a clear mechanistic understanding of the ring’s two-stage behavior (transition from bending to axial stages), while the numerical implementation efficiently captures the system-level response of interconnected rings under complex loading conditions. The parameter calibration methodology presented in Section 2.4 enhances the model’s practical utility for different barrier designs, addressing a fundamental challenge in flexible barrier modeling.

The validation against both component-level tests (single and multi-ring tests) and full-scale impact demonstrations confirms the model’s robustness in predicting barrier performance across different scales. Particularly noteworthy is the model’s ability to reproduce the progressive failure sequence observed in multi-ring tests, where edge rings enter the hardening state before central rings, demonstrating its capacity to capture complex system-level behaviors.

Practical applications and engineering implications

The application case study demonstrates the model’s capability to inform critical engineering decisions. The analysis of energy levels (500 kJ, 3000 kJ, 5000 kJ) provides a quantitative basis for selecting cost-effective barriers, with the 3000 kJ barrier identified as optimal for the specific case study. The investigation of barrier placement reveals crucial design considerations: barriers near the source (Position I) experience reduced energy demands but may face practical installation challenges; barriers at the slope toe (Position II) risk being overtopped or failed due to high-energy impacts and rock bounce phenomena; while barriers at a protected distance (Position III) represent a balanced solution that contains the debris flow.

The model’s ability to simulate rockfall jumping over barriers represents a significant advancement over simplified methods, providing engineers with a more realistic assessment of barrier performance. This capability is particularly valuable for designing protection systems in complex terrain where rockfall trajectories may be unpredictable.

Parameterization strategy and generalization

Addressing the important question of parameter determination for different ring configurations (varying number of rings, diameters, and materials), our study proposes a systematic approach:

(1) For designs with available experimental data, direct calibration using tensile tests provides the most accurate parameter determination.

(2) For preliminary design without test data, the theoretical framework (particularly the δ/D = 0.45 transition criterion) offers reliable initial estimates.

(3) The initial modulus can be approximated based on ring diameter and bending theory, while the hardening modulus can be estimated from the material properties.

(4) Numerical inverse analysis can further refine parameters through iterative simulation.

This hierarchical parameterization strategy ensures the model’s applicability to a wide range of flexible barrier designs while maintaining physical meaningfulness.

Limitations and future research

While the numerical results show satisfactory agreement with experimental data, several theoretical limitations must be acknowledged. The current analytical derivation remains confined to a simplified single-ring impact model and lacks a comprehensive theoretical framework for more complex mechanical states, such as multi-ring tensile interactions and full-scale barrier impact responses. Although numerical simulations effectively capture these behaviors, a corresponding theoretical explanation for system-level performance has not yet been established.

Future work will prioritize the development of theoretical models capable of addressing complex impact scenarios. Efforts will focus on formulating mechanical representations for multi-ring systems under tensile and impact loading, establishing theoretical model for flexible barrier components, and deriving analytically the force-redistribution mechanisms among rings and subsystems. In parallel, modeling capabilities will be enhanced through the integration of detailed connection elements and support systems, improved elastoplastic material constitutive laws, and collaborative frameworks for full-scale physical testing. Further research should also extend the model’s application to more complex conditions including successive impacts and varying rock size distributions.

Conclusion

The DEM-FDM coupled method with the newly developed constitutive model offers a powerful and practical tool for simulating rockslide impacts on flexible barriers. It provides scientific guidance for optimizing barrier energy capacity, installation location, and overall system design, contributing to more resilient infrastructure protection. The model’s ability to accurately capture both component-level behavior and system-level response while maintaining computational efficiency represents a significant advancement in geohazard protection engineering. The methodology and insights gained from this study not only advance the numerical simulation of flexible barriers but also provide a framework for addressing similar challenges in other engineering applications involving complex structure-fluid-solid interactions.

Footnotes

ORCID iD

Ziwei Ge

Author contributions

Ziwei Ge contributed to the conceptualization, methodology, data collection, data analysis, and writing—original draft preparation, validation, and manuscript editing.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: A Project Supported by Scientific Research Fund of Zhejiang Provincial Education Department (Y202454601). This research was supported by the Joint Fund of Zhejiang Provincial Natural Science Foundation of China under Grant No. LGEY25E090019.

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Data Availability Statement

All data generated or analyzed during this study are included in the article.*

References

Albaba

Lambert

Kneib

, et al. (2017) DEM modeling of a flexible barrier impacted by a dry granular flow. Rock Mechanics and Rock Engineering 50(11): 3029–3048.

Asteriou

Tsiambaos

(2018) Effect of impact velocity, block mass and hardness on the coefficients of restitution for rockfall analysis. International Journal of Rock Mechanics and Mining Sciences 106: 41–50.

Boulaud

Douthe

(2022) A comparative assessment of ASM4 rockfall barrier modelling. Engineering Structures 251: 113512.

Cai

Kaiser

Morioka

, et al. (2007) FLAC/PFC coupled numerical simulation of AE in large-scale underground excavations. International Journal of Rock Mechanics and Mining Sciences 44(4): 550–564.

Castanon

Blanco

Castro

-FD

, et al. (2018) Use of explicit FEM models for the structural and parametrical analysis of rockfall protection barriers. Engineering Structures 166: 212–226.

Chen

Zhao

Tomas

, et al. (2023) Retrieving the kinematic process of repeated-mining-induced landslides by fusing SAR/InSAR displacement, logistic model, and probability integral method. Remote Sensing 15(12): 3145.

Coulibaly

Chanut

Lambert

, et al. (2019) Toward a generic computational approach for flexible rockfall barrier modeling. Rock Mechanics and Rock Engineering 52(11): 4475–4496.

Cui

Liang

Pei

, et al. (2023) Structural characteristics of landslide failure boundaries using three-dimensional point clouds: a case study of the Zhaobiyan landslide, China. Bulletin of Engineering Geology and the Environment 82(4): 127.

Dhakal

Bhandary

Yatabe

, et al. (2011) Experimental, numerical and analytical modelling of a newly developed rockfall protective cable-net structure. Natural Hazards and Earth System Sciences 12(12): 3197–3212.

10.

di Prisco

Redaelli

Zerbi

(2023) Inclined block impacts on granular strata: coupled DEM-FDM numerical investigation and rheological modelling. Rock Mechanics and Rock Engineering 56: 6761–6778.

11.

Escallón

Boetticher

Wendeler

, et al. (2015) Mechanics of chain-link wire nets with loose connections. Engineering Structures 101: 68–87.

12.

Gentilini

Govoni

Miranda

, et al. (2012) Three-dimensional numerical modelling of falling rock protection barriers. Computers and Geotechnics 44: 58–72.

13.

Gottardi

Govoni

(2010) Full-scale modelling of falling rock protection barriers. Rock Mechanics and Rock Engineering 43(3): 261–274.

14.

Tang

, et al. (2019) Assessment of debris-flow potential dangers in the Jiuzhaigou Valley following the August 8, 2017, Jiuzhaigou earthquake, Western China. Engineering Geology 256: 57–66.

15.

Kwan

JSH

Sun

Lam

, et al. (2016) Recent advances in landslide risk management measures in Hong Kong. In: 12th International Symposium on Landslides, Napoli, ITALY, 1219–1227.

16.

Leonardi

Wittel

Mendoza

, et al. (2016) Particle-Fluid-Structure interaction for debris flow impact on flexible barriers. Computer-Aided Civil and Infrastructure Engineering 31(5): 323–333.

17.

Gong

Tang

, et al. (2021) Probabilistic analysis of a discrete element modelling of the runout behavior of the Jiweishan landslide. International Journal for Numerical and Analytical Methods in Geomechanics 45(8): 1120–1138.

18.

Sheng

(2023) Seismic performance and failure mechanism of interbedded slopes with steep rock layers. Engineering Geology 326: 107312.

19.

Ling

Sun

, et al. (2021) Characterizing the distribution pattern and geologic and geomorphic controls on earthquake-triggered landslide occurrence during the 2017 Ms7.0 Jiuzhaigou earthquake, Sichuan, China. Landslides 18(4): 1275–1291.

20.

Liu

(2022) Simulation of the earthquake-induced soil-rock mixed accumulation body sliding movement using discrete-continuous coupled approach. Natural Hazards 114(2): 2087–2108.

21.

Mentani

Govoni

Giacomini

, et al. (2018) An equivalent continuum approach to efficiently model the response of steel wire meshes to rockfall impacts. Rock Mechanics and Rock Engineering 51(9): 2825–2838.

22.

Nicot

Cambou

BMG

(2001) Design of rockfall restraining nets from a discrete element modelling. Rock Mechanics and Rock Engineering 34(2): 99–118.

23.

Zhao

, et al. (2018) A new numerical modeling approach for flexible rockfall protection barriers based on failure modes. Advanced Steel Construction 14(3): 479–495.

24.

Sautter

Hofmann

Wendeler

, et al. (2021) Advanced modeling and simulation of rockfall attenuator barriers via partitioned DEM-FEM coupling. Frontiers in Built Environment 7: 659382.

25.

Scuderi

Onyango

Nagle

(2023) A remote sensing and GIS analysis of Rockfall Distributions from the 5 July 2019 Ridgecrest (MW>7.1) and 24 June 2020 Owens Lake (MW>5.8) earthquakes. Remote Sensing 15(8): 1962.

26.

Shao

, et al. (2023) Insight into the characteristics and triggers of loess landslides during the 2013 heavy rainfall event in the tianshui area, China. Remote Sensing 15(17): 4304.

27.

Shi

Zhao

Zhang

, et al. (2021) Coupled discrete-continuum approach for railway ballast track and subgrade macro-meso analysis. International Journal of Pavement Engineering 22(13): 1744–1759.

28.

Tahmasbi

Giacomini

Wendeler

, et al. (2019) On the computational efficiency of the hybrid approach in numerical simulation of rockall flexible chain-link mesh. Rock Mechanics and Rock Engineering 52(10): 3849–3866.

29.

Tang

Zeng

(2024) Rockburst-specific flexible protection net (RFPN) for tunnels: creation, characteristics, and preliminary application. Rock Mechanics and Rock Engineering 57(7): 1–19.

30.

Tomás

Pinheiro

Pinto

, et al. (2023) Preliminary analysis of the mechanisms, characteristics, and causes of a recent catastrophic structurally controlled rock planar slide in Esposende (northern Portugal). Landslides 20(8): 1657–1665.

31.

Verrucci

Forte

, et al. (2023) Instantaneous limit equilibrium back analyses of major rockslides triggered during the 2016-2017 central Italy seismic sequence. Natural Hazard Earth System 23(3): 1177–1190.

32.

Wang

Xiao

Shi

, et al. (2022) Forensic analysis and numerical simulation of a catastrophic landslide of dissolved and fractured rock slope subject to underground mining. Landslides 19(5): 1045–1067.

33.

Gentilini

, et al. (2018) An energy allocation based design approach for flexible rockfall protection barriers. Engineering Structures 173: 831–852.

34.

Yang

Duan

, et al. (2019) A review of flexible protection in rockfall protection. Natural Hazards 99(1): 71–89.

35.

Yang

Tong

, et al. (2024) Deep learning-based automatic rockfall impact force reconstruction for flexible barrier systems in full-scale tests. Automation in Construction 165: 105510.

36.

Liu

Guo

, et al. (2019) Nonlinear numerical modeling of the wire-ring net for flexible barriers. Shock and Vibration 2019: 3040213.

37.

Zhao

(2018) Analysis of Deformation and Energy Dissipation Performance of the Passive Flexible Protection Struction. Doctor Degree Dissertation. South Jiaotong University. in Chinese.

38.

Zhao

, et al. (2020) Coupled numerical simulation of a flexible barrier impacted by debris flow with boulders in front. Landslides 17(12): 2723–2736.