Seepage erosion damage simulation method based on ML-CFD-DEM

Abstract

The invention discloses a seepage erosion damage simulation method based on an ML-CFD-DEM. The seepage erosion damage simulation method comprises the following steps: acquiring a test strength envelope of a soil body in a subsidence area; taking the deviation minimization of the test strength envelope and the predicted strength envelope as a target, constructing a peak strength prediction model based on machine learning, and performing inversion through a preset optimization algorithm to determine mesoscopic parameters; in a DEM solver, constructing a particle phase model according to the mesoscopic parameters, and in a CFD solver, establishing a steady-state seepage field corresponding to the particle phase model; a fluid pressure gradient field obtained by the CFD solver is mapped into equivalent seepage physical force acting on particles in the DEM solver through a one-way coupling interface, so that the particles are driven to move, and the soil collapse process is simulated. According to the method, intelligent and efficient calibration of mesoscopic parameters can be realized, and a high-credibility CFD-DEM one-way coupling model is established on the basis and is used for accurately deducing the whole process from seepage erosion to collapse of the soil body.

Description

Technical Field

[0001] This invention relates to the field of urban road collapse analysis, prediction and remediation, and in particular to a method for simulating seepage erosion damage based on ML-CFD-DEM. Background Technology

[0003] Urban road collapse disasters are characterized by their high degree of concealment, suddenness, and destructiveness. Due to its weak interparticle cohesion, loose structure, and high permeability, sandy soil is highly susceptible to internal erosion (piping) under the seepage action caused by pipeline leakage, which leads to the hollowing out of the soil skeleton and causes rapid and catastrophic road collapse.

[0004] Numerical simulation is a key tool for revealing the mechanisms of such disasters and predicting risks. The finite element method (FEM) is suitable for simulating macroscopic deformation of continuous media, but it struggles to accurately describe discontinuous processes such as particle migration and loss under seepage. The discrete element method (DEM) can accurately simulate the mesoscopic motion of granular systems, while computational fluid dynamics (CFD) can describe fluid flow within pores. By coupling CFD and DEM, the entire process of seepage erosion can be physically and realistically reproduced.

[0005] However, the accuracy of CFD-DEM simulations heavily relies on the proper setting of microscopic parameters in the DEM (such as contact stiffness and friction coefficient). These parameters cannot be directly obtained through physical tests. Traditionally, calibration is performed using a manual trial-and-error method, which involves repeatedly adjusting parameters, running DEM simulations, and comparing macroscopic responses until they match the results of geotechnical tests (such as the strength envelope obtained from direct shear tests). This process is computationally extremely costly; a single complete calibration often requires hundreds of DEM simulations, taking several days or even weeks, which greatly limits the application of this technology in practical engineering.

[0006] Therefore, developing a method that can efficiently and accurately calibrate the micro-parameters of DEM and achieve rapid and accurate simulation of seepage erosion and collapse has become a technical problem that urgently needs to be solved in this field. Summary of the Invention

[0007] To address the aforementioned issues, the present invention aims to provide a seepage erosion damage simulation method based on ML-CFD-DEM, which enables intelligent and efficient calibration of microscopic parameters. Based on this, a highly reliable CFD-DEM unidirectional coupling model is established to accurately simulate the entire process of soil seepage erosion to collapse.

[0008] Based on this, the present invention provides a method for simulating seepage erosion damage based on ML-CFD-DEM, the method comprising:

[0009] Obtain the test strength envelope of the soil in the collapse zone;

[0010] With the goal of minimizing the deviation between the experimental intensity envelope and the predicted intensity envelope, the micro-parameters are determined by inversion using a peak intensity prediction model built based on machine learning and a preset optimization algorithm.

[0011] In the DEM solver, a granular phase model is constructed based on the microscopic parameters, and in the CFD solver, a steady-state seepage field corresponding to the granular phase model is established.

[0012] Through a one-way coupling interface, the fluid pressure gradient field obtained by the CFD solver is mapped into an equivalent seepage force acting on the particles in the DEM solver, thereby driving the particle movement and simulating the soil collapse process.

