Debris flow surface SPH particle adaptive processing method and debris flow real-time simulation method

Abstract

The invention discloses a debris flow surface SPH particle adaptive processing method and a debris flow real-time simulation method. The invention provides a debris flow surface SPH particle adaptive processing method for overcoming the defect that surface particles are easily misjudged in the prior art. A particle density dynamic threshold is used as a surface particle identification standard, and the optimized dynamic threshold dynamically balances and controls debris flow characteristics to make a judgment contribution to the surface particles by using a real-time velocity gradient, a density gradient and a weight. Neighborhood virtual particles are randomly generated around the mass center of the surface particles to support kernel calculation, and the density and speed attributes of the virtual particles are limited. The debris flow real-time simulation method is realized by utilizing a self-adaptive processing technology. According to the optimization scheme, SPH-DEM coupling is adopted, heterogeneous units form a unified SPH framework, hydrodynamic force is calculated through an SPH method, contact force is calculated through a DEM method, and bidirectional coupling force calculation is completed. An SPH method is realized by adopting a Niagara particle system, and a digital twinborn platform foundation of numerical simulation and visual integrated simulation is constructed.

Description

Technical Field

[0001] This invention relates to an adaptive processing method for SPH particles on the surface of debris flows, a real-time simulation method for debris flows, and particularly to an adaptive processing method for surface particles in dynamic simulation of debris flows based on an SPH model, as well as a dynamic simulation method for debris flows. It belongs to the fields of electronic digital data processing and computer-aided design geological hazard simulation technology. Background Technology

[0002] Smoothed particle hydrodynamics (SPH) discretizes the object of study into a system composed of a finite number of particles. These particles carry the material properties of the object and provide computational support to simulate continuous fields. By weighted superposition of the properties of all particles in the neighborhood of the target particle, the field function and derivative are estimated. The motion process of the system is then simulated by solving the corresponding fluid equations. SPH methods are widely used in numerical simulation studies of various fluids.

[0003] Limited by the original problem domain targeted by the method development, the SPH method has inherent defects in boundary condition handling. The core idea of ​​SPH is the particle approximation within the neighborhood, including kernel interpolation approximation and particle approximation. Therefore, in the fluid boundary part, due to the lack of boundary particles, there are insufficient neighborhood particles required for the kernel function, the kernel support domain is truncated, the interpolation integral is incomplete, and the physical quantity calculation is distorted, resulting in abnormal boundary behavior. Fluid boundary problems involve two typical scenarios: fluid-solid boundary scenarios and fluid free surfaces. Since the SPH method is a meshless particle method, when simulating problems involving complex phenomena such as large deformation free surfaces, overlapping and breaking liquid surfaces, and moving boundaries, it is not necessary to introduce complex liquid surface tracking and other computational techniques to capture these complex physical phenomena. It only needs to automatically capture the free surface boundary. Therefore, among various SPH modeling techniques, there are few methods specifically designed to solve the identification of free surface particles. It is mainly suitable for modeling fluids with relatively stable real material components. The existing technology "Numerical Simulation of Bedload Motion Process Based on Lagrange Method" (Zheng Xiaogang, Sichuan University, 2021) discloses two methods commonly used to detect and identify particles on fluid free surfaces. One method is based on particle density, where a particle is considered a free surface particle if its density is less than a set value. The other method is based on the number of particles surrounding a fluid particle, where a particle is considered a free surface particle if the number of particles surrounding it is less than a set value.

[0004] Compared to the aforementioned stable fluids, debris flows, as typical mixed fluids with complex and diverse material components, exhibit frequent and highly varied collisions between particles, resulting in drastic fluctuations in local density and velocity fields. When particle aggregation (such as tumbling and accumulation) occurs in the free surface region, on the one hand, the particle density and quantity may not decrease significantly, easily leading to free surface particles being misidentified as internal fluid particles; on the other hand, short-term decreases in density and quantity caused by internal disturbances may also lead to misidentification as free surface particles. Therefore, existing technologies are not feasible for identifying free surface particles in debris flows based solely on "local particle information anomalies" such as quantity and density, limiting the value of debris flow SPH models in numerical simulation or visualization. The existing technology "Debris Flow Motion Simulation Based on SPH Method" (Wang Yaru, North China University of Water Resources and Electric Power, 2018) discloses a method for meshing the problem domain in a debris flow SPH model to find fluid surface particles. Specifically, for each SPH fluid particle, a hash query is performed on the particles in the surrounding grid. If no particle exists in a certain surrounding grid, it is marked as a particle on the surface of the debris flow slurry; otherwise, it is marked as a particle inside the debris flow slurry. This method belongs to the second commonly used detection and identification method mentioned above.

[0005] The aforementioned existing technologies are more suitable for modeling fluids with relatively stable real-world material components. In contrast, debris flows, as typical mixed fluids with complex and diverse material components, exhibit frequent and varied collisions between particles, resulting in drastic fluctuations in local density and velocity fields. When particle aggregation (such as tumbling and accumulation) occurs in the free surface region, on the one hand, the particle density and quantity may not decrease significantly, easily leading to free surface particles being misidentified as internal fluid particles; on the other hand, short-term decreases in density and quantity caused by internal disturbances may also lead to misidentification as free surface particles. Therefore, existing technologies cannot effectively identify free surface particles in debris flows solely based on "local particle information anomalies" such as quantity and density, limiting the value of debris flow SPH models in numerical simulation or visualization. Summary of the Invention

[0006] The purpose of this invention is to address the shortcomings of existing technologies by providing an adaptive processing method for surface SPH particles in debris flow dynamic simulation, and a debris flow simulation scheme implemented using this method.

[0007] To achieve the above objectives, the present invention first provides the following adaptive treatment method for SPH particles on debris flow surfaces.

