Dynamic process and impact force evolution characteristics of debris flow entraining irregular boulders impacting check dams

In studies on the impact forces of debris flow boulders, most scholars simplify boulders as spheres, thereby neglecting the influence of their irregular shapes. Therefore, this study adopts a coupled SPH-FEM method to establish a numerical flume model for investigating the impact of irregular boulders in debris flows on check dams. The validity of this method is verified through field flume tests. Systematic simulations are conducted to analyze the interactions among the slurry, single boulders, multiple boulders, and check dams. The results indicate that with the increase in slope gradient and boulder mass, the boulder impact forces increase significantly, while the buffering effect of the slurry gradually weakens. During the movement of multiple boulders, collisions and energy exchange between boulders lead to a reduction in their impact forces compared to those of single boulders. Specifically, in the movement of single and multiple boulders, the peak impact force of Boulder 1 decreases by approximately 12.97%–27.50%, and that of Boulder 2 decreases by about 14.60%–30.87%. The research findings provide a theoretical basis for the dynamic response analysis of rigid check dams and the design of engineering protection under debris flow hazard conditions.

1 Introduction

Debris flow is a gravity-driven multiphase flow commonly occurring in mountain gullies, often triggered by earthquakes and rainfall. It is characterized by high velocity, intense mobility, strong destructive power, and unpredictability, capable of transporting boulders (Schöffl et al., 2023; Yongjie et al., 2024). The boulders entrained in debris flows can generate destructive impact forces on rigid check dams (Song et al., 2018; Shen et al., 2018; Zhao et al., 2020; Ng et al., 2021; Luo et al., 2022). Therefore, the impact of these boulders must be considered in the design of rigid check dams. To investigate the fluid-particle-structure interaction, current methodologies are generally categorized into three types: semi-empirical analytical methods, experimental methods, and numerical simulation methods (Liu et al., 2021). Research shows that semi-empirical analytical methods for calculating debris flow impact forces primarily fall into two categories. One category focuses on the fluid phase, with the hydrostatic model, hydrodynamic model, and mixed model being primarily adopted (Vagnon, 2020), where the empirical coefficients in these formulas vary depending on the factors considered. The other category addresses the impact forces exerted by boulders in debris flows, primarily adopting the impulse method, energy method, and Hertz contact theory, with simplified impact force formulas based on Hertz contact theory being predominant (He et al., 2016).

In recent years, numerous scholars have conducted extensive flume experiments to study the impact of large boulders in debris flows. Wang et al. (2019) performed large-scale laboratory tests, and the experimental results showed that large boulders entrained in debris flows often lead to peak impact pressures. Wendeler et al. (2019) combined laboratory tests with full-scale field test results to establish a debris flow load model for flexible barrier structures in mountain torrents, providing a basis for barrier design. Ng et al. (2021) compared the impact dynamics of debris flows with and without large boulders on rigid barriers through large-scale flume tests, proposing a new estimation equation for the load reduction coefficient Kc considering this buffering effect. Lam and Wong (2021) analyzed the impact mechanism and dynamic earth pressure coefficient of bouldery debris flows on reinforced concrete barriers using a 28-m-long flume. Kim and Yune (2025) investigated the influence of arrays of cylindrical baffles in front of a rigid barrier, with different transverse blockage ratios and numbers of rows, on the flow characteristics, dynamic impact force, and boulder travel distance of debris flows containing large boulders through small-scale flume experiments. Wang et al. (2025) explored the influence of boulder size and debris flow density on boulder movement through flume experiments.

With the rapid development of computer technology, numerical simulation has become a powerful analytical tool widely used in various fields. In the numerical simulation of debris flow fluids, an increasing number of researchers are using numerical methods to analyze the flow characteristics, impact forces, and environmental effects of debris flows. For the numerical simulation of debris flows, methods for the fluid domain mainly include mesh-based and particle-based approaches (Jia et al., 2023). Traditional mesh-based methods, such as the Finite Volume Method (FVM) (Mao et al., 2020), have been widely employed in debris flow simulations. For example, Li and Zhao (2018) proposed a debris flow simulation method combining Computational Fluid Dynamics (CFD) and the Discrete Element Method (DEM). However, under conditions of extreme deformation, these mesh-based methods used in such studies can encounter numerous numerical difficulties (such as severe mesh tangling, distortion, and degradation). Although accuracy and stability can be improved through remeshing (Khayyer et al., 2019), this also increases computational difficulty and reduces efficiency (Huaqing et al., 2022). To overcome this issue, Smoothed Particle Hydrodynamics (SPH) provides a more flexible solution within the Lagrangian framework, successfully avoiding mesh dependency issues. SPH is a meshless Lagrangian particle method that discretizes the fluid using particles and approximates the governing equations using kernel functions, making it highly suitable for simulating free-surface flows, large deformations, and complex boundary problems. Due to its ability to approximate mathematical equations through particle interactions, SPH has become an effective analytical tool for solving complex fluid mechanics problems and is widely applied across various engineering fields, demonstrating significant advantages, particularly in simulating complex fluid dynamic behaviors. Numerous researchers (Zhang et al. (2011); Wang (2012)), Hu et al (2015) have applied the coupled SPH-FEM method in the field of impact dynamics. Law et al. (2018) and Liu et al. (2019) applied the SPH method to simulate the fluid flow behavior of debris flows and achieved promising results. Xiong et al. (2023) and Qiao et al. (2023) successfully used the SPH method to simulate the fluid motion characteristics of debris flows. Wei et al. (2016), Dai et al. (2017), and Shi et al. (2022) used the SPH method to simulate debris flows and successfully predicted their propagation processes. Feng et al. (2019), Chen et al. (2019), and Luo et al. (2019) verified the effectiveness of the Finite Element Method (FEM) in simulating the progressive failure of structures induced by debris flow impacts. Li et al. (2025) revealed the buffering effect of debris flow slurry on boulders, identified boulder impact as the primary cause of local damage to frame structures, and clarified the significant influence of initial debris flow velocity on structural damage.