[0013] The process of obtaining the microscopic parameters includes:

[0014] Obtain training samples, which include the microscopic parameters and peak shear strength in the discrete element direct shear simulation;

[0015] Using the microscopic parameters as input and the peak shear strength as output, the peak strength prediction model is trained.

[0016] With the goal of minimizing the deviation between the test strength envelope and the predicted strength envelope, a Bayesian optimization algorithm is used to iteratively search within a preset parameter space to determine the micro-parameter corresponding to the minimum deviation.

[0017] The minimum deviation is defined as the minimum weighted sum of squares of the differences between the predicted peak shear strength and the experimental peak shear strength under at least two different normal stress levels.

[0018] The peak intensity prediction model is the XGBoost model.

[0019] The microscopic parameters include vertical stress, effective modulus, stiffness ratio, friction coefficient, and rolling resistance coefficient.

[0020] The method further includes: the particle phase model constructed by the DEM solver contains the geometry for simulating the soil and pipe, and the steady-state seepage field is established in the CFD solver based on Darcy's law.

[0021] The method further includes setting different simulation conditions by changing the size and / or location of the damaged opening in the geometric structure, so as to analyze the evolution of seepage erosion path and subsidence pit morphology.

[0022] The calculation process of the equivalent seepage force includes:

[0023] ;

[0024] Among them, the This represents the equivalent percolation force acting on a single particle. For particle volume, The density of water, This represents the hydraulic gradient.

[0025] The soil collapse process is divided into three stages: slow development, settlement, and failure. The slope of the collapse pit in the failure stage is related to the friction coefficient and rolling resistance coefficient among the micro-parameters.

[0026] The method further includes: acquiring data on particle migration, porosity changes, and sinkhole development during the soil collapse process, and using this data to quantitatively assess and compare the road collapse risk under different pipeline damage conditions.

[0027] In this invention, machine learning models replace a large number of time-consuming DEM calculations, improving the efficiency of micro-parameter calibration by tens of times and making the engineering application of high-precision CFD-DEM simulation possible. Intelligent algorithms such as Bayesian optimization are used for global automatic optimization, avoiding the subjectivity and local optima problems of manual trial and error, ensuring the accuracy and reliability of parameter calibration. Based on reliable parameters, CFD-DEM coupled simulation can physically and realistically reproduce the complete dynamic process from particle initiation and migration to soil instability, clearly revealing its inherent evolutionary laws. This method can also be directly used for virtual experiments of geological disaster risk assessment and mitigation schemes, providing effective decision support for urban underground space safety. Attached Figure Description

[0028] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0029] Figure 1 This is a flowchart of a seepage erosion damage simulation method based on ML-CFD-DEM provided in an embodiment of the present invention; Figure 2 This is a flowchart of the process for obtaining microscopic parameters provided in the embodiments of the present invention; Figure 3 This is a schematic diagram of the stress-strain curve of sand under direct shear provided in an embodiment of the present invention; Figure 4 This is a schematic diagram of the simulated direct shear strength envelope of sand provided in an embodiment of the present invention; Figure 5This is a schematic diagram of the particulate phase model provided in an embodiment of the present invention; Figure 6 This is a schematic diagram of flow field calculation provided in an embodiment of the present invention; Figure 7 This is a schematic diagram of a collapse pit with a defect size of 13cm provided in an embodiment of the present invention, which is run to 500,000, 1,500,000 and 2,500,000 steps respectively; Figure 8 This is a schematic diagram of a collapse pit with a defect size of 15cm provided in an embodiment of the present invention, which is run to 500,000, 1,500,000 and 2,500,000 steps respectively; Figure 9 This is a schematic diagram of a collapse pit with a defect size of 17cm provided in an embodiment of the present invention, which is run to 500,000, 1,500,000 and 2,500,000 steps respectively; Figure 10 This is a schematic diagram illustrating the leakage development process at different locations of the rupture openings provided in the embodiments of the present invention. Detailed Implementation

[0030] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0031] Figure 1 This is a flowchart of a seepage erosion damage simulation method based on ML-CFD-DEM provided in an embodiment of the present invention. The method includes:

