A discrete element-based in-situ ground stress equilibrium method for rock and soil analysis

Abstract

The present invention discloses a method for in-situ balance of geostress based on discrete elements for rock and soil analysis, and relates to the technical field of numerical simulation of geology and geotechnical engineering. The method comprises: S1, data acquisition step: obtaining corresponding material parameters, control parameters and boundary conditions of the target area, and establishing a particle discrete model; S2, initialization step: initializing the particle discrete model obtained in S1; S3, data analysis step: calculating the contact force between particles and between particles and boundaries, and updating particle velocity and position information; S4, data judgment step: judging whether the number of executed time steps meets the condition of resetting the particle position, and obtaining the in-situ balance result of the geostress of the numerical model. The present invention can realize the initial geostress balance problem in particle discrete element and block discrete element simulation.

Description

Technical Field

[0001] The present invention relates to the technical field of numerical simulation of geology and geotechnical engineering, and in particular to a discrete element-based in-situ ground stress balance method applied to rock and soil body analysis. Background Art

[0002] Geostress is the natural stress that exists in the earth's crust without being disturbed by engineering. In numerical simulations in the field of geology and geotechnical engineering, initial geostress balance is the first step in the simulation solution. The accuracy of the balance result determines the accuracy of subsequent analysis results. The initial geostress balance is to obtain a state where initial stress exists but initial strain does not exist. In general, when modeling, the size of the initially established geometric model and the target geometric model are consistent. At this time, the model is in a state without initial geostress, and gravity and boundary conditions need to be applied for balance calculation. The size of the model after force will be different from the target model, which will have a certain impact on subsequent analysis. Therefore, it is very necessary to balance geostress in some numerical simulations. In current research, numerical simulation methods based on continuous media have become increasingly mature and widely used, such as finite element method, finite volume method, finite difference method, etc.; for the initial geostress balance analysis in the simulation process, many common finite element software on the market involve the solution of the problem, such as Abaqus, Comsol, midas GTS NX, etc., among which Abaqus is the most widely used. The methods provided in this software to balance the initial geostress include automatic balance method, keyword definition initial geostress method, ODB import method and initial geostress extraction method. The automatic balance method can realize the geostress balance function by pre-specifying the allowable displacement change tolerance and automatically selecting the incremental step. This method is simple and easy to operate; the keyword definition initial geostress method is suitable for the situation where the geometric model is relatively regular and the surface is horizontal. It can apply a certain initial stress field by considering the lateral pressure coefficients in the two horizontal directions; the principle of the ODB import method is to calculate an ODB file of initial geostress in advance, and import the file when performing geostress balance; the initial geostress extraction method is to extrapolate the calculated stress to the centroid and then derive the stress components in each direction.

[0003] The above method using the finite element principle can solve the problem of geostress balance based on the continuous simulation process, but when the numerical simulation process involves large deformation and discontinuous deformation damage, such as large-scale landslides, tunnel excavation, ground subsidence, etc., the above method is difficult to be well implemented. The discrete element method based on discontinuous media is more suitable for simulating such phenomena. When using the discrete element method to simulate settlement deformation damage and other problems, the initial geostress balance of the model is also required. At present, the solution to such problems in the discrete element method is usually to directly apply deadweight after the model is established, so that it can settle and balance under the action of gravity. This treatment method will also cause the model size to be different from the original size after gravity deposition, that is, there is a certain size of initial displacement before the simulation. In the analysis of deformation simulation and other problems, such as foundation pit excavation, the initial displacement before the simulation will be accumulated in the subsequent construction simulation, and additional additional displacement will be obtained, making the simulation result unsatisfactory. Therefore, it is crucial to study the initial geostress balance method based on discrete elements.

[0004] Therefore, it is an urgent problem for those skilled in the art to propose an in-situ ground stress balance method based on discrete elements for rock and soil analysis to solve the difficulties existing in the prior art. Summary of the invention

[0005] In view of this, the present invention provides an in-situ geostress balance method based on discrete elements for rock and soil analysis, which can realize the initial geostress balance problem in particle discrete element and block discrete element simulation.

