Simulation method, device and equipment for soil-rock mixed slope sliding

Abstract

The application provides a simulation method, device and equipment for soil-rock mixed slope sliding. In an initialization stage, a soil-rock mixed body model of the slope is constructed. In an execution stage, the following process is cyclically executed until a total time step is reached: soil particle field data and rock particle field data corresponding to each time step are calculated based on a soil particle field model and a rock particle field model constructed for the soil-rock mixed body model of the slope. The soil particle field model and the rock particle field model each include a density field term, a velocity field term, a displacement field term and a stress field term. The soil particle velocity field term includes a soil particle contact acceleration term corresponding to a soil particle contact resultant force. Thus, the method does not require a grid, avoids grid distortion problems caused by large deformation and sliding, can more accurately simulate the non-continuous contact and sliding process of the soil-rock mixed body, and can realize the contact movement and large deformation of the soil-rock interface without coupling other numerical methods, significantly improves the calculation efficiency, and makes parameter calibration more convenient.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of data processing, in particular to a simulation method, device and equipment for soil-rock mixed slope sliding. BACKGROUND

[0002] Soil-rock mixed slope is a special geological body formed by the mixed accumulation of soil and block stones, which is widely distributed in mountainous areas, valleys and artificial filling engineering. Its heterogeneity is significant, and its mechanical behavior is controlled by the coupling of multiple factors such as particle gradation, block stone spatial distribution, and soil-rock interface characteristics, resulting in complex instability mechanism and strong disaster burst. The traditional slope stability analysis method (such as finite element method) has significant limitations in the simulation of soil-rock mixed slope: on the one hand, the finite element method based on the assumption of continuous medium is difficult to depict the discontinuous contact and large deformation sliding process between block stones and soil, and the grid distortion problem leads to calculation failure; on the other hand, although the discrete element method can represent particle motion, it has low calculation efficiency, difficulty in parameter calibration, and insufficient modeling accuracy of complex boundary conditions. SUMMARY

[0003] The present application provides a simulation method, device and equipment for soil-rock mixed slope sliding, which can solve one of the problems in the background art.

[0004] To achieve the above-mentioned purpose, the present application adopts the following technical scheme:

[0005] In a first aspect, a simulation method for soil-rock mixed slope sliding is provided, comprising:

[0006] In the initialization stage, a soil-rock mixed body model of the slope composed of soil particles and rock particles is constructed, and initial particle parameters are assigned; and,

[0007] In the execution stage, before reaching the set total time step, the following process is executed in a loop: based on the soil particle field model and the rock particle field model constructed for the soil-rock mixed body model of the slope, the soil particle field data and the rock particle field data corresponding to each time step are calculated, wherein the soil particle field model includes: soil particle density field term, soil particle velocity field term, soil particle displacement field term and soil particle stress field term constructed based on smoothed particle hydrodynamics (SPH) algorithm, and the rock particle field model includes: rock particle density field term, rock particle velocity field term, rock particle displacement field term and rock particle stress field term constructed based on SPH algorithm, and the soil particle velocity field term includes: soil particle contact acceleration term corresponding to soil particle contact resultant force.

[0008] Based on the technical scheme, in the initialization stage, the slope soil-rock mixture model is constructed, and in the execution stage, the following process is executed in a loop until the total time step is reached: based on the soil particle field model and the rock particle field model constructed for the slope soil-rock mixture model, the soil particle field data and the rock particle field data corresponding to each time step are calculated, the soil particle field model and the rock particle field model both include: a density field term, a velocity field term, a displacement field term, and a stress field term, and the soil particle velocity field term includes a soil particle contact acceleration term corresponding to the soil particle contact force, so that the method does not require a grid, avoids the grid distortion problem caused by large deformation and sliding, can more accurately simulate the non-continuous contact and sliding process of the soil-rock mixture, and without coupling other numerical methods, the contact motion and large deformation of the soil-rock interface can be realized, the calculation efficiency is significantly improved, and the parameter calibration is more convenient.

[0009] In a possible design manner of the first aspect, the soil particle field model is:

[0010]

[0011] wherein i is the serial number of a base particle, j is the serial number of a neighboring particle adjacent to the base particle, N is the number of particles, t is the calculation time, is the particle density, m is the particle mass, v is the particle velocity, W is a smoothing kernel function, x, y, and z are the Cartesian components of the particle displacement, and alpha and beta represent the Cartesian components x, y, and z, is the total stress on the soil particle, is a Dirac function, is an artificial viscosity, is an artificial stress, is an additional force on the base particle, is the particle contact force, P is the isotropic pressure, and S is the anisotropic shear stress.

[0012] In a possible design manner of the first aspect, the particle contact force is:

[0013]

[0014] wherein, is a normal unit vector from the particle j to the particle i, is a normal force, is a tangential friction force, is a damping force, is a unit vector in the direction of the tangential velocity.

[0015] In a possible design manner of the first aspect, the rock particle field model is:

[0016]