[0032] S101. Obtain the test strength envelope of the soil in the collapse zone.

[0033] Direct shear tests were conducted on soil samples, such as representative sand samples, obtained from the road collapse site. Specifically, samples were prepared to the target relative density and sheared under normal stresses of 100 kPa, 200 kPa, and 300 kPa to obtain their peak shear strength. This data was used to plot the Mohr-Coulomb strength envelope, which served as the sole basis for subsequent calibration.

[0034] This embodiment of the direct shear test strictly follows the test specifications. The required data and sampling include: gradation parameters, soil particle specific gravity, maximum / minimum void ratio, and natural moisture content. Based on these parameters, a specimen with the target relative density is prepared. At the same time, key operating conditions such as shear box dimensions, specimen thickness, sample preparation method, normal stress, and shear rate are forcibly recorded. Other site data are not collected.

[0035] S102. With the goal of minimizing the deviation between the test intensity envelope and the predicted intensity envelope, the micro-parameters are determined by inversion using a preset optimization algorithm based on the peak intensity prediction model constructed by machine learning.

[0036] The process of obtaining the microscopic parameters will be described in detail in the second embodiment, and will not be repeated here.

[0037] S103. In the DEM solver, a granular phase model is constructed based on the microscopic parameters, and in the CFD solver, a steady-state seepage field corresponding to the granular phase model is established.

[0038] like Figure 5 As shown, this embodiment uses existing software such as commercial discrete element method software PFC5.0, EDEM, and LIGGGHTS to construct a two-dimensional particle phase model. The model size is 8000 mm × 4000 mm, and a circular hole with a diameter of 1000 mm is set at the center coordinates (4000 mm, 1000 mm) to simulate a buried pipeline. Particle generation adopts layered and servo control technology to ensure that the sample morphology is closer to the actual soil. To balance computational efficiency and accuracy, the soil particle size is enlarged by 40 times (while keeping key dimensionless parameters unchanged).

[0039] In the DEM solver, particle motion is governed by Newton's second law:

[0040] ;

[0041] ;

[0042] Where m is the particle mass, It is the translational velocity of the particle. It is the contact force between particles and between particles and the wall. It is the force exerted by the fluid on the particle, i.e., the equivalent seepage force acting on a single particle, where g is the acceleration due to gravity. It is the particle angular velocity. It is the particle moment of inertia tensor. It is the bending moment experienced by the particle.

[0043] Based on the established particle phase model, a fluid domain of the same size is constructed in the CFD solver. Fluid phase modeling mainly includes domain discretization, boundary condition setting, and solver selection.

[0044] CFD simulations are based on the Navier-Stokes equations, also known as the NS equations, which include the continuity equation, momentum conservation equation, and energy conservation equation. Due to the complexity of solving the NS equations, discretization methods are typically used for numerical solutions. Under the low Reynolds number seepage conditions described in this embodiment, the NS equations can be simplified to Darcy flow form, and their governing equations can be expressed as two-dimensional Laplace equations:

[0045] ;

[0046] in, and Let x and z be the permeability coefficients. and These are the pressure gradients along the x and z directions, respectively, where h is the water head, and x and z refer to the horizontal and vertical directions in the two-dimensional coordinate system, respectively.

[0047] To meet the scale matching and numerical stability requirements of CFD-DEM coupled calculations, fluid mesh generation must satisfy the following criteria:

[0048] Mesh feature length The following resolution requirements should be met:

[0049]

[0050] Among them, the The minimum width of the watershed is represented by the following. Indicates the length of the fluid unit.

[0051] To avoid a single particle being affected by multiple fluid units, a single fluid unit should be able to accommodate multiple particles:

[0052]

[0053] Wherein, r is the average particle radius.

[0054] The specific settings for fluid solution in this embodiment are as follows:

[0055] (A) Fluid domain and mesh;