[0006] In order to achieve the above object, the present invention adopts the following technical solution:

[0007] A discrete element-based in-situ ground stress balance method for rock and soil analysis includes the following steps:

[0008] S1, data acquisition step: obtaining corresponding material parameters, control parameters and boundary conditions of the target area, and establishing a particle discrete model;

[0009] S2, initialization step: initializing the particle discrete model obtained in S1, including: setting the number of time steps for resetting the particle position, detecting and recording the initial position of the particle, and applying gravity conditions to the numerical model;

[0010] S3, data analysis step: based on the initialization results in S2, calculate the contact forces between particles and between particles and boundaries, and update the particle velocity and position information;

[0011] S4, data judgment step: judge whether the number of executed time steps meets the condition of particle position reset. If so, reset the particle position information; if not, continue the calculation and iterate repeatedly until the convergence condition is reached and the calculation is terminated to obtain the in-situ equilibrium result of the ground stress of the numerical model.

[0012] The above method, optionally, the material parameters in S1: including geometric parameters and mechanical parameters;

[0013] Among them, the geometric parameters include particle position coordinate information, particle shape and particle size;

[0014] Mechanical parameters include particle density ρ, normal contact stiffness K n , tangential contact stiffness K s , Maximum contact tension F t , contact adhesion F c Angle of internal friction

[0015] The control parameters include: local damping coefficient, total simulation time step length, and time step reduction factor;

[0016] There are two types of boundary conditions: one is a fixed boundary condition, which limits the horizontal and vertical displacements of the side and bottom boundaries; the other is a viscoelastic boundary condition, which sets normal and tangential springs and dampers on the side and bottom boundaries.

[0017] In the above method, optionally, the particle discrete model content in S1 includes:

[0018] (1) The rock and soil in the study area is discretized into densely packed particles, and the particles in contact with each other are connected by a bonded particle model;

[0019] (2) A polygon is created based on the actual shape and size of the model, and the polygon information is imported into the program as the model boundary. Based on this polygon, a model consisting of a finite number of randomly or regularly closely arranged discrete particles and walls is created.

[0020] In the above method, optionally, the number of time steps for resetting the particle position in S2 includes:

[0021] Use trial and error to determine the number of time steps required to perform particle position resetting during the static analysis.

[0022] In the above method, optionally, detecting and recording the initial position of the particle in S2 includes:

[0023] Contact detection is performed between particles and between particles and boundaries, and the initial positions of all particles in the numerical model are recorded. In contact detection, the PushBufferContact contact buffer strategy is used, and contact detection is restarted when the maximum motion in the particles exceeds the set cutoff distance.

[0024] In the above method, optionally, applying gravity conditions to the numerical model in S2 includes:

[0025] Gravity conditions are applied to the numerical model so that the model can perform conventional static analysis under its own weight. During the static analysis of the numerical model, the loads and boundary conditions applied in each iteration remain consistent.

[0026] In the above method, optionally, the calculation of contact forces between particles and between particles and boundaries in S3 includes:

[0027] Normal force F n With tangential force F s The calculations all follow Hooke's law, as shown in formula (1)(2);

[0028] ΔF n =-K n Δd n (1)

[0029] ΔF s =-K s Δd s (2)

[0030] Where K n , K s are the normal contact stiffness and the tangential contact stiffness respectively; Δd n , Δd s are the changes in the normal and tangential relative displacements of the contact pair in each calculation time step;

[0031] Particle velocity update:

[0032] Based on the updated acceleration, the particle velocity can be updated as shown in formulas (3) to (5);

[0033]

[0034]

[0035] Where i is the particle number, Δt is the time step increment, and m i and I i are the mass and moment of inertia of particle i, is the force on the particle in the x direction at the current moment, It represents the velocity of the particle in the x direction at time t-Δt / 2. The central difference method is used here. The definitions in the y direction and the rotation direction are the same.

[0036] Based on the velocity update, the particle position state update can be obtained, as shown in formulas (6) to (8);

[0037]

[0038] In the formula, x i (t+Δt) and x i (t) They represent the displacement of the particle in the x direction at time t+Δt and time t, respectively. It represents the velocity of the particle in the x direction at time t+Δt / 2. The definitions in the y direction and the rotation direction are similar.