[0017] where i is the serial number of the base particle, j is the serial number of the adjacent particle adjacent to the base particle, N is the number of particles, t is the calculation time, is the particle density, m is the particle mass, v is the particle velocity, W is the smoothing kernel function, x, y, z are the Cartesian components of the particle displacement, and a and b represent the Cartesian components x, y and z, is the total stress on the soil particles, is the Dirac function, is the artificial viscosity, is the artificial stress, is the additional force on the base particle, P is the isotropic pressure, and S is the anisotropic shear stress.

[0018] In a possible design manner of the first aspect, the execution stage further includes:

[0019] Based on the double search architecture, the same property particles are searched, and the interaction relationship between the base particles and the adjacent particles of the same property is established based on SPH:

[0020]

[0021]

[0022] where X is a general term of variables, is the field function of the base particle, is the derivative of the field function of the base particle, is the field function of the adjacent particle, is the smoothing kernel function, is the gradient of the smoothing kernel function, is the smoothing kernel radius of the influence range of the base particle.

[0023] In a possible design manner of the first aspect, the initialization stage includes:

[0024] Discretize the physical space to be simulated into uniform initial particles to form an initial discrete space;

[0025] According to the particle size grading curve of the initial particles, a plurality of particle regions of different particle sizes are divided, particles of each particle size range are generated in the initial discrete space, and each particle region is given an ID identifier for identifying particles of the same property; and

[0026] The initial particle parameters are assigned to the particles.

[0027] In a possible design of the first aspect, the initial particle parameters include: initial positions of the particles, initial velocities of the particles, initial stresses of the particles, particle size gradation, block stone space distribution, material parameters, and boundary conditions.

[0028] In a possible design of the first aspect, the soil particle field model and / or the rock particle field model further include: hydrostatic pressure adjustment conditions and / or shear stress adjustment conditions based on a return mapping algorithm.

[0029] In a second aspect, a device for simulating sliding of a soil-rock mixed slope is provided, and the device includes:

[0030] An initialization unit is configured to construct a soil-rock mixed slope model composed of soil particles and rock particles, and to assign initial particle parameters; and

[0031] An execution unit is configured to cyclically execute the following process before a total time step is reached: based on a soil particle field model and a rock particle field model constructed for the soil-rock mixed slope model, soil particle field data and rock particle field data corresponding to each time step are calculated, wherein the soil particle field model includes a soil particle density field term, a soil particle velocity field term, a soil particle displacement field term, and a soil particle stress field term constructed based on a smoothed particle hydrodynamics (SPH) algorithm, and the rock particle field model includes a rock particle density field term, a rock particle velocity field term, a rock particle displacement field term, and a rock particle stress field term constructed based on the SPH algorithm, and the soil particle velocity field term includes a soil particle contact acceleration term corresponding to a soil particle contact resultant force.

[0032] In a third aspect, an electronic device is provided, and the electronic device includes: a processor, and a memory coupled to the processor, the memory being configured to store a computer program, and the processor being configured to execute the computer program stored in the memory, so that the electronic device performs the method for simulating sliding of a soil-rock mixed slope according to any one of possible implementation manners of the first aspect. BRIEF DESCRIPTION OF DRAWINGS

[0033] In order to more clearly illustrate the technical solutions in the embodiments of the present application, the following will briefly introduce the drawings needed to be used in the embodiments or related technical descriptions. Obviously, the drawings in the following description are only some embodiments of the present application, and other drawings can be obtained by those skilled in the art without any creative labor.

[0034] Figure 1 is a whole flowchart of the method for simulating sliding of a soil-rock mixed slope provided by the embodiments of the present application;

[0035] Figure 2 is a schematic diagram of the particle domain search method of step 1 provided by the embodiments of the present application;

[0036] Figure 3 is a flow chart of the contact pair detection of step 6 provided by the embodiments of the present application;

[0037] Figure 4 is a schematic diagram of the slope model and boundary condition setting provided by the embodiments of the present application;

[0038] Figure 5 is a schematic diagram of the simulation results of the plastic strain distribution field, potential sliding surface, displacement field, etc. of the soil-rock mixed slope provided by the embodiments of the present application, considering different contact stiffnesses of rock mass and soil.

[0039] Figure 6 is a schematic diagram of the simulation results of the plastic strain distribution field, potential sliding surface, displacement field, etc. of the soil-rock mixed slope provided by the embodiments of the present application, considering different internal friction angles of rock mass and soil. DETAILED DESCRIPTION

[0040] In order to make the objects, technical solutions and advantages of the present application clearer, the present application will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present application, and are not used to limit the present application.

[0041] It should be noted that although the functional modules are divided in the device schematic diagram, and the logical order is shown in the flow chart, in some cases, the steps shown or described can be executed in a manner different from the module division in the device or the order in the flow chart. The terms "first", "second", etc. in the specification and claims and the above-described drawings are used to distinguish similar objects, and do not necessarily describe a specific order or sequence.

[0042] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which the present application belongs. The terms used herein are only for the purpose of describing the embodiments of the present application, and are not intended to limit the present application.