[0056] A fluid computational domain is established within a geometric region consistent with the particles, and a triangular computational mesh is generated using Gmsh (e.g., Figure 6 (As shown).

[0057] (B) Governing equations and solutions

[0058] Considering low Reynolds number seepage conditions, the pressure field is solved based on Darcy's law; the governing equations are discretized using the finite volume method and numerically solved using the FiPy solver.

[0059] (C) Boundary conditions

[0060] The upper boundary is set as the pressure inlet, the pipe opening as the pressure outlet, and the left, right and lower boundaries as fixed walls.

[0061] S104. Through a one-way coupling interface, the fluid pressure gradient field obtained by the CFD solver is mapped into an equivalent seepage force acting on the particles in the DEM solver, thereby driving the particle movement and simulating the soil collapse process.

[0062] This embodiment uses a unidirectional fluid-structure interaction method to simulate soil deformation under seepage through the following steps:

[0063] Step 1: Modeling the granular phase (DEM);

[0064] Establish the geometric model of soil and pipe in the particle phase solver; generate the particle system based on the calibrated microscopic parameters, and complete static equilibrium and consolidation; set the contact model and parameters between particles and between particles and the wall.

[0065] Step 2: Fluid phase (CFD) solution;

[0066] Set up the physical model and boundary conditions in the fluid solver; solve the Darcy-type steady-state seepage field to obtain the global pressure gradient distribution; output the pressure gradient field as an interpolable data field for coupling purposes.

[0067] Step 3: Field quantity mapping and update cycle;

[0068] Through the coupling interface, the fluid pressure gradient is interpolated to the centroid position of each particle every N DEM time steps (N=100 in this embodiment);

[0069] When the inlet / outlet boundary conditions change or the operating conditions switch, return to step 2 to recalculate the flow field and update the mapping.

[0070] Step 4: Unidirectional Force Application and Motion Regeneration

[0071] Apply an equivalent percolation force to each particle in the DEM system:

[0072] ;

[0073] Wherein, the For particle volume, The density of water, Indicates the hydraulic gradient;

[0074] The position and velocity of the particles are updated based on the force integral.

[0075] Step 5: Loop and Termination

[0076] Repeat steps 3 and 4 until the set simulation duration is reached.

[0077] The core feature of this unidirectional fluid-structure interaction method is:

[0078] Establish a soil-pipe model in the granular phase and set the contact parameters;

[0079] Solve for the steady-state pressure gradient field in a fluid phase;

[0080] The flow field data is interpolated to the particle centroid at a fixed time step frequency;

[0081] According to the formula Apply physical force to particles and update their motion;

[0082] Interpolation and motion updates are performed cyclically without feeding back particle motion to the fluid phase, i.e., no reverse coupling is performed.

[0083] Through the aforementioned fluid-structure interaction numerical simulation, the system obtains the macroscopic morphological evolution at each time step during the road collapse induced by pipeline leakage, fully reproducing the entire process from soil instability to final failure. Based on this dynamic evolution result, the collapse formation mechanism under seepage-soil interaction can be revealed, and the impact range and development trend of road settlement can be accurately predicted, thus providing key theoretical basis and technical support for the early prevention and engineering treatment of collapse accidents.

[0084] It should be noted that the fluid-structure interaction method used in this embodiment is based on existing technology: a CFD solver calculates the Darcy-type steady-state seepage pressure field, a DEM solver simulates the particle system motion, and the pressure gradient is applied to the particles as a volume force term in a unidirectional coupling manner. This coupling architecture, solver settings, and mesh strategy are all conventional implementations in the field and do not constitute a limitation on the scope of protection of this invention; provided that equivalence is ensured, a bidirectional coupling process can also be used without affecting the substantive content of this invention.

[0085] In this embodiment, the size and location of the pipe rupture are selected as variables, and the groundwater level is assumed to be level with the ground surface and there is no water flow in the pipe to design the simulation conditions, as shown in Table 1 (where the upper right refers to the northeast direction).

[0086] Table 1

[0087] Simulation results show that the soil erosion and collapse process caused by pipeline leakage can be divided into three typical stages:

[0088] The first stage, or slow development stage, involves soil erosion above the pipeline rupture point, forming a localized loose zone, but not yet extending to the surface. During this stage, there is no significant surface subsidence, the overall bearing capacity of the road is not significantly affected, and the risk of collapse is low.

[0089] The second stage is the settlement stage: the loosened area continues to develop upwards and extends to the ground surface, causing visible ground subsidence. The bearing capacity of the foundation gradually decreases, the structural stability weakens, and the risk of collapse increases significantly.

[0090] The third stage, or the failure stage, is characterized by the formation of a distinct inverted triangular sinkhole on the ground surface, complete instability of the foundation, and the occurrence of a collapse accident. The slope of the final sinkhole is close to the internal friction angle of the soil particles.

[0091] Figure 7 This is a schematic diagram of a collapse pit with a defect size of 13cm, as provided in an embodiment of the present invention, after running for 500,000, 1,500,000, and 2,500,000 steps respectively:

[0092] When the operation reached 500,000 steps, only a small amount of particles were lost near the damaged area, and no obvious subsidence was observed on the ground surface.

[0093] When the operation reaches 1.5 million steps, the loosened area extends to the surface, forming a parabolic subsidence zone;

[0094] When the operation reached 2.5 million steps, the sinkhole expanded significantly laterally, and its final shape approached that of an inverted triangle.

[0095] refer to Figure 8 and Figure 9 When the size of the ruptured opening increased to 15 cm and 17 cm, the development process accelerated significantly.

[0096] When the operation reached 500,000 steps, the loosened area had expanded to the ground surface, and significant subsidence and a decrease in the bearing capacity of the foundation had occurred.

[0097] When the circuit reached 2.5 million steps, the collapse was fully formed, the inverted triangular collapse pit was fully developed, and the foundation was completely unstable.

[0098] Depend on Figure 10 It is evident that the leakage development process at different rupture locations still follows the aforementioned three-stage pattern. When the rupture is directly above the pipeline, the collapse develops most rapidly; when the rupture is located to the upper right or directly to the right, the development path of the subsidence pit shifts to the right, indicating that soil failure always develops along the shortest path towards the rupture.

[0099] Figure 2 This is a flowchart of the process for obtaining mesoscopic parameters provided in an embodiment of the present invention. The process for obtaining mesoscopic parameters includes:

[0100] S201. Obtain training samples, which include the microscopic parameters and peak shear strength in the discrete element direct shear simulation.

[0101] The construction process of the peak intensity prediction model includes:

[0102] Conduct direct shear numerical simulations of discrete element soil (such as sand) with preset groups (e.g., 60 groups) to obtain the corresponding peak shear strength data.

[0103] The data for this invention comes from discrete element direct shear numerical simulations. Three levels of vertical stress—100, 200, and 300 kPa—were set in the simulation. Key micro-parameters were sampled within their reasonable engineering ranges, specifically: effective modulus 10–100 MPa, stiffness ratio 1.5–2.5, particle friction coefficient 0.1–0.8, and rolling resistance coefficient 0.1–0.8.

[0104] S202. Using the microscopic parameters as input and the peak shear strength as output, train the peak strength prediction model.

[0105] Using vertical stress, effective modulus, stiffness ratio, friction coefficient, and rolling resistance coefficient from discrete element direct shear simulation as input variables and peak shear strength as output variable, a dataset is constructed for training the peak strength prediction model.

[0106] The dataset is divided into training and test sets according to a preset ratio (e.g., 8:2). The multiple correlation coefficient (R²), root mean square error (RMSE), and mean absolute error (MAE) are used as evaluation metrics to comprehensively assess the prediction accuracy and stability of the model.

[0107] The XGBoost regression model (tree booster) was used for training. The main hyperparameters were set as follows:

[0108] n_estimators = 1000;

[0109] learning_rate=0.1;

[0110] max_depth=2;

[0111] min_child_weight=1;

[0112] gamma = 0.005;

[0113] subsample=0.8;

[0114] colsample_bytree=0.9.

[0115] The remaining parameters are the default values.

[0116] The peak intensity prediction model was evaluated using the multiple correlation coefficient R², root mean square error (RMSE), and mean absolute error (MAE).

[0117] See Figure 3 and Figure 4 In this embodiment, the evaluation results of the peak strength prediction model on the test set are: MAE=0.18, RMSE=0.05, R²=0.95, indicating that the peak strength prediction model has high prediction accuracy and strong generalization ability, and can stably predict peak shear strength.

[0118] S203. With the goal of minimizing the deviation between the test strength envelope and the predicted strength envelope, a Bayesian optimization algorithm is used to iteratively search within a preset parameter space to determine the micro-parameter corresponding to the minimum deviation.

[0119] Based on the trained prediction model, the peak shear strength of soil under different vertical stresses (100kPa, 200kPa, 300kPa) is predicted in direct shear simulation. The strength envelope is fitted and compared with the actual test results to determine a set of optimal micro-parameter values ​​for subsequent discrete element simulation analysis.

[0120] Based on the trained prediction model, the peak shear strength was calculated for vertical stresses of 100, 200, and 300 kPa. The Mohr-Coulomb strength envelope of the sand was then fitted using the three sets of stress-strength data points.

[0121] To determine a set of discrete element micro-parameters that can reproduce the real experimental results, this study uses a Bayesian optimization method for automatic parameter calibration. The specific steps are as follows:

[0122] Define the search space:

[0123] The engineering reasonable ranges of four key microscopic parameters (effective modulus, stiffness ratio, particle friction coefficient, and rolling resistance coefficient) are set as the optimization search space.

[0124] Construct the objective function:

[0125] The optimization objective is to achieve "the best fit between the simulated strength envelope and the experimental strength envelope". Specifically, the sum of squares of the differences between the simulated and experimental peak strengths under multiple normal stresses is calculated (in this embodiment, the weight of each stress level is 1), and this is used as the objective function value to find its minimum value.

[0126] Execute optimization process:

[0127] Latin hypercube sampling is used to generate a preset number (e.g., 200) of initial parameter samples, and a surrogate model is constructed based on a Gaussian process. The expected improvement (EI) is used as the acquisition function to guide the iterative search, which is carried out for a total of 80 rounds to approximate the optimal solution.

[0128] Determine the final parameters:

[0129] Select the parameter combination that performs best from the optimization results (e.g., 4 groups) and perform discrete element verification calculations. Finally, select the group that matches the experimental results the best.

[0130] Through the above process, a set of optimal discrete element micro-parameters for sand were obtained as follows: particle density 1800 kg / m³, porosity 0.16, effective modulus 50 MPa, stiffness ratio 2.0, friction coefficient 0.3, rolling resistance coefficient 0.3, local damping 0.7, and both normal and tangential viscous damping 0.2.

[0131] Verification Results: The stress-strain curves and strength envelopes recalculated using this set of parameters are highly consistent with the indoor test results. Based on this envelope, the equivalent shear strength indices for the sand are: cohesion 2 kPa and internal friction angle 30.65°. This result is very close to the measured values ​​(cohesion 3.67 kPa, internal friction angle 28.9°), verifying the reliability of the parameter calibration process. Therefore, this validated prediction model and mesoscopic parameters can be directly applied to subsequent CFD-DEM coupled calculations for infiltration collapse scenarios.

[0132] Introducing machine learning methods into the micro-parameter calibration process can significantly improve efficiency compared to traditional manual trial and error. For quantitative comparison, the average time to complete one DEM direct shear simulation (single parameter group, single stress level) is set as T.

[0133] Traditional manual trial and error method: It usually requires evaluating no less than 100 sets of micro parameters and verifying them at three vertical stress levels of 100, 200 and 300 kPa. A total of about 300 DEM simulations need to be completed, and the total time is about 300T.

[0134] Machine learning-assisted method (in this embodiment): This method combines Bayesian optimization (BO) with a trained prediction model to optimize parameters. The computational cost of the optimization process is far less than that of DEM simulation and can be ignored. This method only needs to perform DEM verification on a few optimal candidate parameters (4 groups in this example), with a total time consumption of approximately 4T.

[0135] Efficiency Comparison Analysis:

[0136] Ideal speedup: approximately 75x speedup (300T / 4T) in subsequent applications without considering training costs.

[0137] Total time for initial application: If the 60 DEM simulations required to build the initial training set are included, the total time for initial application is approximately 64T. Even so, compared to 300T for manual trial and error, the efficiency is still improved by approximately 4.7 times.

[0138] Long-term benefits: No need to pay for training again in subsequent repetitive calibration tasks, and can achieve an efficiency improvement of about 75 times or more.

[0139] The minimum deviation is defined as the minimum weighted sum of squares of the differences between the predicted peak shear strength and the experimental peak shear strength under at least two different normal stress levels.

[0140] In this invention, machine learning models replace a large number of time-consuming DEM calculations, improving the efficiency of micro-parameter calibration by tens of times and making the engineering application of high-precision CFD-DEM simulation possible. Intelligent algorithms such as Bayesian optimization are used for global automatic optimization, avoiding the subjectivity and local optima problems of manual trial and error, ensuring the accuracy and reliability of parameter calibration. Based on reliable parameters, CFD-DEM coupled simulation can physically and realistically reproduce the complete dynamic process from particle initiation and migration to soil instability, clearly revealing its inherent evolutionary laws. This method can also be directly used for virtual experiments of geological disaster risk assessment and mitigation schemes, providing effective decision support for urban underground space safety.

[0141] The above description is only a preferred embodiment of the present invention. It should be noted that for those skilled in the art, several improvements and substitutions can be made without departing from the technical principles of the present invention, and these improvements and substitutions should also be considered within the scope of protection of the present invention.

Claims

1. A method for simulating seepage erosion damage based on ML-CFD-DEM, characterized in that, include: Obtain the test strength envelope of the soil in the collapse zone; With the goal of minimizing the deviation between the experimental intensity envelope and the predicted intensity envelope, the micro-parameters are determined by inversion using a peak intensity prediction model built based on machine learning and a preset optimization algorithm. In the DEM solver, a granular phase model is constructed based on the microscopic parameters, and in the CFD solver, a steady-state seepage field corresponding to the granular phase model is established. Through a one-way coupling interface, the fluid pressure gradient field obtained by the CFD solver is mapped into an equivalent seepage force acting on the particles in the DEM solver, thereby driving the particle movement and simulating the soil collapse process.

2. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The process of obtaining the microscopic parameters includes: Obtain training samples, which include the microscopic parameters and peak shear strength in the discrete element direct shear simulation; Using the microscopic parameters as input and the peak shear strength as output, the peak strength prediction model is trained. With the goal of minimizing the deviation between the test strength envelope and the predicted strength envelope, a Bayesian optimization algorithm is used to iteratively search within a preset parameter space to determine the micro-parameter corresponding to the minimum deviation.

3. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 2, characterized in that, The minimum deviation is defined as the minimum weighted sum of squares of the differences between the predicted peak shear strength and the experimental peak shear strength at at least two different normal stress levels.

4. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The peak intensity prediction model is the XGBoost model.

5. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The microscopic parameters include vertical stress, effective modulus, stiffness ratio, friction coefficient, and rolling resistance coefficient.

6. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The method further includes: the particle phase model constructed by the DEM solver contains the geometry for simulating the soil and pipe, and the steady-state seepage field is established in the CFD solver based on Darcy's law.

7. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 6, characterized in that, The method further includes setting different simulation conditions by changing the size and / or location of the breach in the geometry, in order to analyze the evolution of seepage erosion paths and subsidence pit morphology.

8. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The calculation process for the equivalent seepage force includes: ; Among them, the This represents the equivalent percolation force acting on a single particle. For particle volume, The density of water, This represents the hydraulic gradient.

9. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The soil collapse process is divided into three stages: slow development, settlement, and failure. In the failure stage, the slope of the collapse pit is related to the friction coefficient and rolling resistance coefficient among the micro-parameters.

10. The seepage erosion damage simulation method based on ML-CFD-DEM as described in claim 1, characterized in that, The method further includes: acquiring data on particle migration, porosity changes, and sinkhole development during the soil collapse process, and using this data to quantitatively assess and compare the road collapse risk under different pipeline damage conditions.

Images