[0008] An adaptive treatment method for SPH particles on debris flow surfaces. The parameters of the SPH model for debris flow were set using the smoothed particle hydrodynamics method, and the model was initialized. For any SPH particle i If the particle density is satisfied ρlocal ( i , t Particle density dynamic threshold ρ thd ( i , t ), determine the particle i If the surface of the debris flow is SPH particles, then the interior of the debris flow is considered to be SPH particles. The particle density dynamic threshold ρ thd ( i , t According to the model in Equation 1, Formula 1 In the formula, ρ 0 – Basic density threshold N -particle i The number of neighboring particles, ρ local ( i , t )-particle i The density.

[0009] The aforementioned SPH particle adaptive processing method for debris flow surfaces uses a dynamic threshold for particle density as the standard for identifying surface particles. The dynamic threshold is determined by the average density of neighboring particles in the current time step. By using the dynamic density to reflect the dynamic characteristics of the debris flow surface, the surface particle identification is linked to the dynamic characteristics of the debris flow surface, thus overcoming the shortcomings of fixed numerical comparison in existing technologies.

[0010] This invention also provides the following optimization scheme for the above-mentioned adaptive processing method. The optimization technique improves the technical solution from multiple technical perspectives, and can be implemented independently while maintaining technical independence, or simultaneously without causing conflicts in the logical relationships between the preceding and following techniques.

[0011] Optimization 1: Construct a dynamic threshold for particle density. ρ thd ( i , t A refined function model of ).

[0012] In Equation 1 model, ρ thd ( i , t ) is limited to the basic density threshold ρ The 0.1.2 times function is an empirical function summarized from theoretical analysis and a large amount of previous experimental data. Based on the data analysis of previous research results, under the condition of satisfying... ρ thd (i , t )≤1.2 ρ Under condition 0, the system can avoid misidentifying non-surface particles in high-density perturbation regions, thus effectively identifying surface particles while ensuring numerical stability. However, using a constant ratio of 1.2 at different time steps under some model conditions can lead to… ρ thd ( i , t The particle fit for specific coordinate positions is not precise, resulting in some distortion of the captured surface particles. The optimization scheme further introduces particle density and velocity-related variables to construct... ρ thd ( i , t The refined function model of ) is specifically expressed according to model 2.

[0013] Formula 2 In the formula, -particle i Real-time density gradient, -particle i Real-time velocity gradient, - Dimensionless particle velocity gradient normalized value, ranging from 0 to 1. - Dimensionless normalized particle density gradient, ranging from 0 to 1. λ v and λ ρ - These are the velocity gradient weighting coefficient and the density gradient weighting coefficient, respectively.

[0014] In Equation 2, the velocity gradient weighting coefficients λ v This is used to measure the impact of local velocity field changes on boundary thresholds. When the velocity gradient is large, the local particle motion state changes rapidly, making them more likely to appear on free surfaces or boundary regions. Increasing their weight in the determination helps improve the sensitivity of dynamic boundary recognition. Density gradient weight coefficient λ ρ This is used to measure the impact of local density field changes on the boundary determination threshold. A large density gradient indicates sparse local particle distribution or drastic density changes, which is usually related to the physical properties of the boundary region. Increasing its weight helps identify surface particles at locations of abrupt density changes. Based on experimental results, the following settings are configured: λ v and λ ρThe value of 0.1 primarily reflects two technical advantages: First, it achieves a better balance between the contributions of velocity gradient and density gradient features to surface determination. On the one hand, it ensures the velocity gradient's sensitivity to rapid changes in the free surface, and on the other hand, it maintains the role of density gradient in surface recognition, thereby improving the overall robustness of the determination. Second, it reduces computational complexity.

[0015] The optimized model in Equation 2 is constructed by simultaneously introducing the normalization conditions of particle velocity gradient and density gradient. ρ thd ( i , t Relative to ρ The dynamic ratio relationship of 0 is used, and the changing characteristics of the velocity field and density field are integrated into the boundary judgment threshold by using weighting coefficients to achieve coordinated control of multi-physics gradients. This can reduce the misjudgment problem of relying on a single index or fixed value, so that the control conditions are dynamically satisfied for the calculation of any particle at any time step in the numerical model process. ρ thd ( i , t )≤1.2 ρ 0, to achieve dynamic equilibrium of the numerical process, as well as the stability and accuracy of surface particle capture.

[0016] Optimization 2: For any SPH particle, calculate the particle density using the properties of its neighboring particles according to Equation 3. ρ local ( i , t ).

[0017] Formula 3 In the formula, N -particle i Neighborhood particles j quantity, m j - Neighborhood particles j quality, h - Smooth kernel radius, and - These are particles i With particles j The position vector.

[0018] Optimization 3: Adaptive dynamic compensation for surface SPH particles after identification.

[0019] At any time step of the iterative cycle, all particles in the SPH model are traversed, and SPH particles on the debris flow surface and those inside the debris flow are identified, and sets of surface SPH particles are constructed respectively. T sur Internal SPH particle collection T int ; Compute set T int The average number of neighboring particles of all particles n int ; For sets T sur For any particle, calculate the number of its neighborhood particles. n real ,Will n int and n real The difference is denoted as the number of neighborhood virtual particles that need to be added within its support domain. n need ; Calculate the centroid coordinates of the support domain Randomly generate a number around the centroid. n need The virtual particle supporting kernel function is calculated, and the neighborhood virtual particle properties satisfy the density... ρ tem ,speed Once the resultant force calculation is complete, the system clears the neighboring virtual particles.

[0020] Formula 4 Formula 5 Formula 6 In the formula, N ( i )-particle i The set of neighborhood particles, -particle j The position vector, ρ local ( j )-particle j density, W ( r ij ,h ) - Smoothing kernel function, -particle j The velocity vector.

[0021] The reason for setting the neighborhood virtual particles used to support the surface particle kernel computation around the centroid is that the offset of the centroid relative to the particle's own position can effectively reveal the main direction of the particle's support domain deficiency. In other words, the centroid position can show the completeness and asymmetry of the current particle support domain. Therefore, the process of supplementing the neighborhood virtual particles around the centroid is actually a benchmark dynamic adjustment of the support domain, effectively and quickly overcoming the impact of support domain defects on kernel computation.

[0022] The aforementioned adaptive SPH particle processing method for debris flow surfaces can be applied to debris flow simulations using SPH-based numerical models. The SPH-based numerical model can be either a debris flow model constructed using the SPH method or a debris flow model constructed using the SPH-DEM coupling method. In this model, the SPH method simulates the continuous fluid phase of the debris flow, while the Discrete Element Method (DEM method) simulates the discrete solid phase.

[0023] Based on the above-mentioned adaptive processing method for SPH particles on the surface of debris flow, this invention also provides a real-time simulation method for debris flow using this technical solution, as follows.

[0024] A real-time debris flow simulation method is provided, which involves constructing a debris flow numerical model based on the SPH method. The SPH-based debris flow numerical model is a debris flow model constructed using the SPH method, and then implemented as follows: Step S100: Define global parameters for the model, create all model micro-units and set basic properties according to computational needs, and register all model micro-units into a unified SPH framework. The model micro-units are SPH particles. Step S200: Construct the global acceleration structure of the SPH framework, perform neighborhood particle search, and complete model initialization; Step S300 includes: Step S310: Calculate the density and pressure of all model micro-units; Step S320: The free interface of the debris flow is treated using the above-mentioned SPH particle adaptive treatment method on the debris flow surface. Step S330: Calculate the resultant force of the model micro-element excluding gravity; Step S340: Calculate the acceleration of all model micro-elements; Step S350: Update the state of all model micro-units; Step S400: Perform collision detection and processing; Step S500: Check whether all model micro-units exceed the range of the acceleration structure. If they do not exceed the range, proceed to the next time step for calculation; otherwise, end the simulation calculation.

[0025] This invention also provides an optimized technical solution for the aforementioned real-time debris flow simulation method. The optimization technique improves the technical solution from multiple technical perspectives, and can be implemented independently while maintaining technical independence, or simultaneously without causing logical conflicts between the preceding and following techniques.

[0026] Optimization 1: Optimization of the debris flow numerical model based on the SPH method.

[0027] A debris flow numerical model based on the SPH method is constructed using the SPH-DEM coupling method. The SPH method simulates the continuous fluid phase in the debris flow, while the DEM method simulates the discrete solid phase.

[0028] In the coupled model, the composition of the mud-rock fluid is simulated by both SPH particles and DEM particles. Therefore, the model micro-units registered in the unified SPH framework are either SPH particles or DEM particles. When calculating the forces between model micro-units, the SPH method is used to calculate the hydrodynamic forces (such as gradient forces and viscous forces) exerted by the fluid phase between model micro-units, while the DEM method is used to calculate the contact forces exerted by the solid phase between model micro-units. These two forces together constitute the resultant force of each model micro-unit, excluding gravity, and are used for acceleration calculations. Simultaneously, the DEM method is also used to calculate the angular acceleration of the DEM particles. Based on the calculation results, the state updates of the model micro-units include the velocity and position of SPH particles, and the velocity, position, and rotation of DEM particles.

[0029] The optimized scheme utilizes the SPH-DEM coupling method to consider the contact forces between debris flow components in the debris flow model. The interaction between components is achieved through both hydrodynamics and contact forces. By uniformly registering the heterogeneous model micro-elements of SPH particles and DEM particles into the global SPH framework, a computational platform for bidirectional coupling forces is obtained. The mechanical implementation problem of the bidirectional coupling model is solved by using the SPH method to calculate hydrodynamics and the DEM method to calculate contact forces. When implementing the SPH particle adaptive processing method on the debris flow surface to process the free interface of the debris flow in step S330, the fluid particles within the global SPH framework (i.e., all SPH particles in the model micro-elements) are used as the target particles. i Each particle is then individually determined to be a surface particle. During the calculation, the target particles... i SPH particles and DEM particles within the support domain are jointly labeled as their neighboring particles. j That is, quantity N It also includes SPH particles and DEM particles within the support domain.

[0030] Optimization 2 involves optimizing the calculation of contact forces between micro-units in the model.

[0031] In the step of calculating the contact force between model micro-units using the DEM method, a calculation concept is provided to "treat the two heterogeneous model micro-units as homogeneous members". The specific scheme is as follows: When calculating the contact force between DEM particles, the conventional DEM method is used; when calculating the contact force between SPH particles, for the target SPH particle, the candidate SPH particles for calculating the overlap amount are all SPH particles within the support domain of the SPH model, and then the conventional DEM method is used to calculate based on the SPH particle attributes; when calculating the contact force between DEM particles and SPH particles, the position coordinates of the target SPH particle in the DEM model space are first calculated, and the DEM particles whose position coordinates are within the search radius of the DEM model are marked as candidate DEM particles for calculating the overlap amount of the target SPH particle, and then the conventional DEM method is used to calculate based on the attributes of the SPH particles and DEM particles; when calculating the contact force between DEM particles and SPH particles, the position coordinates of the target DEM particle in the SPH model space are first calculated, and all SPH particles whose position coordinates are within the support domain of the SPH model are marked as candidate SPH particles for calculating the overlap amount of the target DEM particle, and then the conventional DEM method is used to calculate based on the attributes of the DEM particles and SPH particles.

[0032] Optimization 3 is implemented within the Niagara particle system of Unreal Engine (UE).

[0033] Because the real-time debris flow simulation method employs a global SPH framework with "global particleization properties," it is combined with the Niagara particle system of Unreal Engine. The Niagara system is then used to simulate SPH fluids (independent SPH fluids or SPH-DEM coupled fluids), resulting in a debris flow visualization simulation framework that integrates numerical simulation and real-time vision. This framework possesses high scalability and is an ideal digital twin platform to support data-driven debris flow monitoring and early warning schemes.

[0034] Compared with the prior art, the beneficial effects of the present invention are: (1) The prior art only relies on the basic technical setting of the SPH method for the treatment of debris flow free surface, and adopts the method of comparing the number of particles around the fluid particles to identify free surface particles. These methods are not very adaptable to the typical free surface motion characteristics of debris flow and are prone to misjudgment. The present invention provides an adaptive processing method for SPH particles on debris flow surface to solve two levels of technical problems. The first level is the problem of adaptive identification of debris flow free surface particles. The method uses a dynamic threshold of particle density as the standard for surface particle identification. The dynamic threshold is determined by the average density of neighboring particles in the current time step, thereby linking surface particle identification with the dynamic fluid characteristics of debris flow and overcoming the defect of fixed numerical comparison in the prior art. In the optimization scheme, the method constructs a particle density dynamic threshold model jointly expressed by real-time velocity gradient and real-time density gradient, and based on the analysis of previous experimental data, uses a two-weighted coefficient to realize the collaborative dynamic balance of multi-physics gradient to control the contribution of debris flow characteristics to surface judgment, thereby improving the stability and accuracy of surface particle capture. The second level is the problem of adaptive dynamic compensation of surface SPH particles after identification. For the identified free surface particles, the method randomly generates neighborhood virtual particle support kernels around their centroids to calculate compensation properties. After the resultant force calculation step is completed, the system automatically clears the neighborhood virtual particles. The optimization scheme also specifically designs the density and velocity properties of the neighborhood virtual particles. (2) The debris flow real-time simulation method of the present invention provides a debris flow real-time simulation technology realized by the SPH particle adaptive processing method on the debris flow surface. In the optimization scheme, a debris flow dynamic simulation method based on the SPH-DEM coupling method is adopted. The unified SPH framework composed of SPH particles and DEM particle heterogeneous model micro-units is used as a two-way coupled calculation platform. The SPH method is used to calculate hydrodynamics and the DEM method is used to calculate contact force to solve the mechanical realization problem of the two-way coupled model. Under the coupling framework, the surface SPH particle adaptive processing is performed on the fluid part, and the solid particles are jointly identified as the neighborhood particles of the fluid particles. In the further optimization scheme, for different types of contact forces between SPH particles and DEM particles, a method is designed to find candidate micro-units for contact force calculation in the SPH support domain or within the DEM search radius, respectively, thereby optimizing the coupling effect of mechanical calculation under the unified SPH framework. The method enhances the technical value of the surface adaptive processing method of this invention in debris flow simulation models. (3) It provides a method for realizing SPH fluid (independent SPH fluid or SPH-DEM coupled fluid) by using the Niagara particle system based on Unreal Engine, realizing a real-time simulation scheme that integrates numerical simulation and visualization, which is easy to expand the function and is an ideal debris flow digital twin platform in the field of data-driven debris flow monitoring and early warning. Attached Figure Description

[0035] Figure 1This is a schematic diagram of the main steps in the SPH particle adaptive treatment process for debris flow surfaces.

[0036] Figure 2 This represents the free surface and abnormal flow state of the debris flow SPH model.

[0037] Figure 3 Here are examples of local surface particle determination in the debris flow SPH model: (a) shows the existing fixed threshold method, and (b) shows the adaptive dynamic threshold method of the present invention.

[0038] Figure 4 Here are examples of determining surface particles in the debris flow tumble zone using the SPH model for debris flows. (a) shows the existing fixed threshold method, and (b) shows the adaptive dynamic threshold method of the present invention.

[0039] Figure 5 This is a diagram showing the surface particle compensation effect of the debris flow SPH model.

[0040] Figure 6 The diagrams are of the simulated dam-break experiment setup: (a) a side view of the setup in this example, (b) a front view of the setup in this example, (c) a side view of the setup in Reference 1, and (d) a front view of the setup in Reference 1.

[0041] Figure 7 This is a schematic diagram of the SPH-DEM coupling mode.

[0042] Figure 8 This is a schematic diagram of the main implementation steps of the SPH-DEM coupled numerical simulation of debris flow.

[0043] Figure 9 This is a diagram showing the overall flow state of debris flow at four observation points during the first simulated dam break test of the experimental group.

[0044] Figure 10 This is a comparison of the displacement of the debris flow front in simulated dam break and physical dam break experiments. Detailed Implementation

[0045] The preferred embodiments of the present invention will now be further described with reference to the accompanying drawings.

[0046] Figure 1 This is a schematic diagram of the main steps in the SPH particle adaptive treatment process for debris flow surfaces.

[0047] Example 1

[0048] The method of this invention is used to perform adaptive processing of surface SPH particles in a debris flow dynamic simulation model, and to dynamically and adaptively determine the particle type.

[0049] 1. Debris Flow Numerical Model Setup and Initialization A numerical model of debris flow was constructed using the Smoothed Particle Hydrodynamics (SPH) method. The SPH time integration employed a semi-implicit Euler method, where the velocity was updated first, and then the particle position was calculated using the updated velocity. This method is characterized by its simple computational process, ease of parallel implementation on GPUs, and good numerical stability under relatively small time steps. The physical parameters of the SPH model were comprehensively referenced from existing technologies and previous research experience; the main parameters are shown in Table 1.

[0050] Table 1. Main settings parameters for the debris flow SPH model

[0051] Existing technical document 1: Xiaosong Sun et al , Three-dimensional simulation of asolid-liquid flow by the DEM-SPH method. Journal of Computational Physics ,2013, 248:147–176 (Paper title translation: Three-dimensional simulation of solid-liquid two-phase flow based on DEM-SPH method).

[0052] Existing technical document 2: Chun Liu et al , Consideration of maximum impact forcedesign for a rock shed against dry granular flow. European Journal of Environmental and Civil Engineering , 2022, 9:1–22 (Paper title translation: Design considerations for maximum impact force of dry debris flow rockfall shed).

[0053] Existing technical literature 3: Wang Yaru, Simulation of debris flow motion based on SPH method, North China University of Water Resources and Electric Power, Master's thesis, 2018.

[0054] Complete model initialization.

[0055] In this example, the Niagara particle system of Unreal Engine is used to build a debris flow simulation model, which is then implemented in parallel on a GPU.

[0056] After initializing the SPH model, the standard numerical simulation process is followed, including iterative calculations such as neighborhood search, density calculation, correlation force calculation, time integration, and particle update, as well as boundary treatment and output analysis. Figure 1 This example omits the description of routine steps and focuses on the relevant operational details of the surface SPH particle adaptive treatment method of this invention.

[0057] 2. Determination of SPH Particle Type in Debris Flow In the iterative calculation, at any time step, for any SPH particle i Perform a neighborhood search to establish its neighborhood particle set. Utilize the properties of the neighborhood particles to perform the following calculations.

[0058] The particle density was calculated using the Equation 3 model. ρ local ( i , t ).

[0059] If the particle density is satisfied ρ local ( i , t Particle density dynamic threshold ρ thd ( i , t If the particle is identified as an SPH particle on the surface of the debris flow, it is determined to be an SPH particle inside the debris flow.

[0060] As a preferred implementation scheme, in this example, ρ thd ( i , t The determination of ) adopts an optimization scheme, namely, the calculation of the refined dynamic model expressed in Equation 2.

[0061] To obtain the normalized values ​​of the dimensionless particle velocity gradient and the dimensionless particle density gradient, the particle density gradient must first be calculated separately. velocity gradient .

[0062] Formula 7 Formula 8 In time step t =0.2s particle calculation i 1. i For example, the relevant data for adaptive particle type determination are shown in Table 2.

[0063] Table 2 Examples of Adaptive Particle Type Determination

[0064] We extract a local area from the SPH model simulation domain of debris flow as an example. The particle emission source is placed in the central region of the scene container. Initially, the particles are concentrated at the center, and then diffuse outwards under the action of gradient forces. Upon encountering the container boundary, they collide and bounce, forming a surge and generating backflow, ultimately causing the particles in the central region to re-aggregate, resulting in a significant increase in local density.

[0065] Figure 2 This represents the free surface and abnormal flow state of the debris flow SPH model. Figure 2 The diagram shows that due to the lack of particle density on the free surface, particles at the bottom surge upwards under the influence of gradient forces, while particles at the center diffuse outwards. The areas highlighted in the image indicate locations where the debris flow exhibits abnormal flow patterns.

[0066] Figure 3 Here are examples of local surface particle determination in the debris flow SPH model: (a) shows the existing fixed threshold method, and (b) shows the adaptive dynamic threshold method of the present invention. Figure 3 The pink particles are surface SPH particles, and the white particles are internal SPH particles. Figure 3 (a) shows the result of the prior art fixed threshold method (for specific operation of the prior art fixed threshold, refer to formula (3-21) in prior art document 4): ρ i <α ρ 0, ρ i For any i Particle density, ρ 0 represents the initial density, and α is a fixed threshold coefficient typically set between 0.8 and 0.96. In this example, the fixed threshold coefficient is set to 0.95. This shows that the fixed threshold method can effectively identify surface particles in the diffusion stage, but when central particles aggregate, because the density does not fall below the fixed threshold, the actual surface particles are not mistakenly identified as internal SPH particles. The results of using the adaptive dynamic threshold method of this invention... Figure 3 As shown in Figure (b), the use of a local dynamic adaptive threshold based on real-time density gradient and velocity gradient effectively improves the accuracy of determining surface particles in particle aggregation regions and avoids simulation distortion.

[0067] Existing technical literature 4: Zheng Xiaogang, Numerical simulation of bedload motion process based on Lagrange method, Sichuan University, Master's thesis, 2021.

[0068] Figure 4 Here are examples of determining surface particles in the debris flow tumble zone using the SPH model for debris flows. (a) shows the existing fixed threshold method, and (b) shows the adaptive dynamic threshold method of the present invention. Figure 4 The meanings of pink and white particles Figure 3 The same applies. In the model, when particles collide with the container boundary and cause turbulence, particle aggregation is likely to occur in the central turbulent region. In this case: Figure 4 As shown in (a), the existing technology cannot accurately identify the debris flow as surface particles because the density of particles on the surface of the debris flow in the turbulent area does not decrease significantly. Figure 4As shown in Figure (b), even when the surface particle density is still high, the threshold dynamic balance control method of the present invention can still accurately determine that it is a surface particle.

[0069] Example 2 The method of this invention is used to adaptively process surface SPH particles in a debris flow dynamic simulation model. Specifically, it is based on Example 1, but with adaptive dynamic compensation for surface SPH particles. The similarities to Example 1 will not be repeated, and the description of the conventional steps of the SPH model will be omitted. The focus is on describing the relevant operational content of the adaptive surface SPH particle processing method of this invention.

[0070] 1. Construct a collection of type particles At any given time step, traverse all particles in the SPH model and construct a type particle set based on the particle determination results from Example 1, including the surface SPH particle set. T sur Internal SPH particle collection T int .

[0071] 2. Randomly generate neighborhood virtual particles At any given time step, execute: Compute set T int The average number of neighboring particles of all particles n int .

[0072] Formula 9 For sets T sur For any particle, calculate the number of its neighborhood particles. n real ,Will n int and n real The difference is denoted as the number of neighborhood virtual particles that need to be added within its support domain. n need ; Calculate the centroid coordinates of the support domain (Equation 4) Randomly generate a number around the centroid. n need Virtual particle support kernel function calculation; neighborhood virtual particle properties satisfy density ρ tem (Equation 5) Speed (Formula 6).

[0073] Once the resultant force calculation is complete, the system clears the neighboring virtual particles.

[0074] In time step tTaking a calculation of 0.2s as an example, the SPH global total consists of 3928 particles, and the set T sur With sets T int The number of particles are 236 and 3692 respectively. For the set... T int , n int =21.3. (Based on sets) T sur SPH particles on any surface i Taking 3 as an example, its neighborhood particle number n real =14, and a number of random additions need to be made around its centroid coordinates. n need =21.3-14=7.3≈7 (rounded) The neighborhood virtual particles are used as the supporting kernel function for calculation. i The attribute data of the 14 neighboring particles of 3 are shown in Table 3.

[0075] Table 3 Set T sur Neutron i 3's 14 neighborhood particle attribute data

[0076] Utilizing neighboring particles j Attribute data interpolation calculation determined, particles i 3. Centroid coordinates ≈(4.970, 4.894, 4.975), the properties of randomly generated neighborhood virtual particles around the centroid should satisfy density. ρ tem =1497.9kg / m 3 ,speed = (2.01, 2.01, 2.02) m / s.

[0077] Figure 5 This is a diagram showing the surface particle compensation effect of the debris flow SPH model. Figure 5 and Figure 2 The comparison shows that, after compensating for the lack of SPH particles on the surface, the fluid avoids unstable phenomena such as central particles surging to the surroundings and bottom particles surging to the top, and the system can effectively simulate the stable state of the fluid.

[0078] Example 3 A debris flow simulation model was built using the Niagara particle system in Unreal Engine to simulate the debris flow dam breakage process. Two tasks were accomplished: First, the SPH particle adaptive processing method of the debris flow surface was applied to the simulation model of the complete debris flow movement process, and the effectiveness of the application scheme was verified through experimental data. Second, real-time rendering was superimposed on the debris flow numerical simulation model in the Niagara particle system to achieve real-time visual effects in the complete debris flow movement.

[0079] 1. Simulation Testing System 1.1 Simulated Dam Failure Experimental Device Existing technology document 1 discloses a solid-liquid two-phase debris flow dam-break experimental device built in a real environment, and uses experimental data to verify the effectiveness of the constructed debris flow numerical model. In this example, the experimental research system (hardware and testing process) of document 1 is reproduced in a virtual environment using Unreal Engine, resulting in a virtual testing system for debris flow dam-break experiments (hereinafter referred to as the simulated dam-break system). Then, the same tests as in document 1 are carried out in the simulated dam-break system, and the simulation test results are compared with its experimental data to verify the accuracy and reliability of the simulated dam-break system.

[0080] A simulated dam-break experimental device was built in the Niagara virtual system based on the content disclosed in existing technical document 1. Figure 6 These are diagrams of a simulated dam-break experiment setup: (a) a side view of the setup in this example, (b) a front view of the setup in this example, (c) a side view of the setup in Reference 1, and (d) a front view of the setup in Reference 1. In Reference 1, the setup uses a = 5 cm as the unit, and the complete domain size is set to 4a × 3a × 2a. The simulated dam-break experiment setup built in this example maps this size using the number of grid cells as the proportional unit of the virtual space.

[0081] 1.2 Setting up the SPH-DEM numerical model for debris flow As a preferred embodiment, a coupled SPH-DEM numerical model of debris flow is constructed by combining the Smoothed Particle Hydrodynamics (SPH) method with the Discrete Element Method (DEM). The SPH method simulates the continuous fluid phase in debris flow, characterizing the flow and wave motion of the fluid. The DEM method simulates the discrete solid phase in debris flow, and is suitable for solid particles with obvious particle characteristics and independent movement, such as boulders and gravel, in debris flow.

[0082] Figure 7This is a schematic diagram of the SPH-DEM coupling mode. In the left figure, the curved arrows indicate the connecting fluid phase, and the solid circles indicate the discrete solid phase. The coupling calculation of SPH and DEM is achieved through the bidirectional mechanical interaction between the continuous fluid phase and the solid phase. DEM particles are introduced into the SPH particle system, forming a unified SPH framework system composed of heterogeneous model micro-units (DEM particles and SPH particles). The SPH method is used to calculate the hydrodynamic forces between model micro-units, and the DEM method is used to calculate the contact forces between model micro-units. That is, SPH particles and DEM particles jointly participate in the hydrodynamic forces (gradient forces and viscous forces) calculated by the SPH kernel function, and also jointly participate in the contact forces calculated by the DEM method.

[0083] For the calculation of contact forces between micro-units in the model, the following distinctions can be made based on different properties of the contact forces.

[0084] (1) Calculate the contact force between solid particles f SScol At that time, the contact force between DEM particles was calculated using the conventional DEM method.

[0085] (2) Calculate the contact force between fluid particles f FFcol When calculating the overlap of a target SPH particle, the candidate SPH particles are all SPH particles in the SPH model support domain, and then the overlap is calculated using the conventional DEM method based on the SPH particle properties.

[0086] (3) Calculate the contact force between the fluid particles and the solid particles. f SFcol First, the position coordinates of the target SPH particle in the DEM model space are calculated, and the DEM particles whose position coordinates are within the search radius of the DEM model are recorded as candidate DEM particles for calculating the overlap of the target SPH particles. Then, the conventional DEM method is used to calculate the overlap based on the properties of the SPH particles and DEM particles.

[0087] (4) Calculate the contact force between the solid particles and the fluid particles. f FScol First, the position coordinates of the target DEM particles in the SPH model space are calculated, and all SPH particles in the SPH model support domain with these position coordinates are marked as candidate SPH particles for calculating the overlap of the target DEM particles. Then, the conventional DEM method is used to calculate the overlap based on the properties of the DEM particles and SPH particles.

[0088] Considering the mechanical properties of the debris flow material components, this embodiment does not consider the contact force between fluid particles and solid particles. The force calculation for the model micro-element is as follows: F f The net force acting on the fluid particles (SPH particles) includes the gradient force between the fluid particles.f FFpre Viscosity f FFvis Contact force f FFcol and gradient forces from solid particles. f SFpre Viscosity f SFvis Contact force f SFcol (Equation 10); F s The resultant force on the solid particles (DEM particles) includes the gradient force from the continuous fluid. f FSpre Viscosity f FSvis and the contact force between solid particles f SScol (Equation 11).

[0089] Formula 10 Formula 11 Contact force includes normal contact force F n tangential contact force F t Its size | F n |、| F t Calculate according to Equations 12 and 13 respectively, where, δ , k , η , v These represent contact overlap, stiffness coefficient, damping coefficient, and relative velocity, respectively, with subscripts... n and t These represent the normal and tangential directions, respectively. μ The coefficient of friction is tangential.

[0090] Formula 12 Formula 13 The physical parameters of the SPH-DEM model are comprehensively referenced from existing technologies and previous research experience. The SPH model parameters are set the same as in Example 1, while the DEM model uses a spring-damped model. The main parameter settings for the model are shown in Table 4.

[0091] Table 4. Main settings parameters for the debris flow SPH-DEM model

[0092] 2. Numerical Realization of Debris Flow SPH-DEM Coupled Figure 8 This is a schematic diagram of the main implementation steps of the SPH-DEM coupled numerical simulation of debris flow.

[0093] DEM particles are introduced into the SPH particle system to complete the SPH-DEM model initialization, start the time-step cycle, and execute numerical simulation calculations. The following description omits the routine steps of the SPH-DEM coupled model simulation, focusing only on the main steps. Since the two types of model micro-units participate in numerical calculations with equal status after the DEM particles are introduced into the SPH particle system, the three types of computational objects are referred to as SPH particles, DEM particles, and collectively as model micro-units.

[0094] For ease of description, the main steps of coupled numerical implementation are represented as four stages: Phase I: Model Definition Define global parameters of the model, such as computational domain, smoothness length, and time step; Create all model micro-units and set basic properties according to computational needs; Register all model micro-units into a unified SPH framework.

[0095] Phase II: Model Initialization Constructing a global acceleration structure for the SPH framework: An acceleration structure is built based on the current positions of all model micro-units. To integrate the SPH method into the Niagara particle system, a spatial grid-type global acceleration structure is used. In this example, a spatial hash grid is chosen. Perform a neighborhood search on all model micro-units to complete model initialization.

[0096] Phase III: Calculation of Mechanical and Particle Properties Calculate the density and pressure of the micro-units in the model; For the debris flow fluid portion (i.e., taking any SPH particle as the target particle), the method of Example 2 is used to adaptively identify particles at the debris flow free interface and perform density compensation. In the calculation, SPH particles within the support domain and DEM particles are jointly labeled as target SPH particles. i Neighborhood particles j That is, the number of neighboring particles N Includes SPH particles and DEM particles within the supported domain; Calculate the resultant force of the micro-element in the model excluding gravity. F f and F s This includes: using the SPH method to calculate the gradient force and viscous force between model micro-elements, and using the DEM method to calculate the normal contact force and tangential contact force between model micro-elements; Calculate the acceleration of the micro-element in the model; Based on the tangential contact force between DEM particles, update the tangential components of its local coordinate system and calculate the angular acceleration; Update the state of the model micro-units, including: updating the velocity and position of SPH particles, and updating the velocity, position, and rotation of DEM particles; Collision Detection and Handling: If the updated position of a model micro-unit is within the scene object content, a collision is confirmed to have occurred within this time step, requiring collision handling. In this example, a restitution coefficient is introduced to correct the velocity of the model micro-unit after the collision, and based on previous research results, the restitution coefficient is set to 0.35–0.45, preferably 0.4. That is, the magnitude of the normal velocity of the colliding particle becomes 40% of the original normal velocity, the direction is reversed, while the horizontal velocity remains unchanged.

[0097] Phase IV: Simulated Domain Detection The system checks whether the micro-units of the detection model exceed the range of the acceleration structure. If they do not exceed the range, it proceeds to the next time step for calculation; otherwise, it terminates the numerical simulation.

[0098] 3. Simulated dam failure test and results Experimental and control groups were set up. The experimental group underwent adaptive processing (including identification and compensation) of surface SPH particles according to the method in Example 2, while the control group did not undergo any special processing of surface SPH particles. The following description of the routine steps of the SPH-DEM coupled model simulation, as well as the content that is the same as in Examples 1 and 2, is omitted.

[0099] The experimental group and the control group each conducted three simulated dam-break tests, and recorded the displacement values ​​of the mud and fluid front at 12 observation times from t=0.0s to t=0.2s each time (i.e., 12 observation values).

[0100] Figure 9 This is a diagram showing the overall flow state of the debris flow at four observation times during the first simulated dam break test conducted by the experimental group: t=0.05s, t=0.10s, t=0.15s, and t=0.20s. The white and blue particles in the diagram represent SPH particles and DEM particles, respectively.

[0101] The average value of each observation in the three test results of the experimental group and the control group was calculated separately. The results and average values ​​of the experimental group and the control group were then compared with the snapshot of the physical dam break test experiment in existing technical literature 1. The results are shown in [reference needed]. Figure 10 Table 5.

[0102] Figure 10This is a comparison of the displacement at the mudflow front between the simulated dam-break and physical dam-break experiments. The curves in the figure are named as follows: SR1-ADC, SR2-ADC, SR3-ADC, and ASR-ADC represent the 1st, 2nd, and 3rd simulation observations and their average values ​​after density compensation in the experimental group, respectively; SR1-BDC, SR2-BDC, SR3-BDC, and ASR-BDC represent the 1st, 2nd, and 3rd simulation observations and their average values ​​before density compensation in the experimental group, respectively; the experimental snapshot represents the control group data. Figure 10 The data is shown in Table 5.

[0103] Table 5. Debris flow front data from snapshots of simulated and physical dam-break experiments in the experimental group.

[0104] The simulation error was evaluated using the mean absolute error (MAPE) index, based on the average results of the experimental and control groups respectively. MAPE is calculated according to Equation 14, where... n The total number of samples, E i For the first snapshot of the experiment i One observation value, S i For the first i The average of the observations.

[0105] Formula 14 MAPE assessment showed that the control group had a MAPE of 12.3507% and the experimental group had a MAPE of 8.552%.

[0106] Figure 10 The data in Table 5 and the MAPE results show that: (1) The adaptive processing method for surface SPH particles in the debris flow dynamic simulation model of the present invention can be well applied to the debris flow numerical model based on the SPH method, supporting the dynamic simulation of the entire debris flow dam break movement; (2) In the debris flow numerical model based on the SPH method, by adopting the adaptive processing method for surface SPH particles in the debris flow dynamic simulation model of the present invention, the simulated dam break debris flow front displacement at each time node maintains a higher consistency with the experimental snapshot of the physical dam break. Compared with the prior art, the adaptive method of the present invention can improve the accuracy of numerical simulation and the stability of the system.

[0107] The specific embodiments of this invention employ a Niagara particle system based on Unreal Engine to realize SPH fluid (independent SPH fluid or SPH-DEM coupled fluid), achieving real-time debris flow simulation integrating numerical simulation and visualization. This technical framework is easily expandable with front-end real-time monitoring data access and overall closed-loop control functions. Combined with current real-time data processing, simulation analysis, and real-time visual presentation capabilities, it represents an ideal digital twin platform for debris flow monitoring and early warning in the data-driven field.

Claims

1. An adaptive treatment method for SPH particles on the surface of debris flows, characterized in that: The parameters of the SPH model for debris flow were set using the smoothed particle hydrodynamics method, and the model was initialized. For any SPH particle i If the particle density is satisfied ρ local ( i , t Particle density dynamic threshold ρ thd ( i , t ), determine the particle i If the surface of the debris flow is SPH particles, then the interior of the debris flow is considered to be SPH particles. The particle density dynamic threshold ρ thd ( i , t According to the model in Equation 1, Formula 1 In the formula, ρ 0 – Basic density threshold N -particle i The number of neighboring particles, ρ local ( i , t )-particle i The density.

2. The adaptive treatment method for SPH particles on the debris flow surface according to claim 1, characterized in that: The ρ thd ( i , t According to the model in Equation 2, Formula 2 In the formula, -particle i Real-time density gradient, -particle i Real-time velocity gradient, - Dimensionless particle velocity gradient normalized value, ranging from 0 to 1. - Dimensionless normalized particle density gradient, ranging from 0 to 1. λ v and λ ρ - These are the velocity gradient weighting coefficient and the density gradient weighting coefficient, respectively.

3. The adaptive treatment method for SPH particles on the debris flow surface according to claim 1, characterized in that: Particle density ρ local ( i , t ) Calculate according to Equation 3, Formula 3 In the formula, h - Smooth kernel radius, m j -particle i Neighborhood particles j quality and - These are particles i With particles j The position vector.

4. The adaptive treatment method for SPH particles on the debris flow surface according to any one of claims 1 to 3, characterized in that: At any time step of the loop iteration, execute: Iterate through all particles in the SPH model, identify SPH particles on the debris flow surface and SPH particles inside the debris flow, and construct separate sets of surface SPH particles. T sur Internal SPH particle collection T int ; Compute set T int The average number of neighboring particles of all particles n int ; For sets T sur For any particle, calculate the number of its neighborhood particles. n real ,Will n int and n real The difference is denoted as the number of neighborhood virtual particles that need to be added within its support domain. n need Calculate the centroid coordinates of the support domain. Randomly generate a number around the centroid. n need The virtual particle supporting kernel function is calculated, and the neighborhood virtual particle properties satisfy the density... ρ tem With speed , Formula 4 Formula 5 Formula 6 In the formula, N ( i )-particle i The set of neighborhood particles, -particle j The position vector, ρ local ( j )-particle j density, W ( r ij ,h - Smoothing kernel function, -particle j The velocity vector; Once the resultant force calculation is complete, the system clears the neighboring virtual particles.

5. The adaptive treatment method for SPH particles on the debris flow surface according to claim 4, characterized in that: It is applied to debris flow simulation in the numerical model of debris flow based on the SPH method.

6. A real-time simulation method for debris flow, characterized in that: A debris flow numerical model based on the SPH method is constructed, wherein the debris flow numerical model based on the SPH method is a debris flow model constructed using the SPH method, and then implemented as follows: Step S100: Define global parameters for the model, create all model micro-units and set basic properties according to computational needs, and register all model micro-units into a unified SPH framework. The model micro-units are SPH particles. Step S200: Construct the global acceleration structure of the SPH framework, perform neighborhood particle search, and complete model initialization; Step S300 includes: Step S310: Calculate the density and pressure of all model micro-units; Step S320: The debris flow free interface is treated using the SPH particle adaptive treatment method for debris flow surface as described in claim 4. Step S330: Calculate the resultant force of the model micro-element excluding gravity; Step S340: Calculate the acceleration of all model micro-elements; Step S350: Update the state of all model micro-units; Step S400: Perform collision detection and processing; Step S500: Check whether all model micro-units exceed the range of the acceleration structure. If they do not exceed the range, proceed to the next time step for calculation; otherwise, end the simulation calculation.

7. The real-time simulation method for debris flow according to claim 6, characterized in that: The debris flow numerical model based on the SPH method is a debris flow model constructed using the SPH-DEM coupling method, wherein the SPH method simulates the continuous fluid phase part of the debris flow, and the DEM method simulates the discrete solid phase part of the debris flow. In step S100, the model micro-unit is an SPH particle or a DEM particle. In step S320, debris flow free interface processing calculation is performed with any SPH particle as the object particle. In the calculation, SPH particles in the support domain of the SPH particle and DEM particles are jointly marked as their neighboring particles. In step S330, the gradient force and viscous force between model micro-units are calculated using the SPH method, and the contact force between model micro-units is calculated using the DEM method to obtain the resultant force of the model micro-units excluding gravity. In step S340, calculating the acceleration of all model micro-units includes: updating the tangential components of the local coordinate system and calculating the angular acceleration based on the tangential contact force between DEM particles; In step S350, updating the state of all model micro-units includes: updating the velocity and position of SPH particles, and updating the velocity, position, and rotation of DEM particles.

8. The real-time simulation method for debris flow according to claim 7, characterized in that: Step S330, in the calculation of the contact force between micro-elements of the model using the DEM method, When calculating the interparticle contact force in DEM, the conventional DEM method is used. When calculating the contact force between SPH particles, for the target SPH particle, the candidate SPH particles for calculating the overlap amount are all SPH particles in the SPH model support domain of the SPH body particle, and then the conventional DEM method is used to calculate based on the SPH particle properties. When calculating the contact force between SPH particles and DEM particles, first calculate the position coordinates of the target SPH particle in the DEM model space, and record the DEM particles whose position coordinates are within the search radius of the DEM model as candidate DEM particles for calculating the overlap of the target SPH particles. Then, calculate the force using the conventional DEM method based on the properties of SPH particles and DEM particles. When calculating the contact force of SPH particles on DEM particles, the position coordinates of the target DEM particles in the SPH model space are first calculated, and all SPH particles in the SPH model support domain with these position coordinates are marked as candidate SPH particles for calculating the overlap of the target DEM particles. Then, the conventional DEM method is used to calculate the overlap based on the properties of DEM particles and SPH particles.

9. The real-time simulation method for debris flow according to claim 8, characterized in that: Step S330, the resultant forces of the model micro-units other than gravity are: The net forces acting on SPH particles, besides gravity, include gradient forces, viscous forces, and contact forces between SPH particles, as well as gradient forces, viscous forces, and contact forces from DEM particles. The combined forces acting on DEM particles, besides gravity, include gradient forces from the fluid, viscous forces, and contact forces between DEM particles.

10. The real-time simulation method for debris flow according to any one of claims 6 to 9, characterized in that: Implemented within the Niagara particle system of Unreal Engine, the global acceleration structure is a spatial grid-type global acceleration structure.

Images