[0039] In the above method, optionally, the specific contents of S4 include:

[0040] Based on the number of time steps required to reset the particle positions during the static analysis process set in S2, it is determined whether the currently executed time step is an integer multiple of the number of time steps required to reset the particle positions;

[0041] If the current time step meets the preset conditions, it is necessary to reset the positions of all particles in the current step to the recorded initial positions of the particles just after updating the force, acceleration and velocity of the current time step, keep the contact force between the particles unchanged, and continue with the next iterative calculation.

[0042] The above method is optional. The sign of judging the end of calculation is whether the convergence condition is met. The convergence condition is as follows:

[0043]

[0044] In the formula, is the velocity of particle i at the current time step, where v ix t and v iy t are the x-direction velocity and u-direction velocity of particle i in the current time step, Δt is the time step increment, R i is the radius of particle i, and the recommended value of ε is 10 -8 .

[0045] If the current calculation step does not meet the convergence conditions, it is necessary to continue the iterative calculation. If it meets the conditions, the calculation ends and the initial geostress equilibrium model of the numerical model is obtained.

[0046] It can be seen from the above technical solution that, compared with the prior art, the in-situ ground stress balance method based on discrete elements applied to rock and soil analysis of the present invention has the following beneficial effects:

[0047] Based on the discrete element method, the initial geostress balance is performed in situ. The model obtained after the balance can achieve consistency between the simulated geometric model and the actual geometric model, that is, the initial geostress field is applied to the numerical model in the state of "zero" displacement field to avoid affecting the subsequent deformation simulation analysis. It is applied to particle discrete element and block discrete element simulation, and has significant advantages in non-continuous and large deformation failure problems in the field of geotechnical engineering. Among them, in the in-situ geostress balance method, since the topological relationship between particles basically does not change and the relative displacement generated is very small, there is basically no need to update the contact retrieval when performing the balance, which greatly shortens the balance calculation time. BRIEF DESCRIPTION OF THE DRAWINGS

[0048] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on the provided drawings without paying creative work.

[0049] Figure 1 A flow chart of a discrete element based in-situ ground stress balance method for rock and soil analysis provided by the present invention;

[0050] Figure 2 The present invention provides a schematic diagram of the particle contact model structure, wherein 2a is a schematic diagram of a normal contact spring, and 2b is a schematic diagram of a tangential contact spring;

[0051] Figure 3 Schematic diagram of the bonding contact constitutive model provided by the present invention, wherein 3a is the normal bonding contact constitutive model, and 3b is the tangential bonding contact constitutive model;

[0052] Figure 4 Schematic diagram of the Mohr-Coulomb constitutive model provided by the present invention, wherein 4a is the normal Mohr-Coulomb constitutive model, and 4b is the tangential Mohr-Coulomb constitutive model;

[0053] Figure 5 A schematic diagram of an initial slope model provided by the present invention;

[0054] Figure 6 The slope direction displacement map provided by the present invention;

[0055] Figure 7 A stress curve diagram of selected observation points developed over time provided by the present invention;

[0056] Figure 8 A graph of directional stress variation over time provided by the present invention;

[0057] Fig. 9 A comparison chart of the calculated solution and the analytical solution provided by the present invention;

[0058] Fig.10 The contact force chain diagram provided by the present invention;

[0059] Fig.11 A slope direction displacement diagram after applying gravity provided by the present invention;

[0060] Fig.12 A stress curve diagram of the selected observation point developing over time after applying gravity provided by the present invention;

[0061] Fig.13 A graph showing the change of directional stress over time after gravity is applied provided by the present invention;

[0062] Fig.14 A comparison diagram of the calculated solution and the analytical solution after applying gravity provided by the present invention;

[0063] Fig.15 The contact force chain diagram after gravity is applied provided by the present invention. DETAILED DESCRIPTION