Both in flume model tests and numerical simulations, scholars such as Yongjie et al. (2024), Ng et al. (2021), Luo et al. (2022), Liu et al. (2019), and Lu et al. (2025) have uniformly adopted spherical blocks. However, boulders under natural conditions are usually irregular in shape, thus the influence of their shape on the impact force is neglected.

In this study, the geometries of irregular boulders were acquired through 3D scanning, and the SPH-FEM coupled analysis method was adopted. The numerical model was validated through field physical tests to ensure its reliability. The focus of this study is to investigate the impact forces exerted by irregular boulders in debris flows on rigid check dams, thereby providing a more comprehensive and scientific reference for the design and risk assessment of rigid check dams against debris flow hazards.

2 Numerical methodology

2.1 SPH method

SPH is a meshless particle method based on the Lagrangian description. It effectively addresses the problem of mesh distortion inherent in traditional methods and offers notable advantages in handling problems such as free-surface flows and deformable boundaries. The core computational principle of SPH is based on kernel approximation and particle approximation theories. Its essence lies in solving continuum mechanics equations in a meshless manner using discrete particle groups.

In the Smoothed Particle Hydrodynamics (SPH) method, the first step is to employ a meshless kernel function (Yao et al., 2024), which is defined by the following equations:

$〈f\left(x\right)〉=\int f\left(x^{\prime }\right)W\left(x-x^{\prime },h\right)\mathrm{d}{x}^{\prime }$

$〈\nabla f\left(x\right)〉=\int \nabla f\left(x^{\prime }\right)W\left(x-x^{\prime },h\right)\mathrm{d}{x}^{\prime }$

where $W\left(\mathrm{x}-{\mathrm{x}}^{\prime },h\right)$ is the kernel function, and $h$ is the smoothing length. $〈f\left(\mathrm{x}\right)〉$ and $〈\nabla f\left(\mathrm{x}\right)〉$ are the approximations of the function $f\left(\mathrm{x}\right)$ and its derivative, respectively.

This study employs the cubic spline kernel function (Yongjie et al., 2024) (Liu et al., 2021) to balance computational accuracy and efficiency.