[0043] As Figure 1As shown, a simulation method for soil-rock mixed slope sliding failure is provided to overcome the shortcomings of the prior art. The method includes mainly involving: first, obtaining the particle size distribution curve of the research slope and the soil-rock mixed information and the strength parameters. Then, a two-dimensional soil-rock mixed slope model is established by using the improved soil-rock mixture-smoothed particle hydrodynamics method (Soil and Rock Mixture-SPH, SRM-SPH). Next, the elastoplastic constitutive model based on the non-associated flow rule of the Drucker-Prager (DP) yield criterion and the improved kernel function linear elastic constitutive model are used to describe the deformation inside the soil and rock, and the point-point contact criterion is coupled to realize the contact movement between the soil-rock interface. The plastic strain distribution field, potential sliding surface, and displacement field of the slope are simulated and calculated by using the soil-rock mixed model.

[0044] As shown in Figure 1 The simulation method for soil-rock mixed slope sliding of the present embodiment includes the following steps:

[0045] Step 1: Assign initial particle parameters to the slope model: including the position of the particle, the initial velocity of the particle, the initial stress of the particle, the particle size distribution, the spatial distribution of the block stone, the material strength parameters and the boundary conditions;

[0046] Step 2: Rock particles and rock particles interact in pairs under the double search architecture, and soil particles and soil particles interact in pairs;

[0047] Step 3: Calculate the change of soil density, construct the stress-strain relationship of the Drucker-Prager ideal elastoplastic soil material model under the non-associated flow rule, and the changes of the remaining internal forces and external forces acting on the soil;

[0048] Step 4: Regression analysis is performed on the soil that satisfies the yield criterion, and the stress state is pulled back to the yield surface from the state of exceeding the plastic yield surface;

[0049] Step 5: Calculate the change of rock density, construct the rock material model and the changes of the remaining internal forces and external forces acting on the rock;

[0050] Step 6: Calculate the contact force between the rock and the soil and the contact force between the rocks;

[0051] Step 7: Calculate the stress field and displacement field during the failure process of the soil-rock mixed slope;

[0052] Step 8: If the loop is not completed, return to Step 2 and continue the loop until the predetermined total time step is reached. Output the final velocity field, block stone distribution field, and displacement field of the particles.

[0053] Specifically, step 1: the soil and rock mixture with particle size distribution and block stone space distribution is discretized into particles, each particle is assigned with position information, initial velocity, initial stress, strength parameters, and boundary conditions, etc. which will be specifically embodied in the following cases. The block stone distribution is detailed as follows:

[0054] 1.1 Initial spatial discretization: discretize the physical space to be simulated into uniform initial particles to form an initial discrete space, as shown in the upper left of Fig. Figure 2 .

[0055] 1.2 Import the particle size distribution curve: obtain and import the particle size distribution curve of the target particle system, which reflects the corresponding relationship between the radius and volume percentage of the debris particles, as shown in the upper right of Fig. Figure 2 .

[0056] 1.3 Generate initial range of particles and assign different IDs: according to the particle size distribution curve, divide multiple particle regions of different particle sizes, generate particles of each particle size range in the initial discrete space, and assign a unique ID to each particle region, as shown in the lower right of Fig. Figure 2 , different particle regions with IDs are generated.

[0057] 1.4 Search for particles within different ID numbers: search and identify the particle set corresponding to each ID mark through the adjacent particle algorithm, complete the particle marking of different particle size distribution regions, as shown in the lower left of Fig. Figure 2 , the accurate identification of particles in each ID region is realized.

[0058] Step 2: based on the double search architecture, i.e. same ID search each other, different ID does not search, the basic particle and its influence domain within the adjacent same ID particle pair, through the link list search method, wherein the basic particle is acted on the adjacent same ID particle by the following formula:

[0059] (1)

[0060] (2)

[0061] In the formula: is the serial number of the basic particle, j is the serial number of the adjacent particle, N is the total number of particles; Xᵢ is the spatial position vector of the i-th particle (target particle, i.e. the particle whose field variable f or its gradient ∇⋅f needs to be calculated); Xⱼ is the spatial position vector of the j-th particle (adjacent particle, i.e. the surrounding particle that contributes to the field variable of the target particle i); is the field function of the basic particle i; is the field function of the adjacent particle j; is the field function derivative of the underlying particle i; is the particle density, m is the particle mass; is the smoothing kernel radius of the influence range of the underlying particle i; is the smoothing kernel function. is the gradient of the smoothing kernel function. A B-spline function is adopted as the kernel function to reduce the effect of instability during particle compression, which is in the form of:

[0062] (3)

[0063] where , denotes the normalized distance, is the distance between particle i and j, and λ is a normalization constant (dimensionless), whose core role is to ensure that the kernel function satisfies the "normalization condition" (i.e., the sum of the weights of all neighboring particles on the target particle is 1), thereby ensuring the conservation of physical quantities such as density and mass. The value of λ is determined by the spatial dimension: in two-dimensional space ,

[0064] Step 3: Calculate the density change, calculate the artificial viscosity , artificial stress and the total stress and external force change, which can be expressed as:

[0065] (4)