[0064] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0065] In this application, relational terms such as first and second, etc. are used only to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any such actual relationship or order between these entities or operations. The terms "include", "comprise" or any other variants thereof are intended to cover non-exclusive inclusion, so that a process, method, article or device including a series of elements includes not only those elements, but also other elements not explicitly listed, or also includes elements inherent to such process, method, article or device. In the absence of more restrictions, the elements defined by the sentence "comprise one..." do not exclude the existence of other identical elements in the process, method, article or device including the elements.

[0066] The present invention can be used in many general or special computing device environments or configurations, such as personal computers, server computers, handheld or portable devices, tablet devices, multi-processor devices, distributed computing environments including any of the above devices or devices, etc.

[0067] Reference Figure 1 As shown, the present invention discloses a discrete element based in-situ ground stress balance method for rock and soil analysis, comprising the following steps:

[0068] S1, data acquisition step: obtaining corresponding material parameters, control parameters and boundary conditions of the target area, and establishing a particle discrete model;

[0069] S2, initialization step: initializing the particle discrete model obtained in S1, including: setting the number of time steps for resetting the particle position, detecting and recording the initial position of the particle, and applying gravity conditions to the numerical model;

[0070] S3, data analysis step: based on the initialization results in S2, calculate the contact forces between particles and between particles and boundaries, and update the particle velocity and position information;

[0071] S4, data judgment step: judge whether the number of executed time steps meets the condition of particle position reset. If so, reset the particle position information; if not, continue the calculation and iterate repeatedly until the convergence condition is reached and the calculation is terminated to obtain the in-situ equilibrium result of the ground stress of the numerical model.

[0072] Furthermore, the material parameters in S1 include geometric parameters and mechanical parameters;

[0073] Among them, the geometric parameters include particle position coordinate information, particle shape and particle size;

[0074] Mechanical parameters include particle density ρ, normal contact stiffness K n , tangential contact stiffness K s , Maximum contact tension F t , contact adhesion F c Angle of internal friction

[0075] The control parameters include: local damping coefficient, total simulation time step length, and time step reduction factor;

[0076] There are two types of boundary conditions: one is a fixed boundary condition, which limits the horizontal and vertical displacements of the side and bottom boundaries; the other is a viscoelastic boundary condition, which sets normal and tangential springs and dampers on the side and bottom boundaries.

[0077] Furthermore, the particle discrete model in S1 includes:

[0078] (1) The rock and soil in the study area is discretized into densely packed particles, and the particles in contact with each other are connected by a bonded particle model;

[0079] (2) A polygon is created based on the actual shape and size of the model, and the polygon information is imported into the program as the model boundary. Based on this polygon, a model consisting of a finite number of randomly or regularly closely arranged discrete particles and walls is created.

[0080] Specifically, the numerical model of the simulation is to discretize the rock and soil into individual particles, and then stack and cement the particles, which can well simulate the discontinuity, heterogeneity and large deformation failure problems of the rock and soil. The particle model is to establish polygons according to the actual shape and size of the model, and import the polygon information into the program as the model boundary, and based on this polygon, establish a model composed of a finite number of randomly or regularly closely arranged discrete particles and walls. The model adopts the following assumptions:

[0081] 1) In the discrete element modeling process of continuous media, the study area is discretized into densely packed particles. The particles that are in contact with each other are connected by a bonded particle model.

[0082] 2) The model includes a wall model and a particle model. Normal and tangential springs are introduced between particles and between particles and the wall to generate contact forces, such as Figure 2 As shown in the figure, 2a is a schematic diagram of the normal contact spring, and 2b is a schematic diagram of the tangential contact spring. The normal force and tangential force between particles can be simulated by the normal spring and the tangential spring. Both the particles and the wall are simplified as rigid bodies, without deformation, and contain state attributes such as position, force, acceleration and velocity.

[0083] Furthermore, the number of time steps for resetting the particle position in S2 includes:

[0084] Use trial and error to determine the number of time steps required to perform particle position resetting during the static analysis.

[0085] Specifically, the number of time steps required to reset the particle positions during the static analysis is set, that is, the number of iterations required to reset the positions. The number of time steps must be small enough to ensure that the contact topological relationship between particles does not change during the deposition process under the action of gravity; the number of time steps cannot be too small to affect the calculation efficiency, which can be determined by trial and error.