$W\left(r,h\right)=C\left\{\begin{array}{cc}1-\frac{3}{2}{R}^2+\frac{3}{4}{R}^3& 0\le R<1\\ \frac{1}{4}{\left(2-R\right)}^3& 1\le R<2\\ 0& R\ge 2\end{array}\right.$

where $W\left(r,h\right)$ is the centrally-peaked function ($r=\mid \mathbf{x}-{\mathbf{x}}^{\prime }\mid$), $R=r_{ij}/h$, $C$ is a normalization constant whose value depends on the spatial dimensions (in three-dimensional problems, $C=1/\left(\pi h^3\right)$).

The second step in the SPH method is particle approximation. This step discretizes the computational domain into particles that carry physical properties such as mass, density, velocity, and stress. By superimposing the contributions of these particles via the kernel function, the approximate expressions for the function $f\left({\mathbf{x}}_i\right)$ and its derivative $\nabla f\left({\mathbf{x}}_i\right)$ at particle $i$ are obtained as follows:

$〈f\left(x_i\right)〉={\displaystyle \sum _{j=1}^N}\frac{m_j}{{\rho }_j}f\left(x_j\right)·W_{ij},$

$〈\nabla ·f\left(x_i\right)〉=-{\displaystyle \sum _{j=1}^N}\frac{m_j}{{\rho }_j}f\left(x_j\right)·{\nabla }_j{W}_{ij},$

where $i,j$ are particle indices, $N$ is the total number of particles within the support domain, ${\rho }_j$ is the density of particle $j$, $m_j$ is the mass of particle $j$, and ${\nabla }_j{W}_{ij}$ is the derivative of the kernel function $W_{ij}$ at particle $j$ (${\nabla }_j{W}_{ij}=-{\nabla }_i{W}_{ij}$).

In the SPH simulation of debris flow, treated as a weakly compressible fluid. The pressure term is explicitly updated through an equation of state (EOS) to avoid the computational cost of implicit solutions. The EOS-Gruneisen adopted in this study is defined as follows (Yongjie et al., 2024; Chen et al., 2019; Yao et al., 2024):

$p=\frac{{\rho }_0{C}^2\mu \left[1+\left(1\frac{{\gamma }_0}{2}\right)\mu -\frac{a}{2}{\mu }^2\right]}{{\left[1-\left(S_1-1\right)\mu -S_2\frac{{\mu }^2}{\mu +1}-S_3\frac{{\mu }^3}{{\left(\mu +1\right)}^2}\right]}^2}+\left({\gamma }_0+a\mu \right)E_0$

where $C$ is the intercept of the adiabatic speed of sound curve, $S_1$, $S_2$, and $S_3$ are fitting coefficients, ${\gamma }_0$ is the Gruneisen coefficient, $a$ is the volume correction coefficient for ${\gamma }_0$, $E_0$ is the internal energy per initial volume, $\mu =\rho /{\rho }_0-1$, and ${\rho }_0$ and $\rho$ are the initial density and the changed density, respectively.

2.2 SPH-FEM coupling mechanism

To fully leverage the advantages of both SPH and FEM, coupling of the two methods is necessary (Figure 1). The SPH-FEM coupling algorithm follows a primary-secondary approach, where SPH particles are designated as secondary nodes and FEM elements are designated as primary nodes. This contact algorithm continuously monitors potential penetration of SPH particles into the segmented surfaces of FEM elements. The “contact thickness” parameter defines the outward extension from the contact surface where physical contact is evaluated. When any SPH particle enters this detection zone, it experiences a repulsive force that scales proportionally with either the predicted penetration depth or the contact spring stiffness.

Figure 1

Figure 1. Schematic diagram of SPH-FEM coupling.

The formula for the normal contact force ${\mathbf{f}}_{sf,n}$ is as follows:

$f_{sf,n}=\left(k_{sf,n}{\delta }_{sf}+c_{sf,n}{\dot{\delta }}_{sf}\right)n$

where $k_{sf,n}$ is the normal spring stiffness, ${\dot{\delta }}_{sf}$ is the relative normal velocity, $c_{sf,n}$ is the normal damping coefficient, and $\mathbf{n}$ is the unit normal displacement vector.

The formula for the normal spring stiffness $k_{sf,n}$ is (Jiang et al., 2025):

$k_{\text{sf},\mathrm{n}}=\mathit{max}\left(SOFSCL\frac{m}{2\Delta t^2},k_1\frac{K{s}^2}{V}\right)$

where $SOFSCL=0.1$ (soft constraint penalty factor, default value), $k_1=0.1$ (penalty scale factor, default value), $K$ is the bulk modulus of the material, $s$ is the segment area, and $V$ is the element volume.

The formula for the SPH-FEM overlap ${\delta }_{sf}$ is:

${\delta }_{\text{sf}}=h_i-d_{\mathit{min}}$

where $h_i$ is the smoothing length of the SPH particle, and $d_{\mathit{min}}$ is the distance from the center of the SPH particle to the contact surface.

The formula for the tangential contact force ${\mathbf{f}}_{sf,t}$ is (Lu et al., 2025):

$f_{sf,t}=\left\{\begin{array}{l}\left(k_{sf,t}{\delta }_{sf,t}+c_{sf,t}{\dot{\delta }}_{sf,t}\right)\left|f_{sf,n}\right|\mu >\left|k_{sf,t}{\delta }_{sf,t}+c_{sf,t}{\dot{\delta }}_{sf,t}\right|\\ \frac{\left(k_{sf,t}{\delta }_{sf,t}+c_{sf,t}{\dot{\delta }}_{sf,t}\right)}{\left|k_{sf,t}{\delta }_{sf,t}+c_{sf,t}{\dot{\delta }}_{sf,t}\right|}\left|f_{sf,n}\right|\mu \left|f_{sf,n}\right|\mu \le \left|k_{sf,t}{\delta }_{sf,t}+c_{sf,t}{\dot{\delta }}_{sf,t}\right|\end{array}\right.$

where $k_{sf,t}$ is the tangential spring stiffness, ${\delta }_{sf,t}$ is the incremental tangential displacement, ${\dot{\delta }}_{sf,t}$ is the rate of change of the tangential displacement increment (relative tangential velocity), $\mu$ is the friction coefficient, and $c_{sf,t}$ is the tangential damping coefficient.

3 Flume model construction

3.1 Field model construction

To investigate the impact dynamics of viscous debris flow on rigid check dams, this study constructed a large-scale on-site debris flow flume test facility (Figure 2). The experimental apparatus comprises three core components: a slurry preparation pond, a debris flow channel, and a signal acquisition zone. The key advantage of this on-site facility lies in the use of locally sourced soil and boulders, which maximizes the authenticity of their physical and mechanical properties, thereby making the generated debris flow phenomena more representative of actual conditions. The particle size distribution curve of the debris flow slurry is presented in Figure 3. While previous researchers including Ng et al. and Wang et al. conducted physical debris flow experiments using 6 m (Ng et al., 2021) and 15 m (Wang et al., 2025) long flumes respectively, this study constructed an 11 m long debris flow channel with a cross-section of 1 m × 1 m and a slope of 16° based on specific site conditions. The channel sides were constructed using geotextile soil bags arranged in a staggered pattern, providing not only excellent structural stability but also effective simulation of the rough interface and permeability characteristics of natural soil slopes. The channel slope of 16° was determined with reference to the gradient of the initiation zone in the Liandeng Gully debris flow, thereby directly replicating the steep terrain conditions found in nature. The signal acquisition zone was equipped with high-precision impact force sensors (Figure 4), consisting of two main parts: smaller sensors in the lower section for measuring the impact forces of debris flow slurry, and larger sensors in the upper section specifically designed to capture the impact forces of boulders. This configuration enabled the simultaneous separation and real-time acquisition of impact forces from both the slurry and boulders in the debris flow. All objects in the experiment started with zero initial velocity and were set in motion solely by gravitational force. To comprehensively observe the dynamic movement process of the debris flow, multiple cameras were strategically positioned at key locations along the large-scale flume. These cameras documented the entire process, recording the initiation, transportation, and deposition characteristics of the debris flow.

Figure 2

Figure 2. Field experimental setup.

Figure 3

Figure 3. Particle size distribution curve of the slurry.

Figure 4

Figure 4. Impact force sensor installation.

3.2 Numerical model construction

This study developed three numerical models of flume tests (Figure 5) using the LS-DYNA platform, representing three distinct working conditions. Figures 5a,b show Boulder1 and Boulder 2 placed at different locations with different masses, while Figure 5c depicts both boulders simultaneously positioned in the debris flow channel. The debris flow slurry was simulated using the Smoothed Particle Hydrodynamics (SPH) method with the MAT_NULL material model (Chen et al., 2019) (Lu et al., 2025). A total of 23,125 SPH particles were generated with a particle spacing of 0.04 m. The EOS_GRUNEISEN equation of state was employed for fluid pressure calculation (Yongjie et al., 2024; Chen et al., 2019). The interaction between the debris flow slurry and both the channel and rigid check dam was defined using the CONTACT_AUTOMATIC_NODES_TO_SURFACE algorithm, while the CONTACT_AUTOMATIC_SURFACE_TO_SURFACE algorithm was applied to characterize interactions between boulders and the channel and check dam. Furthermore, Form16 was selected as the control equation in CONTROL_SPH, and a density reinitialization scheme was implemented to mitigate pressure oscillations caused by boundary effects in conventional SPH methods. The check dam, irregular debris flow boulders, and channel were all simulated using the MAT_RIGID material model (Shen et al., 2018; Li et al., 2020), with the check dam and channel subjected to fixed constraints. The numerical model of the debris flow channel was built at a 1:1 scale relative to the physical experiment, featuring identical dimensions of 11 m length, 1 m width, and 1 m height, with a 16° slope. The irregular boulders were geometrically reconstructed based on 3D scanning data (Figure 6). Boulder 1 and Boulder 2 were selected as representatives of typical debris flow boulders: Boulder 1 exhibits significant elongation along its long axis, while Boulder 2 has approximately equal dimensions along all three axes. The input parameters for the numerical model are summarized in Table 1. All motions in the numerical model started with zero initial velocity and were driven solely by gravitational force.

Figure 5

Figure 5. Numerical models: (a) Single boulder 1, (b) Single boulder 2, (c) Multiple boulders.

Figure 6

Figure 6. 3D scans of irregular boulders: (a) Scan of boulder 1, (b) Scan of boulder 2.

Table 1

Table 1. Material parameters used in the model.

4 Results and analysis

4.1 Model validation

To validate the accuracy and reliability of the established numerical model, this section systematically presents a comparative analysis between numerical simulations and field physical experiments. During the numerical model construction, the debris flow flume slope was set at 16°, consistent with the field test conditions. Figure 7 shows the dynamic process captured during the field physical experiment. After t = 1.70 s, the slurry and boulders became fully mixed, making it difficult to observe the detailed movement of individual boulders. However, footage from alternative camera angles confirmed boulder impacts on the sensors, demonstrating that the impact forces recorded by these sensors were indeed from the boulders. Figure 8 presents the time-history curves of impact forces for both debris flow slurry and boulders, comparing numerical simulations with field physical experiments. The field tests recorded a peak impact force of 2614.3 N for the slurry and 7.50 kN for the boulder. The numerical simulations yielded corresponding values of 2764.63 N for the slurry and 7.41 kN for the boulder, representing relative errors of approximately 5.59% for the slurry and 1.2% for the boulder——both within acceptable margins. The discrepancies are mainly attributed to the inherent non-uniformity of the materials, the simplifications inherent in numerical simulations, and the inevitable differences between the idealized simulations and the actual on-site conditions. Beyond the close agreement in terms of peak impact forces between numerical simulations and field experiments, the simulated time-history curves of impact forces for both the slurry and the boulders also exhibit consistent trends with field observations. These results demonstrate that the numerical model developed in this study can effectively simulate the interaction between boulder-laden debris flow and rigid check dam.

Figure 7

Figure 7. Dynamic process of the physical experiment.

Figure 8

Figure 8. Impact Force Evolution of Debris Flow Slurry and Boulder: (a) Slurry, (b) Boulder (The flume slope is 16°).

The Froude number (Fr), serving as a key dimensionless parameter representing the relative magnitude of inertial force to gravitational force in a fluid, directly governs the dynamic characteristics of the flow system. The corresponding Froude number can be calculated using the following modified formula (Jia et al., 2023):

$Fr=\frac{v}{\sqrt{gh\cos \theta }}$

where:$v$ is the flow velocity in front of the barrier, $g$ is the gravitational acceleration, $h$ is the flow depth, $\theta$ is the channel slope.

In this study, the Froude number measured in the field physical experiment was 4.19, while the value obtained through numerical simulation was 4.13 (Table 2). The relative deviation between the two is only about 1.43%, indicating good numerical agreement. This demonstrates that both the experimental and numerical methods exhibit high reliability in obtaining this key dynamic parameter—the Froude number. Both values fall within the natural Froude number range of 1–10 (Lam and Wong, 2021; Wang et al., 2025; Zhou et al., 2018; Choi et al., 2020). The flow processes reproduced by both the field physical experiment and the numerical simulation are essentially consistent with those of natural debris flows in terms of dynamic characteristics, thus effectively reflecting the motion dominated by the combined effects of inertial and gravitational forces in natural debris flows. This provides a solid similarity basis for subsequent mechanistic analysis and derivation of behavioral mechanisms.

Table 2

Table 2. Calculation of froude number.

4.2 Coupling analysis of a single boulder

This section elucidates the dynamic characteristics of a single boulder within a debris flow impacting a rigid check dam. Taking the case of Boulder 1 on a flume slope of 25° as an example, Figure 9 illustrates the evolution of slurry and boulder movement during the impact process from t = 0 s to t = 2.52 s. Motion initiates at t = 0 s. At t = 0.22 s, the slurry begins to interact with the boulder. From t = 0.22 s to t = 2.52 s, the boulder and slurry move together until the boulder collides with the rigid check dam. At t = 1.79 s, the leading edge of the slurry impacts the dam, and at t = 2.52 s, the boulder collides with the dam. It can be observed that the slurry consistently remains at the leading edge of the debris flow movement, indicating continuous interaction between the slurry and the boulder after t = 0.22 s. When the slurry front reaches the dam, it collides with the structure, leading to a decrease in the velocity of the debris flow slurry. A portion of the slurry accumulates in front of the dam and interacts with the still-descending boulder, further reducing the boulder’s velocity. Consequently, the boulder does not impact the dam at its maximum velocity, demonstrating that the accumulated slurry exerts a significant cushioning effect on the impact force of the boulder.

Figure 9

Figure 9. Evolution of Boulder 1 Impacting the Rigid Check Dam (T = 0 s–2.52 s).

Figure 10 displays the time-history curves of impact forces exerted by a single boulder on the rigid check dam under gradually increasing flume slope gradients. Analysis reveals that the variation trends of impact forces for both Boulder 1 and Boulder 2 are generally consistent: the forces increase rapidly, reach a peak impact force, then decline gradually before eventually stabilizing. When the slope exceeds 20°, two distinct peaks are clearly observed in the impact force time-history curves, with the amplitude of the second peak consistently being lower than that of the first. This phenomenon is primarily attributed to the rebound of the boulder after its initial impact with the dam, leading to a secondary impact. Furthermore, as the slope continues to increase, the timing of the maximum impact force advances progressively, reflecting that steeper slopes enhance both the intensity of boulder motion and the abruptness of the impact.

Figure 10

Figure 10. Process of single boulder impacting the rigid check dam: (a) Boulder 1, (b) Boulder 2.

Figure 11 illustrates the relationship between flume slope and boulder motion velocities (impact velocity and maximum velocity). As shown in Figure 11a for Boulder 1, as the slope increases from 20° to 40°, both the maximum velocity and the impact velocity generally show an increasing trend, with the maximum velocity consistently exceeding the impact velocity. However, the difference between the two velocities gradually narrows as the slope gradient increases. From Figure 11b for Boulder 2, it is observed that when the slope reaches approximately 35°, the impact velocity and maximum velocity coincide. The discrepancy in their velocity variation patterns between Boulder 1 and Boulder 2 primarily stems from the coupled effects of boulder mass and slope gradient, where Boulder 2 has a mass of 229.39 kg, significantly greater than Boulder 1’s mass of 138.56 kg. Analysis indicates that during debris flow movement, as the slope gradient increases and the boulder mass increases, the buffering effect of the slurry on the boulders significantly weakens.

Figure 11

Figure 11. Relationship between flume slope and velocity: (a) Boulder 1, (b) Boulder 2.

4.3 Coupling analysis of multiple boulders

Figure 12 illustrates the evolution characteristics of the joint movement of the slurry and boulders during the impact of multiple boulders on a rigid check dam at a flume slope of 20° with both Boulder 1 and Boulder 2, from t = 0.23 s to t = 3.29 s. The specific process is as follows: at t = 0.23 s, the slurry collides with Boulder 1; subsequently, at t = 0.65 s, the slurry collides with Boulder 2; at t = 0.91 s, Boulder 1 and Boulder 2 undergo mutual collision. Thereafter, the slurry moves together with both Boulder 1 and Boulder 2 and they sequentially interact with the rigid check dam: at t = 2.00 s, the leading edge of the slurry impacts the dam first; at t = 3.16 s, Boulder 2 collides with the dam; finally, at t = 3.29 s, Boulder 1 impacts the dam. It is noteworthy that Boulder 2 collides with the dam approximately 0.13 s earlier than Boulder 1. The collision with Boulder 2 exerted a retarding effect on the movement of Boulder 1.

Figure 12

Figure 12. Evolution process of multiple boulders impacting the rigid check dam (T = 0.23s–3.29 s).

Figure 13 presents the time-history curves of impact forces exerted by Boulder 1 and Boulder 2 on the rigid check dam during the movement of multiple boulders. Overall, the variation pattern of the impact forces is consistent with that observed in single-boulder impacts: upon contact with the dam, the impact force rapidly increases to a peak value, then gradually attenuates and stabilizes. However, compared to single-boulder impacts, the impact process under multiple boulders conditions exhibits certain differences. In single-boulder impacts, multiple peaks may appear in the impact force time-history curves due to strong initial collisions with the dam, often followed by rebound or secondary contact. In contrast, under multiple boulders conditions, in addition to direct interactions between the boulders and the dam, mutual collisions and interference among the boulders also occur. These interactions dissipate part of the impact energy, thereby significantly reducing the occurrence of distinct peaks in the impact force curves. Therefore, the multi-boulder impact process not only reflects the direct effects on the dam but also highlights the influence of dynamic coupling effects between the boulders on the characteristics of the impact force process.

Figure 13

Figure 13. Process of multiple boulders impacting the rigid check dam: (a) Boulder 1, (b) Boulder 2.

Figure 14 presents a comparative analysis of the peak impact forces under single-boulder versus multiple-boulder conditions. Overall, the peak impact forces under multiple-boulder conditions are reduced compared to those under single-boulder conditions. Specifically, for Boulder 1, the peak impact force decreases by approximately 13%–27% in multiple-boulder scenarios. For Boulder 2, the most significant reduction in peak impact forces occurs at a slope of 20°, with a decrease of about 31%, while reductions under other working conditions remain within the range of 15%–26%. These results indicate that interactions between the boulders significantly influence the impact effects. This suggests that during the process of mutual interference between the multiple boulders, part of the impact energy is dissipated, leading to a reduction in the impact intensity of the boulders on the rigid check dam.

Figure 14

Figure 14. Comparison of peak impact forces: (a) Boulder 1, (b) Boulder 2.

Taking the 20° slope case as an example, Figures 15, 16 comparatively present the normal velocity and resultant velocity curves during the movement of the single and multiple boulders, respectively. The general consistency observed between the normal and resultant velocity trends indicates that the boulder impact force is predominantly governed by the normal velocity component. During the initial motion stage, both scenarios exhibit similar behavior characterized by gradually increasing velocities over time, demonstrating that the boulders’ movement is primarily controlled by gravitational forces and slurry interactions. However, as the motion progresses, inter-boulder collisions in the multiple boulders scenario lead to distinctive velocity characteristics compared to the single boulder case: during collision instances, Boulder 1 experiences velocity reduction while Boulder 2 shows corresponding increase, clearly reflecting mutual energy exchange between the two boulders. Following collisions, the overall velocities in multiple boulders cases remain lower than those in single boulder cases, accompanied by significant fluctuations primarily due to repeated inter-boulder collisions. After collision with the check dam, all boulders exhibit rapid velocity decay approaching zero, regardless of their being single or multiple. These findings demonstrate that compared to single boulder movement, multiple boulder motion involves greater complexity where velocity evolution not only reflects external forces but more importantly reveals the crucial role of inter-boulder interactions in the dynamics of multi-boulder systems.

Figure 15

Figure 15. Resultant velocity comparison at 20° slope: (a) Boulder 1, (b) Boulder 2.

Figure 16

Figure 16. Normal velocity comparison at 20° slope: (a) Boulder 1, (b) Boulder 2.

5 Discussion

Three-dimensional scanning technology was employed to reconstruct the real geometric morphologies of irregular boulders, overcoming the limitations of simplifying boulders as spherical particles in conventional studies. Validation through field experiments shows relative errors of only 5.59% for the slurry impact force and 1.2% for the boulder impact force, with a Froude number deviation of 1.43%, demonstrating the model’s reliability for simulating complex debris flow-check dam interactions.

From a dynamic mechanism perspective, the slurry exerts a buffering effect on the impact forces of boulders. This phenomenon is consistent with the findings reported by Yongjie et al. (2024), Ng et al. (2021) and Luo et al. (2022), but the present study further quantifies the influences of slope gradient and boulder mass on this buffering effect. Specifically, as the slope gradient increases, the slurry’s buffering effect on boulder velocity decreases (Figure 17); meanwhile, a larger boulder mass results in a weaker buffering effect of the slurry (Figure 18). All boulders in the figures were modeled with the geometric shape of Boulder 2 to minimize the influence of other variables on the simulation results. These results provide crucial theoretical support for the design of check dams in high-slope debris flow watersheds.

Figure 17

Figure 17. Comparison Diagram of Impact Force Processes and Speeds for Different Slopes (a) Impact Force Processes (b) Speed Comparison.

Figure 18

Figure 18. Comparison Diagram of Impact Force Processes and Speeds for Different Masses (a) Impact Force Processes (b) Speed Comparison.

Comparative analysis of single and multiple boulder impacts reveals that the total impact force in multi-boulder scenarios is not a simple superposition of individual boulder forces. Kinetic energy is partially dissipated through inter-boulder collisions and friction, reducing the overall impact. This finding suggests that considering only single boulder impacts may overestimate the destructive potential, while multi-boulder conditions better represent natural scenarios.

This study has several limitations. First, both the boulders and the check dam were modeled as rigid bodies in the simulation, without accounting for the effects of boulder breakage or local deformation of the dam under actual conditions, which may lead to an overestimation of the peak impact force. Second, only two typical shapes of irregular boulders were selected, which does not cover the full diversity of boulder shapes encountered in natural debris flow environments. As a result, the influence of irregular boulder geometries on impact forces has not been fully and accurately quantified.

6 Conclusion

Based on the coupled SPH-FEM method, this study establishes a three-dimensional numerical model to simulate the impact of irregular boulder-laden debris flows on rigid check dams, and the model is validated via large-scale on-site flume tests. Focusing on the analysis of impact characteristics and dynamic mechanisms under single-boulder and multiple-boulder impact scenarios, the main conclusions are as follows:

1. During the movement of irregular boulder-laden debris flows, the debris flow slurry can exert a buffering effect on the boulders’ impact force by dissipating the boulders’ energy. Since the slurry contacts the check dam first at the front of the flow, a portion of the slurry accumulates in front of the check dam, reducing the boulders’ movement velocity. Consequently, the boulders do not reach their maximum kinetic energy when impacting the dam, thereby lowering the peak impact force.

2. As the flume slope gradient increases, the impact force of the boulders also increases, and the timing of the peak impact force advances. This indicates that the boulder impact is more intense and abrupt under steeper slope conditions. Increasing the flume slope gradient accelerates the overall movement velocity of the debris flow, making the slurry’s buffering effect on boulder velocity increasingly weak. Meanwhile, Increasing the boulder mass enhances the boulders’ own inertial force, further reducing the buffering capacity of the debris flow slurry. Therefore, under the coupled effects of slope gradient and boulder mass, the impact intensity and abruptness of the boulders are significantly enhanced.

3. During the movement of multiple boulders, mutual collisions between boulders lead to the boulders’ movement velocity remaining consistently in a fluctuating state. Both the overall velocity and peak impact force of the boulders are attenuated compared to the single-boulder scenario, reflecting the energy dissipation effect during mutual boulder collisions. Specifically, compared with the single-boulder scenario, the peak impact force of Boulder 1 decreases by approximately 12.97%–27.50%, and that of Boulder 2 decreases by approximately 14.60%–30.87%. This indicates that the energy dissipation caused by mutual collisions between multiple boulders weakens the impact force of individual boulders on the rigid check dam.

Data availability statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author contributions

YW: Data curation, Software, Validation, Writing – original draft. HG: Writing – review and editing. GC: Writing – original draft, Writing – review and editing. SS: Project administration, Supervision, Writing – review and editing. LM: Writing – review and editing.

Funding

The author(s) declared that financial support was not received for this work and/or its publication.

Conflict of interest

Author SS was employed by Kunming Prospecting Design Institute of China Nonferrous Metals Industry Co., Ltd.

The remaining author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Chen, H., Li, J., Feng, S., Gao, H. Y., and Zhang, D. M. (2019). Simulation of interactions between debris flow and check dams on three-dimensional terrain. Eng. Geol. 251, 48–62. doi:10.1016/j.enggeo.2019.02.001

CrossRef Full Text | Google Scholar

Choi, C. E., Ng, C. W. W., Liu, H., and Wang, Y. (2020). Interaction between dry granular flow and rigid barrier with basal clearance: analytical and physical modelling. Can. Geotech. J. 57 (2), 236–245. doi:10.1139/cgj-2018-0622

CrossRef Full Text | Google Scholar

Dai, Z., Huang, Y., Cheng, H., and Xu, Q. (2017). SPH model for fluid–structure interaction and its application to debris flow impact estimation. Landslides 14 (3), 917–928. doi:10.1007/s10346-016-0777-4

CrossRef Full Text | Google Scholar

Feng, S., Gao, H., Gao, L., Zhang, L. M., and Chen, H. X. (2019). Numerical modeling of interactions between a flow slide and buildings considering the destruction process. Landslides 16 (10), 1903–1919. doi:10.1007/s10346-019-01220-9

CrossRef Full Text | Google Scholar

He, S., Liu, W., and Li, X. (2016). Prediction of impact force of debris flows based on distribution and size of particles. Environ. Earth Sci. 75 (4), 298. doi:10.1007/s12665-015-5180-2

CrossRef Full Text | Google Scholar

Hu, Y. G., Lu, W. B., and Chen, M. (2015). Implementation and verification of SPH-FEM coupled analysis method for blast-induced damage. Chin. J. Rock Mech. Eng. 34 (S1), 2740–2748. doi:10.13722/j.cnki.jrme.2014.0104

CrossRef Full Text | Google Scholar

Huaqing, M., Lianyong, Z., Zihan, L., Chen, M., Xia, X., and Zhao, Y. (2022). A review of recent development for the CFD-DEM investigations of non-spherical particles. Powder Technol. 412, 117972. doi:10.1016/j.powtec.2022.117972

CrossRef Full Text | Google Scholar

Jia, Z. S., Li, Z., Nicolin, G., Martínez-Estévez, I., Crespo, A. J., Sun, Z., et al. (2023). A resolved SPH-DEM coupling method for analysing the interaction of polyhedral granular materials with fluid. Ocean. Eng. 287 (P2), 115938. doi:10.1016/j.oceaneng.2023.115938

CrossRef Full Text | Google Scholar

Jiang, L., Xiao, W., Lai, Z., and Mou, B. (2025). Dynamic characteristics and impact load properties of high-speed railway piers under debris flow impact. Transp. Geotech. 52, 101562. doi:10.1016/j.trgeo.2025.101562

CrossRef Full Text | Google Scholar

Khayyer, A., Tsuruta, N., Shimizu, Y., and Gotoh, H. (2019). Multi-resolution MPS for incompressible fluid-elastic structure interactions in ocean engineering. Appl. Ocean Res. 82, 397–414. doi:10.1016/j.apor.2018.10.020

CrossRef Full Text | Google Scholar

Kim, J. B., and Yune, Y. C. (2025). Flume investigation of debris flow entrained boulders with cylindrical baffles and a rigid barrier. Eng. Geol. 344, 344 107836–107836. doi:10.1016/j.enggeo.2024.107836

CrossRef Full Text | Google Scholar

Lam, H. W. K., and Wong, A. L. (2021). Experimental and numerical study of dynamic soil debris impact load on reinforced concrete debris-resisting barriers. Landslides 18 (3), 955–966. doi:10.1007/s10346-020-01529-w

CrossRef Full Text | Google Scholar

Law, R. P. H., Kwan, J. S. H., and Ko, F. W. Y. (2018). “Three dimensional debris mobility modelling coupling smoothed particles hydrodynamics and ArcGIS,” in 19th International Conference on Soil Mechanics and Geotechnical Engineering, 3501–3504.

Google Scholar

Li, X., and Zhao, J. (2018). A unified CFD-DEM approach for modeling of debris flow impacts on flexible barriers. Int. J. Numer. Anal. Methods Geomechanics 42 (14), 1643–1670. doi:10.1002/nag.2806

CrossRef Full Text | Google Scholar

Li, B., Wang, C., and Li, Y. (2020). Dynamic response analysis of retaining dam under the impact of solid-liquid two-phase debris flow based on the coupled SPH-DEM-FEM method. Geofluids, 2020. doi:10.1155/2020/6635378

CrossRef Full Text | Google Scholar

Li, X., Li, B., Qie, L., Li, Y., and Wang, H. (2025). Dynamic response of a frame structure impacted by debris flow containing large boulders. Structures 74, 108543. doi:10.1016/j.istruc.2025.108543

CrossRef Full Text | Google Scholar

Liu, C., Yu, Z. X., and Luo, L. R. (2019). Study on dynamic behavior of concrete retaining dam under impact of debris flow containing large boulders. J. Vib. Shock 38 (14), 161–168+238. doi:10.13465/j.cnki.jvs.2019.14.023

CrossRef Full Text | Google Scholar

Liu, C., Yu, Z., and Zhao, S. (2021). A coupled SPH-DEM-FEM model for fluid-particle-structure interaction and a case study of wenjia gully debris flow impact estimation. Landslides 18 (7), 1–23. doi:10.1007/s10346-021-01640-6

CrossRef Full Text | Google Scholar

Lu, Y. C., Hao, M. J., Zhao, D. B., Liu, B., Han, Y., and He, L. (2025). Dynamic analysis of debris flow impact on double-column bridge piers based on SPH-FEM coupling. J. Shenzhen Univ. Sci. Eng. 42 (04), 419–427. doi:10.3724/sp.j.1249.2025.04419

CrossRef Full Text | Google Scholar

Luo, Y. H., Shen, P., and Zhang, M. L. (2019). How does a cluster of buildings affect landslide mobility: a case study of the shenzhen landslide. Landslides J. Int. Consortium Landslides 16 (3), 2421–2431. doi:10.1007/s10346-019-01239-y

CrossRef Full Text | Google Scholar

Luo, G., Zhao, Y., Shen, W., and Wu, M. (2022). An analytical method for the impact force of a cubic rock boulder colliding onto a rigid barrier. Nat. Hazards 112 (1), 1–16. doi:10.1007/s11069-021-05196-5

CrossRef Full Text | Google Scholar

Mao, J., Zhao, L., Di, Y., Liu, X., and Xu, W. (2020). A resolved CFD–DEM approach for the simulation of landslides and impulse waves. Comput. Methods Appl. Mech. Eng. 359, 112750. doi:10.1016/j.cma.2019.112750

CrossRef Full Text | Google Scholar

Ng, C. W. W., Liu, H., Choi, C. E., Kwan, J. S. H., and Pun, W. K. (2021). Impact dynamics of boulder-enriched debris flow on a rigid barrier. J. Geotechnical Geoenvironmental Eng. 147 (3), 04021004. doi:10.1061/(asce)gt.1943-5606.0002485

CrossRef Full Text | Google Scholar

Qiao, Z., Shen, W., Berti, M., and Li, T. (2023). An advanced SPH model for protective constructions of debris flows adopting the modified HBP constitutive law. Landslides 20 (11), 2437–2453. doi:10.1007/s10346-023-02123-6

CrossRef Full Text | Google Scholar

Schöffl, T., Nagl, G., Koschuch, R., Schreiber, H., Hübl, J., and Kaitna, R. (2023). A perspective of surge dynamics in natural debris flows through pulse-doppler radar observations. J. Geophys. Res. Earth Surf. 128 (9), e2023JF007171. doi:10.1029/2023jf007171

CrossRef Full Text | Google Scholar

Shen, W., Zhao, T., Zhao, J., Dai, F., and Zhou, G. G. (2018). Quantifying the impact of dry debris flow against a rigid barrier by DEM analyses. Eng. Geol. 241, 86–96. doi:10.1016/j.enggeo.2018.05.011

CrossRef Full Text | Google Scholar

Shi, H., Huang, Y., and Feng, D. (2022). Numerical investigation on the role of check dams with bottom outlets in debris flow mobility by 2d SPH. Sci. Rep. 12, 20456. doi:10.1038/s41598-022-24962-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Song, D., Choi, E. C., Zhou, D. G., Kwan, J. S. H., and Sze, H. Y. (2018). Impulse load characteristics of bouldery debris flow impact. Géotechnique Lett. 8 (2), 111–117. doi:10.1680/jgele.17.00159

CrossRef Full Text | Google Scholar

Vagnon, F. (2020). Design of active debris flow mitigation measures: a comprehensive analysis of existing impact models. Landslides 17 (2), 313–333. doi:10.1007/s10346-019-01278-5

CrossRef Full Text | Google Scholar

Wang, S. (2012). Study on failure modes of concrete gravity dams under underwater explosion impact loads based on FEM/SPH method. Tianjin University. Master’s Thesis.

Google Scholar

Wang, D. P., Chen, Z., and He, S. M. (2019). Physical model test on dynamic interaction between debris flow and bridge pier. Rock Soil Mech. 40 (09), 3363–3372. doi:10.16285/j.rsm.2018.1011

CrossRef Full Text | Google Scholar

Wang, F., Wang, J., Chen, X., Qiu, H., Lou, C., and Li, Y. (2025). Experimental study on the coupling motion mechanism of boulders in debris flow. Landslides 22 (5), 1–23. doi:10.1007/s10346-024-02450-2

CrossRef Full Text | Google Scholar

Wei, W., Guang qi, C., Zheng, H., Zhou, S., Zhang, H., and Jing, P. (2016). 3D numerical simulation of debris-flow motion using SPH method incorporating non-Newtonian fluid behavior. Nat. Hazards 81 (3), 1981–1998. doi:10.1007/s11069-016-2171-x

CrossRef Full Text | Google Scholar

Wendeler, C., Volkwein, A., McArdell, B. W., and Bartelt, P. (2019). Load model for designing flexible steel barriers for debris flow mitigation. Can. Geotechnical J. 56 (6), 893–910. doi:10.1139/cgj-2016-0157

CrossRef Full Text | Google Scholar

Xiong, H., Hao, M., Zhao, D., Qiu, Y., and Chen, X. (2023). Study of the dynamics of water-enriched debris flow and its impact on slit-type barriers by a modified SPH–DEM coupling approach. Acta Geotech. 19 (2), 1019–1045. doi:10.1007/s11440-023-02106-w

CrossRef Full Text | Google Scholar

Yao, C., Xu, C., Zhou, X., Liu, Q., and Qiang, B. (2024). Study on the destruction process of piers by debris flow impact using SPH-FEM adaptive coupling method. KSCE J. Civ. Eng. 28 (8), 3162–3175. doi:10.1007/s12205-024-1233-y

CrossRef Full Text | Google Scholar

Yongjie, Z., Yuangang, M., Gang, L., Shen, W., Gao, G., Zhao, M., et al. (2024). Influence of barrier shape on impact dynamics of debris flow entraining a boulder onto rigid barriers. J. Mt. Sci. 21 (12), 3971–3985. doi:10.1007/s11629-024-8733-6

CrossRef Full Text | Google Scholar

Zhang, T. L., Shi, H. H., Wang, C., Dong, R. L., Jia, H. X., Zhang, X. N., et al. (2011). Aerodynamic characteristics of solid particles’ acceleration by shock waves. Shock Waves 21 (3), 243–252. doi:10.1007/s00193-011-0317-z

CrossRef Full Text | Google Scholar

Zhao, L., He, W. J., Yu, X. Z., Liu, Y. P., Zhou, Z. H., and Chan, S. L. (2020). Coupled numerical simulation of a flexible barrier impacted by debris flow with boulders in front. Landslides 17 (12), 1–14. doi:10.1007/s10346-020-01463-x

CrossRef Full Text | Google Scholar

Zhou, G. G. D., Song, D., Choi, C. E., Pasuto, A., Sun, Q. C., and Dai, D. F. (2018). Surge impact behavior of granular flows: effects of water content. Landslides 15, 695–709. doi:10.1007/s10346-017-0908-6

CrossRef Full Text | Google Scholar