[0066] where the Greek superscripts α and β represent the Cartesian components x, y, and z, and the Einstein convention applies to repeated indices; m and ρ represent the mass and density of the particle, respectively, t is the calculation time, is the particle velocity; δ is the Dirac function, which represents the additional force on the underlying particle, such as gravity. and are the artificial viscosity and artificial stress, respectively, used to reduce non-physical oscillations caused by zero-energy modes and stretching instability; is the artificial viscosity, used to reduce non-physical oscillations caused by zero-energy modes. The artificial viscosity is defined by the following equation:

[0067] (5)

[0068] (6)

[0069] (7)

[0070] (8)

[0071] where c is the speed of shear wave, describing the propagation speed of shear wave in the material, determined by the material elastic modulus (reflecting the material stiffness, the ability to resist deformation) and the reference density (the initial density of the material), is the average density of particles; a and b are artificial viscosity coefficients, determined by numerical experiments, used to adjust the "strength" of artificial viscosity, balancing numerical stability and physical reality (too large coefficients will over-dissipate, too small will not suppress vibration); u is the intermediate calculation of artificial viscosity, reflecting the "severity" of relative motion between particles, which is the core variable to adjust the local strength of artificial viscosity, calculated by equation (6) and is and coordinate values, is the elastic modulus, is the initial density, is the average density.

[0072] This oscillation can be effectively suppressed by adding artificial viscosity in the momentum equation, which will eliminate the unwanted motion in the undisturbed area, thus helping to produce a more realistic stress distribution and results in this area and ultimately the sliding distance, is the artificial stress, used to handle tensile instability in solid simulation; obtained by the following equation:

[0073] (9)

[0074] (10)

[0075] where: represents the particle spacing the value of the smoothing kernel function calculated at the particle spacing; is the initial particle spacing, represents the smoothing kernel function calculated at the initial particle spacing, which is a constant for normalization; and are related to the stress state of the rock-soil material, where is the tensor related to the stress state of particle i, is the tensor related to the stress state of particle j, used to construct the force resisting tensile instability. Its definition is based on the principal stress coordinate system of particle i, which can be calculated by the following equation:

[0076] (11)

[0077] (12)

[0078] where: is the principal stress direction angle of particle i, i.e. the angle between the maximum principal stress direction and the X-axis; , is the normal stress component of particle i in its principal stress coordinate system, is the reference response term sensitive to effective tensile stress (xx direction), is the reference response term sensitive to effective tensile stress (yy direction), is the actual response term in the xx direction in the global coordinate system, is the actual response term in the yy direction in the global coordinate system, is the tangential coupling response term in the xy direction in the global coordinate system; the optimal values of e = 0.3 and n = 4 are obtained by Gray et al. from the dispersion equation, where , and are calculated according to the following equations:

[0079] (13)

[0080] (14)

[0081] where, is the instantaneous normal stress component in the x-direction in the global coordinate system, is the instantaneous normal stress component in the y-direction in the global coordinate system, is the instantaneous shear stress component in the global coordinate system.

[0082] In addition to adding an artificial viscosity term, the gradient calculation of the traditional SPH relies on the derivative of the kernel function, but in the case of uneven particle distribution or high gradient area, numerical instability and energy non-conservation phenomena are prone to occur. For example, the initial kernel function used in equation (14) does not satisfy the zero-order or first-order completeness. In order to improve the consistency and accuracy of the SPH method, the kernel gradient correction (KGC) technique is used in SPH to improve the approximation accuracy. The KGC technique adjusts the gradient term of the kernel function by introducing a correction matrix or a weighting factor. which can be adjusted as follows:

[0083] (15)

[0084] (16)

[0085] In the formula: , is the relative position vector component of particle j with respect to particle i, i.e. , ; represents the SPH particle volume, and its value is . The correction matrix of the particle i, which physically compensates for the gradient calculation error caused by the uneven distribution of particles, its inverse matrix To ensure that the nuclear gradient approximation satisfies the first-order completeness, thereby significantly improving the calculation accuracy; the KGC-based gradient SPH particle approximation scheme has a second-order accuracy, which can ensure the symmetry and conservation of the SPH nuclear gradient calculation.

[0086] Total stress tensor Composed of isotropic pressure P and anisotropic shear stress Two parts:

[0087] (17)

[0088] Where, is the kronecker symbol.

[0089] The isotropic pressure P is not obtained from the state equation, but is directly obtained from the constitutive equation according to its standard definition, and the formula is as follows:

[0090] (18)

[0091] After various transformations, the stress-strain relationship of the Drucker-Prager ideal elastic-plastic rock and soil material model under the non-associated flow rule can be expressed as:

[0092] (19)

[0093] (20)