[0086] Furthermore, detecting and recording the initial position of the particle in S2 includes:

[0087] Contact detection is performed between particles and between particles and boundaries, and the initial positions of all particles in the numerical model are recorded. In contact detection, the PushBufferContact contact buffer strategy is used, and contact detection is restarted when the maximum motion in the particles exceeds the set cutoff distance.

[0088] Specifically, after contact detection, force can only be transferred between particles that meet the contact conditions and between particles and boundaries. In contact detection, the PushBufferContact contact buffer strategy is used. It is not necessary to start contact detection at each time step. Contact detection is restarted only when the maximum movement in the particles exceeds the set cutoff distance. In particular, in this in situ stress balance method, since the topological relationship between particles basically does not change, the relative displacement generated is very small. This phenomenon will be explained later. Therefore, there is no need to update the contact retrieval when performing balance, which greatly shortens the balance calculation time.

[0089] The contact model is divided into the bonding contact model and the Mohr-Coulomb contact model, such as Figure 3 , Figure 4 As shown, 3a is the normal bonding contact constitutive model, 3b is the tangential bonding contact constitutive model, 4a is the normal Mohr-Coulomb constitutive model, and 4b is the tangential Mohr-Coulomb constitutive model. The one less than the contact threshold is the bonding model, which has a certain tensile strength and shear strength and can transmit contact force. This bonding model is not necessary when the ground stress is balanced, but may be needed in subsequent simulations. The contact buffer strategy uses a larger potential contact threshold to find potential contact pairs, namely Mohr-Coulomb contact, which will not be assigned strength and cannot withstand tension.

[0090] Furthermore, the gravity conditions imposed on the numerical model in S2 include:

[0091] Gravity conditions are applied to the numerical model so that the model can perform conventional static analysis under its own weight. During the static analysis of the numerical model, the loads and boundary conditions applied in each iteration remain consistent.

[0092] Furthermore, the calculation of contact forces between particles and between particles and boundaries in S3 includes:

[0093] Normal force F n With tangential force F s The calculations all follow Hooke's law, as shown in formula (1)(2);

[0094] ΔF n =-K n Δd n (1)

[0095] ΔF s =-Ks Δd s (2)

[0096] Where K n , K s are the normal contact stiffness and the tangential contact stiffness respectively; Δd n , Δd s are the changes in the normal and tangential relative displacements of the contact pair in each calculation time step;

[0097] Particle velocity update:

[0098] Based on the updated acceleration, the particle velocity can be updated as shown in formulas (3) to (5);

[0099]

[0100] Where i is the particle number, Δt is the time step increment, and m i and I i are the mass and moment of inertia of particle i, is the force on the particle in the x direction at the current moment, It represents the velocity of the particle in the x direction at time t-Δt / 2. The central difference method is used here. The definitions in the y direction and the rotation direction are the same.

[0101] Based on the velocity update, the particle position state update can be obtained, as shown in formulas (6) to (8);

[0102]

[0103]

[0104]

[0105] In the formula, x i (t+Δt) and x i (t) They represent the displacement of the particle in the x direction at time t+Δt and time t, respectively. It represents the velocity of the particle in the x direction at time t+Δt / 2. The definitions in the y direction and the rotation direction are similar.

[0106] Specifically, the movement of particles follows Newton's second law. During the entire calculation process, forces act between contacting particles and between particles and boundaries. Due to factors such as movement, expansion, overlap, and contact destruction of particles after being subjected to gravity, the contact force of particles is constantly changing, and the unbalanced force on the particles will cause the particles to accelerate.

[0107] Furthermore, the specific contents of S4 include:

[0108] Based on the number of time steps required to reset the particle positions during the static analysis process set in S2, it is determined whether the currently executed time step is an integer multiple of the number of time steps required to reset the particle positions;

[0109] If the current time step meets the preset conditions, it is necessary to reset the positions of all particles in the current step to the recorded initial positions of the particles just after updating the force, acceleration and velocity of the current time step, keep the contact force between the particles unchanged, and continue with the next iterative calculation.