[0094] In the formula, denotes the rotation rate tensor, and denotes the rigid rotation of the material unit, which is defined as , where the first two terms are the action terms of Jaunman stress rate, which are used to ensure the objectivity of the stress tensor under material rotation, and are independent of the coordinate system; γ is a dummy index in tensor operation, G is the shear modulus, and K is the bulk modulus; is the partial strain rate tensor, which represents the rate of shape change; is the volume strain rate, which is the trace of the strain rate tensor (i.e. + In the two-dimensional case); The plastic multiplier represents a scalar that determines the magnitude of plastic deformation. When the material is in an elastic state, = 0, when plastic deformation occurs, > 0; the third and fourth terms are related to elastic action, and the last term is related to plastic action. Among them is the kronecker symbol. Among them , and It means that the superscripts a and β of are simultaneously taken as x, y, z. The key parameters are solved as follows:

[0095] (21)

[0096] (22)

[0097] (23)

[0098] (24)

[0099] In the formula: and J2are the first and second invariants of stress tensor, respectively. The former and are Drucker-Prager constants, which are related to the material constants of cohesion c and internal friction angle .

[0100] (25)

[0101] Step 4: Due to the characteristics of the ideal elastic-plastic constitutive model, when the material is in plastic deformation or in tension state during the calculation process, numerical calculation may produce errors, and appropriate stress adjustment is needed to meet the calculation requirements. The core idea of the method is: the first case, if the material particle is in tension state, the so-called tensile crack will occur, then the return mapping algorithm is used to adjust the hydrostatic pressure to remove the tensile cracking. When the stress state of the material particle reaches the tension state in the n-step calculation, it satisfies the following conditions:

[0102] (26)

[0103] wherein, is the total normal stress component in the x-direction of the global coordinate system at the adjusted time step n, is the total normal stress component in the y-direction of the global coordinate system at the adjusted time step n.

[0104] The second case, when plastic deformation occurs, calculation errors may cause the stress state of the material particle to be beyond the plastic yield surface. That is to say, the stress state of the particle must always remain on the plastic yield surface during plastic loading, so the return mapping algorithm is used to adjust the shear stress to keep the stress state on the yield surface. When the stress state of the material particle exceeds the yield surface in the n-step calculation, it satisfies the following conditions:

[0105] (27)

[0106] Here the scaling factor is introduced to adjust the magnitude of shear stress. For the Drucker-Prager yield criterion, the scaling factor at time step n is given by:

[0107] (28)

[0108] The shear stress component is adjusted by the scaling factor according to the relation:

[0109] (29)

[0110] where is the deviatoric normal stress tensor component in x-direction at time step n in global coordinate system, is the deviatoric normal stress tensor component in y-direction at time step n in global coordinate system, is the deviatoric shear stress tensor component in x-y direction at time step n in global coordinate system, is the deviatoric shear stress tensor component in x-z direction at time step n in global coordinate system, is the deviatoric shear stress tensor component in y-z direction at time step n in global coordinate system.

[0111] Step 5. Calculate the momentum and mass equations for soil and rock.

[0112] (30)

[0113] Step 6. Calculate the contact forces between rock and soil, and between blocks and blocks, mainly through particle contact forces, which are calculated by the Hertz contact model. First, a new link-list search is built, and then all possible pairs of particles that may come into contact are globally detected: for any two particles i and j, if their types are different and the distance is less than the contact detection distance they are considered as a contact pair. The distance between particles is defined as:

[0114] (31)

[0115] where represents the coordinates of particle i in the d-th dimension, represents the coordinates of particle j in the d-th dimension, and dim represents the total dimension.

[0116] For each contact pair, the normal force , the tangential friction force and the damping force are calculated. The normal force is given by the Hertz model:

[0117] (32)

[0118] where, is the normal stiffness coefficient, which is set to:

[0119] (33)

[0120] The magnitude of the tangential friction force is determined by the product of the normal force and the friction coefficient, which is chosen to be the static friction coefficient or the kinetic friction coefficient :

[0121] (34)

[0122] where the relative tangential velocity is obtained by subtracting the normal component from the relative velocity :

[0123] (35)

[0124] Here, is the unit normal vector pointing from particle j to particle i. The damping force is proportional to the magnitude of the tangential velocity:

[0125] (36)

[0126] where, is the damping coefficient. The total contact force is:

[0127] (37)

[0128] where, is the unit vector in the direction of the tangential velocity, i.e. when . Finally, the contact forces act in opposite directions on the two particles, i.e. particle i experiences a force and particle j experiences a force . Applying these contact forces to the soil mass yields the updated particle accelerations:

[0129] (38)

[0130] where, and are the masses of particles i and j, respectively.

[0131] Step 7: Calculate the stress, velocity, and displacement fields during the failure process of the soil-rock mixed slope. For rock:

[0132] (39)

[0133] For soil:

[0134] (40)

[0135] Step 8: The frog leapfrog algorithm is a commonly used time integration algorithm in SPH due to its high computational efficiency and strong global search capability. It integrates the kernel function in discrete form. In the frog leapfrog algorithm, the particle density, velocity, and position progress over time.

[0136] (41)

[0137] in, It is the initial time. This refers to the time step. The time step for particles within the computational domain is chosen to satisfy the Courant-Friedrichs-Lewy (CFL) condition.

[0138] Specific Case Analysis