[0110] Specifically, a method is innovatively proposed to apply the initial geostress field in the state of "zero" displacement field to the numerical model to avoid affecting the subsequent deformation simulation analysis. The basis for the realization of this method is as follows:

[0111] 1) The requirement for balancing the initial geostress is that the model obtained by applying gravity and performing equilibrium iteration under self-weight conditions has initial stress but the displacement is infinitely close to zero, that is, the size of the obtained model is consistent with the actual model size and there is geostress;

[0112] 2) The discrete element method used in this method follows Newton's second law. The calculation of force, acceleration, velocity, and displacement is relatively simple. The update of the normal contact force and tangential contact force between particles and between particles and the wall is in the form of increments, as shown in formulas (9) and (10).

[0113] F n t =F n t-Δt +ΔF n (9)

[0114] F s t =F s t-Δt +ΔF s (10)

[0115] In the formula, F n t and F n t-Δt are the normal contact forces at the current time step and the previous time step, ΔF n is the normal contact force increment in this time step, and the same is true for the shear direction. The calculation of the normal contact force increment and the tangential contact force increment is shown in the above formulas (1) and (2). The normal displacement increment Δd n and the tangential displacement increment Δd s The further calculation of is shown in formula (11)(12).

[0116] Δd n =vn Δt (11)

[0117] Δd s =v s Δt (12)

[0118] In the formula, v n and v s are the normal relative velocity and the tangential relative velocity respectively, and their calculations are shown in formulas (13) and (14).

[0119] v n =v x(ij) n x +v y(ij) n y (13)

[0120] v s =-v x(ij) n y +v y(ij) n x -θ i R i -θ j R j (14)

[0121] Where v x(ij) and v y(ij) are the relative velocities of particles i and j in the x and y directions, respectively, n x and n y are the unit vectors of particles i and j in the x and y directions respectively, θ i and θ j are the angular velocities of particles i and j, R i and R j are the radii of particles i and j respectively.

[0122] It can be seen from formulas (11) and (12) that during the iterative calculation of the contact force, only the current time step speed needs to be retained, which can be used to update the contact force increment of the next time step, and then update the total unbalanced force. After the force balance iteration is completed, a model with initial ground stress can be obtained for subsequent simulation analysis. In the whole process, the update of the contact force is determined by the relative speed updated in the previous time step. It can be seen from formulas (13) and (14) that the relative position of the particles has a certain influence on the speed update. Therefore, the number of time steps nstep set in S2 for resetting the particle position should not be too large to avoid changes in the particle contact topology relationship. If the relative position changes of the topological particles are small enough, the influence of the above problems can be ignored.

[0123] The specific operations are as follows:

[0124] By setting the number of time steps required to reset the particle position in the static analysis process in advance in S2, that is, the number of iterations required to reset the position, when performing subsequent iterative calculations, first check whether the number of time steps currently executed meets the following conditions:

[0125] if stepI / nstep is int? (15)

[0126] Where stepI represents the current time step, nstep is the number of time steps required to reset the position after each execution, and the recommended value is 10. int represents an integer.

[0127] If the current time step meets the preset conditions, it is necessary to reset the positions of all particles in the current step to the initial positions of the particles recorded in S4 just after the force, acceleration and velocity of the current time step are updated, keep the contact force between the particles unchanged, and continue with the next iterative calculation.

[0128] Furthermore, the sign of the end of the calculation is whether the convergence condition is met. The convergence condition is as follows:

[0129]

[0130] In the formula, is the velocity of particle i at the current time step, where v ix t and v iy t are the x-direction velocity and y-direction velocity of particle i in the current time step, Δt is the time step increment, R i is the radius of particle i, and the recommended value of ε is 10 -8 .

[0131] If the current calculation step does not meet the convergence conditions, it is necessary to continue the iterative calculation. If it meets the conditions, the calculation ends and the initial geostress equilibrium model of the numerical model is obtained.

[0132] In a specific embodiment, a simple two-dimensional slope model is used to verify the effectiveness of the initial geostress balance method. In the numerical calculation of the slope, the initial geostress balance under the action of deadweight is very important and is also a key step before the simulation calculation. If the initial geostress balance is not performed, the subsequent simulation analysis will be affected by the gravity, resulting in unsatisfactory results. Figure 5 Shown is an initial slope model before gravity is applied.