[0139] Based on the plane strain assumption, a soil-rock mixed slope model and various boundary conditions are established, such as... Figure 4 As shown. To fully understand the impact of different soil and rock particles on the overall slope failure, the soil and rock mixture distributed on the upper part of the slope is represented in the SPH model. The total number of SPH solid particles is 27,000, and the number of virtual particles is approximately 1,100 (particles at the 3-layer displacement boundary are used to prevent boundary truncation). The initial particle spacing is set to 0.2 m. To prevent false viscosity caused by an excessively large smooth core radius, the smooth core radius is set to approximately 1 times the initial particle spacing (h = 0.2 m) for stress calculation. Point-to-point contact algorithms are used between different types of particles, such as soil and rock, and rock and rock. However, the SPH integral solution algorithm is used within particles, such as within soil or rock blocks. Considering the stability of the calculation, the time step is set according to CFL. The physical properties of rock blocks and strength parameters are shown in Table 1. SPH-related parameters, such as artificial viscosity and artificial stress, are also shown in Table 1. All strength physical properties can be obtained on-site. Specific SPH micro-parameters need to be adjusted according to actual working conditions.

[0140] Regarding the boundary settings, it can be seen that this simulation has two types of displacement boundaries (lateral constraints on the left and right sides, and full constraints at the bottom) and corresponding stress boundaries. Natural gravity is applied on top of this, and the boundary truncation problem is handled through three layers of virtual particles to ensure accurate transmission of stress and displacement fields. All boundary conditions can be solved by superimposing two different sets of smooth kernel functions to adapt to the numerical processing requirements of continuous and discontinuous regions in the soil-rock mixture.

[0141] For the simulation setting, firstly, an initial stress balance calculation of 1000 steps is performed, at this time, the strength parameters of the soil body such as the cohesion and the internal friction angle are set to be maximum values, so that the model reaches a static balance state under the action of gravity and without plastic deformation. After the stress field is stable, the real strength parameters of the soil body can be further applied to trigger the slope instability process, and the whole process from the initial balance state to the sliding failure of the soil-rock mixture is simulated.

[0142] Table 1 Mechanical and physical parameters of the soil-rock slope

[0143]