[0133] The microscopic parameters of the model are shown in Table 1.

[0134] Table 1 Model micromechanical parameters

[0135]

[0136] Apply downward gravity to all units of the model without performing initial ground stress balance operation. Perform 10,000 iterations until the calculation converges and all units reach a stable state. Figure 6 From the model diagram of the slope displacement in the y direction, it can be seen that the model diagram that only performs self-weight stress balance produces obvious vertical displacement and the geometric dimensions change;

[0137] During the simulation, 29 particle points with the same x-direction coordinates and equal spacing in the y-direction were selected to observe the stress in the y-direction, such as Figure 7 As shown in the figure, the corresponding stress curve developing with time is obtained, as shown in Figure 8 As shown, from Figure 8 It can also be seen that the calculation has reached convergence and the deadweight balance has been completed. The changes of the y-direction stress of the 29 observation points with the slope height are shown and compared with the analytical solution. Fig. 9 The resulting contact force chain is shown in Fig.10 shown.

[0138] The above results are all the results without initial geostress balance. From the above results, it can be seen that the displacement after applying gravity balance does not meet the conditions, so the initial geostress balance operation is performed next.

[0139] Similarly, using the above slope model, the downward gravity is applied to all units of the model, and the initial ground stress balance operation is performed. 10,000 steps of iteration are performed until the calculation converges and all units reach a stable state. The corresponding results are shown in the figure below. Figures 11 to 15 shown.

[0140] It can be seen from the above figure that after applying gravity for initial ground stress balance, the displacement change of the model meets the conditions. The stress of the 29 selected observation points gradually increases from top to bottom. The obtained stress simulation results are consistent with the analytical results. The contact force chain diagram is also consistent with the actual conditions, which shows the effectiveness of the in-situ ground stress balance method based on particle discrete elements.

[0141] Each embodiment in this specification is described in a progressive manner, and the same or similar parts between the embodiments can refer to each other, and each embodiment focuses on the differences from other embodiments. In particular, for the system or system embodiment, since it is basically similar to the method embodiment, the description is relatively simple, and the relevant parts can refer to the partial description of the method embodiment. The system and system embodiments described above are merely schematic, wherein the units described as separate components may or may not be physically separated, and the components displayed as units may or may not be physical units, that is, they may be located in one place, or they may be distributed on multiple network units. Some or all of the modules may be selected according to actual needs to achieve the purpose of the scheme of this embodiment. Ordinary technicians in this field can understand and implement it without creative work.

[0142] Those skilled in the art may further appreciate that the units and algorithm steps of each example described in conjunction with the embodiments disclosed herein may be implemented by electronic hardware, computer software, or a combination of both.

[0143] In order to clearly illustrate the interchangeability of hardware and software, the above description has generally described the components and steps of each example according to function. Whether these functions are performed in hardware or software depends on the specific application and design constraints of the technical solution. Professional and technical personnel can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of the present invention.

[0144] The above description of the disclosed embodiments enables one skilled in the art to implement or use the present invention. Various modifications to these embodiments will be apparent to one skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the present invention. Therefore, the present invention will not be limited to the embodiments shown herein, but rather to the widest scope consistent with the principles and novel features disclosed herein.

Claims

1. A discrete element based in situ ground stress balance method for rock and soil analysis, characterized in that: The following steps are involved: S1, data acquisition step: obtaining corresponding material parameters, control parameters and boundary conditions of the target area, and establishing a particle discrete model; S2, initialization step: initializing the particle discrete model obtained in S1, including: setting the number of time steps for resetting the particle position, detecting and recording the initial position of the particle, and applying gravity conditions to the numerical model; S3, data analysis step: based on the initialization results in S2, calculate the contact forces between particles and between particles and boundaries, and update the particle velocity and position information; S4, data judgment step: judge whether the number of executed time steps meets the condition of particle position reset. If so, reset the particle position information; if not, continue the calculation and iterate repeatedly until the convergence condition is reached and the calculation is terminated to obtain the in-situ equilibrium result of the ground stress of the numerical model.

2. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: Material parameters in S1: including geometric parameters and mechanical parameters; Among them, the geometric parameters include particle position coordinate information, particle shape and particle size; Mechanical parameters include particle density ρ, normal contact stiffness K n , tangential contact stiffness K s , Maximum contact tension F t , contact adhesion F c Angle of internal friction The control parameters include: local damping coefficient, total simulation time step length, and time step reduction factor; There are two types of boundary conditions: one is a fixed boundary condition, which limits the horizontal and vertical displacements of the side and bottom boundaries; the other is a viscoelastic boundary condition, which sets normal and tangential springs and dampers on the side and bottom boundaries.

3. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 2, characterized in that: The particle discrete model in S1 includes: (1) The rock and soil in the study area is discretized into densely packed particles, and the particles in contact with each other are connected by a bonded particle model; (2) A polygon is created based on the actual shape and size of the model, and the polygon information is imported into the program as the model boundary. Based on this polygon, a model consisting of a finite number of randomly or regularly closely arranged discrete particles and walls is created.

4. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: The number of time steps for resetting particle positions in S2 includes: Use trial and error to determine the number of time steps required to perform particle position resetting during the static analysis.

5. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: The initial position of the particle detected and recorded in S2 includes: Contact detection is performed between particles and between particles and boundaries, and the initial positions of all particles in the numerical model are recorded. In contact detection, the PushBufferContact contact buffer strategy is used, and contact detection is restarted when the maximum motion in the particles exceeds the set cutoff distance.

6. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: The gravity conditions imposed on the numerical model in S2 include: Gravity conditions are applied to the numerical model so that the model can perform conventional static analysis under its own weight. During the static analysis of the numerical model, the loads and boundary conditions applied in each iteration remain consistent.

7. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: The contact force calculation between particles and between particles and boundaries in S3 includes: Normal force F n With tangential force F s The calculations all follow Hooke's law, as shown in formulas (1) and (2); ΔF n =-K n Δd n (1) ΔF s =-K s Δd s (2) Where K n , K s are the normal contact stiffness and the tangential contact stiffness respectively; Δd n , Δd s are the changes in the normal and tangential relative displacements of the contact pair in each calculation time step; Particle velocity update: Based on the updated acceleration, the particle velocity can be updated as shown in formulas (3) to (5); Where i is the particle number, Δt is the time step increment, and m i and I i are the mass and moment of inertia of particle i, is the force on the particle in the x direction at the current moment, It represents the velocity of the particle in the x direction at time t-Δt / 2. The central difference method is used here. The definitions in the y direction and the rotation direction are the same. Based on the velocity update, the particle position state update can be obtained, as shown in formulas (6) to (8); In the formula, x i (t+Δt) and x i (t) They represent the displacement of the particle in the x direction at time t+Δt and time t, respectively. It represents the velocity of the particle in the x direction at time t+Δt / 2. The definitions in the y direction and the rotation direction are similar.

8. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: The specific contents of S4 include: Based on the number of time steps required to reset the particle positions during the static analysis process set in S2, it is determined whether the currently executed time step is an integer multiple of the number of time steps required to reset the particle positions; If the current time step meets the preset conditions, it is necessary to reset the positions of all particles in the current step to the recorded initial positions of the particles just after updating the force, acceleration and velocity of the current time step, keep the contact force between the particles unchanged, and continue with the next iterative calculation.

9. The in-situ ground stress balance method based on discrete element for rock and soil analysis according to claim 1, characterized in that: The sign of judging the end of calculation is whether the convergence condition is met. The convergence condition is as follows: In the formula, is the velocity of particle i at the current time step, where v ix t and v iy t are the x-direction velocity and y-direction velocity of particle i in the current time step, Δt is the time step increment, R i is the radius of particle i, and ε is 10 -8 ; If the current calculation step does not meet the convergence conditions, it is necessary to continue the iterative calculation. If it meets the conditions, the calculation ends and the initial geostress equilibrium model of the numerical model is obtained.

Images