[0144] The initial parameters of the soil-rock mixture slope model include material parameters (strength parameters, elastic modulus, Poisson's ratio, cohesion, internal friction angle, particle size distribution, block stone spatial distribution, etc. of the soil body and the rock body) and boundary conditions (displacement boundary, stress boundary, etc.);

[0145] The key parameter calculation of the model particles in a single time step includes the pairing between particles based on the linklist search method (1)-(2), the updating of the smoothing kernel function (3), the change of the density rate (4), the change of the internal force (4), the calculation of the stress-strain relationship of the soil body and the rock body (19)-(20), the judgment of the plastic zone and the tensile zone of the soil body according to the Drucker-Prager criterion (25), the calculation of the stress state of the particles after the adjustment of the plastic zone and the tensile zone (26)-(29), the calculation of the contact force between the rock body and the soil body and the contact force between the blocks (31)-(38), and finally the characteristic parameter values of the soil body and the rock body in a single time step (39)-(40);

[0146] The time integration of the key parameters in a single time step is performed to obtain the characteristic parameter values of the particles in each time step (41);

[0147] It is judged whether the calculation step reaches the set calculation step, if not, steps (2)-(3) are repeated; if the maximum calculation step is reached, the calculation is terminated.

[0148] Through the above steps, the simulation results such as the plastic strain distribution field, the potential sliding surface and the displacement field of the soil-rock mixture slope are obtained, as shown in Figure 5 and Figure 6 .

[0149] The embodiment also provides a simulation device for the sliding of a soil-rock mixture slope, which comprises:

[0150] An initialization unit is configured to construct a slope soil-rock mixture model composed of soil particles and rock particles, and to assign initial particle parameters; and

[0151] The execution unit is configured to cyclically execute the following process before reaching the set total time step: based on the soil particle field model and the rock particle field model constructed for the soil-rock mixture slope model, the soil particle field data and the rock particle field data corresponding to each time step are calculated, wherein the soil particle field model includes a soil particle density field term, a soil particle velocity field term, a soil particle displacement field term and a soil particle stress field term constructed based on a smoothed particle hydrodynamics (SPH) algorithm, and the rock particle field model includes a rock particle density field term, a rock particle velocity field term, a rock particle displacement field term and a rock particle stress field term constructed based on the SPH algorithm, and the soil particle velocity field term includes a soil particle contact acceleration term corresponding to a soil particle contact force.

[0152] The embodiment relates to the technical field of geological engineering, in particular to a simulation method and device for failure and damage of a soil-rock mixed slope. The embodiment is based on a smoothed particle hydrodynamics (SPH) theoretical framework, and first proposes an SRM-SPH method for simulating a rock contact soil-rock mixed slope. The method can accurately describe the interaction between soil and rock and the large deformation sliding process through a double search architecture, a coupling point contact criterion and an elastic constitutive model and an elastoplastic constitutive model. Based on the double search architecture, rock and soil particles search internally, and false viscosity that may occur in a traditional SPH architecture is overcome. A contact model is proposed at a soil-rock interface to effectively simulate the contact force at the soil-rock interface, and simultaneously solve the contact between soil and rock and the interface contact problem. In addition, the method has simple parameter calibration and high calculation efficiency, and has higher soil-rock mixed modeling efficiency and calculation efficiency than MPM-DEM and MPM-DDA coupling, is suitable for stability analysis and disaster prediction of large-scale slope engineering, and has significant technical advantages and engineering application prospects.

[0153] The improved SRM-SPH method is applied to the sliding simulation of soil-rock mixed slopes for the first time, and the whole process analysis from large deformation of soil, plastic zone expansion to sliding failure is realized. Compared with the traditional finite element method (FEM), the SPH method does not need mesh, avoids the mesh distortion problem caused by large deformation and sliding, and can more accurately simulate the non-continuous contact and sliding process of soil-rock mixture. Compared with the SPH-DEM and SPH-DDA hybrid method, the embodiment is based on a unified SPH framework, and does not need to couple other numerical methods to realize the contact motion and large deformation of the soil-rock interface, the calculation efficiency is significantly improved, and the parameter calibration is more simple. In addition, the method of the embodiment effectively handles the false cohesion phenomenon that may occur between soil and rock through the point-point contact criterion and the double search framework, and provides more accurate plastic zone and displacement field distribution. The embodiment can be widely applied to the stability analysis, landslide disaster evolution prediction, disaster prevention scheme optimization of soil-rock mixed slope, tailings dam, roadbed filling and other engineering scenes, and is especially suitable for large-scale and long-term engineering simulation requirements, and has significant technical substitution value and engineering application prospect compared with traditional numerical methods.

[0154] The electronic device can be a desktop computer, a notebook computer, a palm computer, a cloud server, and the like. The electronic device can include, but is not limited to, a processor and a memory.

[0155] The electronic device can be a desktop computer, a notebook computer, a palm computer, a cloud server, and the like. The electronic device can include, but is not limited to, a processor and a memory.

[0156] The processor can be a central processing unit (CPU), and can also be other general-purpose processors, digital signal processors (DSP), application specific integrated circuits (ASIC), field-programmable gate arrays (FPGA) or other programmable logic devices, discrete gates or transistor logic devices, discrete hardware components, and the like. The general-purpose processor can be a microprocessor or the processor can also be any conventional processor. The processor is a control center of the electronic device, and connects all parts of the electronic device through various interfaces and lines.

[0157] The memory can be used to store the computer program, and the processor realizes various functions of the electronic device by running or executing the computer program stored in the memory and calling data stored in the memory.

[0158] The memory can mainly include a program storage area and a data storage area. The program storage area can store an operating system, application programs required by at least one function, and the like. The data storage area can store data created according to the use of the mobile phone, and the like. In addition, the memory can include a high-speed random access memory, and can also include a non-volatile memory, such as a hard disk, a memory, a plug-in hard disk, a smart media card (SMC), a secure digital (SD) card, a flash card, at least one disk storage device, a flash memory device, or other volatile solid-state memory devices.

[0159] The embodiment of the present application further provides a storage medium, which is a computer readable storage medium, and the computer program is stored in the computer readable storage medium. When the computer program is executed by a processor, the steps of each method embodiment described above can be realized. The computer program includes computer program code, which can be in the form of source code, object code, an executable file, or some intermediate form, etc. The computer readable medium can include any entity or device capable of carrying the computer program code, a recording medium, a U disk, a mobile hard disk, a magnetic disk, an optical disk, a computer memory, a read-only memory (ROM), a random access memory (RAM), an electrical carrier signal, a telecommunication signal, a software distribution medium, etc.

[0160] The embodiment of the present application further provides a computer program product, which includes a computer program or instructions, and when the computer program or instructions are run on a computer, the computer program or instructions make the computer execute the method of any possible implementation manner described above.

[0161] The above is the preferred embodiment of the present application. It should be noted that, for those skilled in the art, without departing from the principle of the present application, a number of improvements and refinements can be made, which are also considered to be within the protection scope of the present application.

Claims

1. A simulation method for sliding of soil-rock mixed slopes, characterized in that, include: In the initialization phase, a slope soil-rock mixture model composed of soil particles and rock particles is constructed, and initial particle parameters are assigned. as well as, During the execution phase, before reaching the set total time step, the following process is executed cyclically: Based on the soil particle field model and rock particle field model constructed for the slope soil-rock mixture model, the soil particle field data and rock particle field data corresponding to each time step are calculated. The soil particle field model includes soil particle density, soil particle velocity, soil particle displacement, and soil particle stress terms constructed based on the Smooth Particle Hydrodynamics (SPH) algorithm. The rock particle field model includes rock particle density, rock particle velocity, rock particle displacement, and rock particle stress terms constructed based on the SPH algorithm. The soil particle velocity term includes the soil particle contact acceleration term corresponding to the resultant contact force of soil particles. The soil particle field model is as follows: Where i is the index of the fundamental particle, j is the index of the nearest neighbor particle, N is the number of particles, t is the computation time, ρ is the particle density, m is the particle mass, v is the particle velocity, W is the smoothing kernel function, x, y, and z are the Cartesian components of the particle displacement, and α and β represent the Cartesian components x, y, and z. It is the total stress on soil particles. It is the Dirac function. It's artificial viscosity. It is artificial stress. It is an additional force acting on fundamental particles. It represents the net force of particle contact, where P is the isotropic pressure and S is the anisotropic shear stress. The rock mass particle field model is as follows: Where i is the index of the fundamental particle, j is the index of the nearest neighbor particle, N is the number of particles, t is the computation time, ρ is the particle density, m is the particle mass, v is the particle velocity, W is the smoothing kernel function, x, y, and z are the Cartesian components of the particle displacement, and α and β represent the Cartesian components x, y, and z. It is the total stress on soil particles. It is the Dirac function. It's artificial viscosity. It is artificial stress. It is the additional force on the fundamental particle, P is the isotropic pressure, and S is the anisotropic shear stress.

2. The simulation method for sliding of soil-rock mixed slope as described in claim 1, characterized in that, The net force of particle contact is: in, Let be the normal unit vector pointing from particle j to particle i. For normal force, For tangential friction, For damping force, It is a unit vector in the direction of tangential velocity.

3. The simulation method for sliding of soil-rock mixed slopes as described in claim 1 or 2, characterized in that, The execution phase also includes: Based on a dual-search architecture, particles with similar properties are searched, and based on SPH, the interaction relationships between fundamental particles and neighboring particles with similar properties are established: Where X is the collective term for variables. It is the field function of the fundamental particles. It is the derivative of the field function of the fundamental particle. It is the field function of the neighboring particle. It is a smooth kernel function. The gradient of the smooth kernel function, It is the smooth nuclear radius of the influence range of the fundamental particle.

4. The simulation method for sliding of soil-rock mixed slope as described in claim 3, characterized in that, The initialization phase includes: The physical space to be simulated is discretized into uniform initial particles to form an initial discrete space. Based on the initial particle size distribution curve, several particle regions with different particle size ranges are divided. Particles of each particle size range are generated in the initial discrete space, and each particle region is assigned an ID identifier to identify particles with the same properties; and, Assign initial particle parameters to the particles.

5. The simulation method for sliding of soil-rock mixed slope as described in claim 4, characterized in that, The initial particle parameters include: initial position of the particle, initial velocity of the particle, initial stress of the particle, particle size distribution, spatial distribution of the boulders, material parameters, and boundary conditions.

6. The simulation method for sliding of soil-rock mixed slope as described in claim 1, characterized in that, The soil particle field model and / or the rock particle field model further include: hydrostatic pressure adjustment conditions and / or shear stress adjustment conditions based on the return mapping algorithm.

7. A simulation device for sliding of soil-rock mixed slope, characterized in that, include: The initialization unit is used to construct a slope soil-rock mixture model composed of soil particles and rock particles, and to assign initial particle parameters. as well as, The execution unit is used to cyclically execute the following process before reaching the set total time step: Based on the soil particle field model and rock particle field model constructed for the slope soil-rock mixture model, calculate the soil particle field data and rock particle field data corresponding to each time step. The soil particle field model includes: soil particle density field terms, soil particle velocity field terms, soil particle displacement field terms, and soil particle stress field terms constructed based on the Smoothed Particle Hydrodynamics (SPH) algorithm. The rock particle field model includes: rock particle density field terms, rock particle velocity field terms, rock particle displacement field terms, and rock particle stress field terms constructed based on the SPH algorithm. The soil particle velocity field term includes: the soil particle contact acceleration term corresponding to the resultant contact force of soil particles. The soil particle field model is as follows: Where i is the index of the fundamental particle, j is the index of the nearest neighbor particle, N is the number of particles, t is the computation time, ρ is the particle density, m is the particle mass, v is the particle velocity, W is the smoothing kernel function, x, y, and z are the Cartesian components of the particle displacement, and α and β represent the Cartesian components x, y, and z. It is the total stress on soil particles. It is the Dirac function. It's artificial viscosity. It is artificial stress. It is an additional force acting on fundamental particles. It represents the net force of particle contact, where P is the isotropic pressure and S is the anisotropic shear stress. The rock mass particle field model is as follows: Where i is the index of the fundamental particle, j is the index of the nearest neighbor particle, N is the number of particles, t is the computation time, ρ is the particle density, m is the particle mass, v is the particle velocity, W is the smoothing kernel function, x, y, and z are the Cartesian components of the particle displacement, and α and β represent the Cartesian components x, y, and z. It is the total stress on soil particles. It is the Dirac function. It's artificial viscosity. It is artificial stress. It is the additional force on the fundamental particle, P is the isotropic pressure, and S is the anisotropic shear stress.

8. An electronic device, characterized in that, The electronic device includes: a processor, and a memory coupled to the processor. The memory is used to store computer programs; and The processor is configured to execute the computer program stored in the memory, so that the electronic device performs the simulation method for sliding of soil-rock mixed slope as described in any one of claims 1-6.

Images