Karst area rock sample penetration damage simulation method and system based on discrete elements

Abstract

The invention relates to a karst area rock sample penetration damage simulation method and system based on discrete elements, and belongs to the field of geotechnical engineering. The invention discloses a karst area rock sample penetration damage simulation method based on discrete elements. The karst area rock sample penetration damage simulation method comprises the following steps: S1, constructing a fracture damage constitutive model reflecting a karst area rock discontinuous surface position; s2, inputting parameters and conditions based on a discrete element method and the fracture damage constitutive model reflecting the position of the rock discontinuous surface in the karst area, and calculating particle stress information; s3, calculating the particle speed and displacement according to the particle stress information to obtain a calculation result; s4, outputting the calculation result; wherein in the step S2, according to the position relation between the particles and the discontinuous surface, the damage coefficient of the discontinuous surface is introduced, the particle stress information is dynamically corrected, and the penetration mechanical response under different discontinuous surface positions is accurately simulated. According to the method, the position and the crushing damage degree of the discontinuous surface of the rock are considered, a new rock sample model is constructed, and more accurate simulation is realized.

Description

Technical Field

[0001] This invention relates to the field of geotechnical engineering, and in particular to a method and system for simulating penetration damage of rock samples in karst areas based on discrete element method. Background Technology

[0002] Rocks in karst regions often contain numerous discontinuities (such as fissures), and their mechanical behavior is significantly affected by the location of these discontinuities. Existing discrete element method (DEM) simulations are mostly based on homogeneous rock models and do not consider the influence of discontinuity location on damage evolution, leading to discrepancies between simulation results and actual experiments. For example, prior art document 1 (CN119294217 B) discloses a DEM-based method and system for simulating rock penetration damage, which corrects particle stress through a staged constitutive model and damage values, but does not address the directional influence of discontinuity location on damage in karst regions.

[0003] In karst rock penetration tests, samples with the same GSI value may exhibit different peak penetration forces, resulting in significant dispersion in the test results and making simulation difficult. This is because karst rocks contain numerous discontinuities or fissures. The distance between the penetration point and these discontinuities or fissures varies, leading to substantial differences in the peak penetration force.

[0004] Therefore, there is an urgent need for a damage simulation method that can reflect the location characteristics of rock discontinuities in karst areas in order to achieve more accurate simulation. Summary of the Invention

[0005] Based on this, in order to solve the technical problem that "existing methods for simulating rock sample penetration tests in karst areas use purely linear constitutive models that are difficult to simulate actual mechanical states, and the calculated results of the constitutive models that do not consider the location of discontinuities differ from the actual results", this invention provides a method and system for simulating rock sample penetration damage in karst areas based on discrete element method. Specifically, it involves a method for constructing constitutive models that reflect the degree of rock damage and the distance between the penetration point and the discontinuity during penetration tests of rock samples in karst areas. This method can provide a reference for the mechanical characteristics of simulated rock sample penetration tests in karst areas, reflect the damage of the penetration test model under different degrees of damage and locations of discontinuities, and ensure production safety during construction.

[0006] To achieve the above objectives, the present invention adopts the following technical solution:

[0007] This invention provides a method for simulating penetration damage of rock samples in karst areas based on the discrete element method, comprising the following steps:

[0008] S1. Construct a constitutive model of fracture damage that reflects the location of rock discontinuities in karst areas;

[0009] S2. Based on the discrete element method and the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, input parameters and conditions to calculate particle stress information; wherein, S2 includes the following steps:

[0010] S21. Initialize the fracture damage constitutive model that reflects the location of rock discontinuities in karst areas;

[0011] S22, Input material parameters;

[0012] S23. Determine the coordinates of discrete particles;

[0013] S24. Number the discrete particles respectively and input the corresponding data information;

[0014] S25. Apply initial boundary conditions and determine the time step;

[0015] S26. Calculate the particle stress information according to the calculation rules of the fracture damage constitutive model that reflects the position of rock discontinuities in karst areas;

[0016] S27. Correct the forces acting on the particles. The formulas for calculating the normal and tangential forces of the particles are as follows:

[0017]

[0018] Among them, F nmodified For the corrected normal force, F smodified For the corrected tangential force, F n,pire For the uncorrected normal force, F s,pure The uncorrected tangential force is given by D, where D is the damage value of the particle. f The damage coefficient for the discontinuity is between 0 and 1, where R is the distance between the discontinuity and the probe point. b Where d is the particle radius BC The distance between the particle and the probe point;

[0019] S3. Calculate the particle velocity and displacement based on the particle force information to obtain the calculation results;

[0020] S4. Output the calculation results.

[0021] The discrete element method (DEM) is well-suited for problems involving discontinuous media and is mature in studying large displacement deformation and nonlinear relationships. It has been widely applied in geotechnical engineering, structural geology, geophysics, mining engineering, and other fields, yielding numerous significant results. For rock samples from karst areas, penetration tests reveal a trend where the penetration force varies with the degree of fragmentation and damage. Furthermore, for samples with similar degrees of fragmentation, the magnitude of the penetration force is directly related to the location of the discontinuity. Therefore, to better represent the mechanical properties and failure patterns of rock samples from karst areas during penetration tests, this invention considers the location and degree of fragmentation and damage of the rock discontinuity, constructing a new rock sample model to achieve more accurate simulation.

[0022] The present invention provides a discrete element method for simulating penetration damage of rock samples in karst areas. In step S1, a fracture damage constitutive model reflecting the location of discontinuities in the rock in the karst area is constructed to describe the stress state of particles in the penetration test model when subjected to penetration damage by a spherical probe. In step S27 of step S2, a damage coefficient for the discontinuity is introduced based on the positional relationship between the particles and the discontinuity, dynamically correcting the particle stress information and accurately simulating the penetration mechanical response at different discontinuity locations.

[0023] In some embodiments, the discrete element method for simulating rock sample penetration damage in karst areas according to the present invention further includes the following steps:

[0024] S5. In conjunction with the indoor penetration test scheme, the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas is verified.

[0025] This invention is based on the simulation experiment of the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, and compares the results with those obtained from actual experiments to verify the reliability of the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas.

[0026] This invention also proposes a discrete element method-based system for simulating the penetration damage of rock samples in karst areas, used to execute any of the discrete element method-based methods for simulating the penetration damage of rock samples in karst areas described above.

[0027] Furthermore, the karst area rock sample penetration damage simulation system of the present invention includes:

[0028] The constitutive model construction module is used to construct a fracture damage constitutive model that reflects the location of rock discontinuities in karst areas;

[0029] The first calculation module, based on the discrete element method and the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, takes input parameters and conditions to calculate the force information of particles.

[0030] The second calculation module is used to calculate the particle velocity and displacement based on the particle force information and obtain the calculation results.

[0031] The result output module is used to output the calculation results.

[0032] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0033] This invention provides a method and system for simulating the penetration damage of rock samples in karst areas based on discrete element method (DEM). It constructs a constitutive model of fracture damage reflecting the location of discontinuities in karst rocks, taking into account the impact of rock sample damage degree and discontinuity location on the failure state of the penetration test model. This model is then applied to DEM software to form a method describing the failure characteristics of rock penetration tests in karst areas. This method is particularly suitable for simulating the differences in mechanical response of rock samples in karst areas due to different discontinuity locations during penetration tests. It can provide a reference for the mechanical characteristics of simulated penetration tests of rock samples in karst areas, reflecting the failure of penetration test models under different damage degrees and discontinuity locations, thus ensuring production safety during construction. Attached Figure Description

[0034] Figure 1 The diagram shows the variation of normal and tangential forces with relative displacement increments in a constitutive model for fracture damage at the location of rock discontinuities in karst areas.

[0035] Figure 2 To reflect the influence of the damage value of the constitutive model on the normal and tangential forces at the location of rock discontinuities in karst areas;

[0036] Figure 3 A schematic diagram showing the application of the method to the location of rock discontinuities in karst areas;

[0037] Figure 4 Flowchart for calculating the constitutive model of fracture damage to reflect the location of rock discontinuities in karst areas;

[0038] Figure 5 A schematic diagram of the cutting tool used in the penetration test;

[0039] Figure 6 This is a simulation result diagram showing the distance between the discontinuity and the probe point when the distance is relatively small.

[0040] Figure 7 This is a simulation result diagram when the distance between the discontinuity and the probe point is moderate;

[0041] Figure 8 This is a simulation result diagram when the distance between the discontinuity and the probe point is large. Detailed Implementation

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

[0043] It should be noted that, for ease of description, only the parts relevant to the invention are shown in the accompanying drawings. Unless otherwise specified, the embodiments and features described herein can be combined with each other.

[0044] It should be understood that the terms "system," "apparatus," "unit," and / or "module" used in this application are a method of distinguishing different components, elements, parts, sections, or assemblies at different levels. However, if other terms can achieve the same purpose, they may be replaced by other expressions.

[0045] As indicated in this application and claims, unless the context clearly indicates otherwise, the words "a," "an," "a," and / or "the" are not specifically singular and may include the plural. Generally, the terms "comprising" and "including" only indicate the inclusion of expressly identified steps and elements, which do not constitute an exclusive list, and the method or apparatus may also include other steps or elements. An element defined by the phrase "comprising an..." does not exclude the presence of other identical elements in the process, method, product, or apparatus that includes the element.

[0046] In the description of the embodiments of this application, "a plurality of" refers to two or more. The terms "first" and "second" are used for descriptive purposes only and should not be construed as indicating or implying relative importance or implicitly specifying the number of indicated technical features. Thus, a feature defined as "first" or "second" may explicitly or implicitly include one or more of that feature.

[0047] Furthermore, flowcharts are used in this application to illustrate the operations performed by the system according to embodiments of this application. It should be understood that the preceding or following operations are not necessarily performed precisely in sequence. Instead, the steps can be processed in reverse order or simultaneously. Additionally, other operations can be added to these processes, or one or more steps can be removed from them.

[0048] The technical solution of the present invention will be described in full and clear below with reference to the accompanying drawings in the embodiments of the present invention.

[0049] Example 1

[0050] Please see Figure 1This is a schematic flowchart of an optional example of the discrete element method for simulating rock sample penetration damage in karst areas according to the present invention. The method is applicable to computer equipment, and the simulation method proposed in this embodiment may include, but is not limited to, the following steps:

[0051] S1. Construct a constitutive model of fracture damage that reflects the location of rock discontinuities in karst areas;

[0052] S2. Based on the discrete element method and the fracture damage constitutive model that reflects the location of rock discontinuities in karst areas, input parameters and conditions to calculate particle stress information;

[0053] S3. Calculate the particle velocity and displacement based on the particle force information to obtain the calculation results;

[0054] S4. Output the calculation results;

[0055] S5. In conjunction with the indoor penetration test scheme, the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas is verified.

[0056] In some feasible embodiments, S1 specifically includes:

[0057] A damage and fracturing constitutive model reflecting the location of discontinuities in rocks in karst areas is constructed to describe the mechanical state of rock particles when they are subjected to localized damage.

[0058] Rocks in karst regions, being brittle materials, exhibit different stages of failure rather than simple linear elastic failure. This invention defines a constitutive model for the fracture damage of rocks in karst regions. This model divides the process from stress to cohesive failure into an initial stage, a critical stage, and a high-strain stage. The normal and tangential stiffness are inconsistent in each stage to demonstrate the mechanical state of rocks in karst regions at different stages.

[0059] Secondly, during penetration tests, rock samples from karst regions exhibit localized damage, which leads to a decrease in the overall peak penetration strength. Therefore, adding the concept of damage to the constitutive model of karst rock samples can effectively control the attenuation of interparticle forces due to changes in the degree of damage, thus reflecting the overall strength reduction of the karst rock sample. The damage value is determined by the initial bond quantity between the particle and its attachments and the number of bonds that break during penetration. A decrease in the number of bonds on the particle increases the damage value, thus reducing the actual force experienced by the particle and reflecting the localized fracturing process of the rock sample.

[0060] In constructing the constitutive model, in the initial stage when strain is relatively small, rock materials in karst regions typically exhibit linear elastic behavior. For example... Figure 1As shown, for normal tensile force, it is assumed that the interparticle force in this stage is proportional to the increase in relative displacement between particles, and the slope (stiffness) is relatively large, reflecting the initial elastic modulus of the rock material. Let A be the cross-sectional area, and the normal contact stiffness coefficient in the initial stage be... F n For the normal force of the previous time step, Δδ n The magnitude of the normal cohesion force between particles during this stage is calculated as the relative normal displacement increment, according to the following formula:

[0061]

[0062] Among them, F n1 This represents the normal force at the current time step in the initial stage. The parallel bond is used because the mechanical properties of this model after bond failure are similar to those of the parallel bond model. Therefore, the parallel bond form is also used in the definition, hence the horizontal line represents the force and other parameters on the parallel bond before the bond breaks.

[0063] As the relative displacement continues to increase, exceeding a certain critical value, microscopic damage begins to occur within the rock, such as the propagation of microcracks; this process is called the critical stage. The rate of increase in normal force during this stage is typically faster than in the initial stage, manifested as a steeper curve slope. The changes in normal force during this stage are as follows:

[0064]

[0065] Among them, F n2 k represents the normal force at the current time step during the critical phase. n,2 δ is the normal contact stiffness coefficient during the critical stage. n,TH1 The first relative normal displacement threshold, This represents the normal force offset during the critical phase, used to ensure the continuity of the curve.

[0066] Upon entering the high-strain region, the rock material may undergo a certain degree of stress release or recompression; this process is known as the high-strain stage. The slope of this stage is typically smaller than that of the first two stages, indicating that the increase in the material's normal force gradually slows down.

[0067]

[0068] Among them, F n3 k represents the normal force at the current time step during the high-strain phase. n,3 δ is the normal contact stiffness coefficient during the high strain stage. n,TH2 This is the second relative normal displacement threshold. This represents the normal force offset during the high strain stage, and is also used to ensure the continuity of the curve.

[0069] For normal pressure, the relationship between normal force and relative displacement increment is considered linear:

[0070]

[0071] In this embodiment, the calculations are mainly performed using the three formulas corresponding to the initial stage, critical stage, and high strain stage. The normal pressure is selected by those skilled in the art based on actual needs. For example, CN119294217B omits this formula, which is a conventional technique in the field.

[0072] For tangential forces, such as Figure 1 As shown, calculations are performed separately for the initial stage, critical stage, and high strain stage, with the following formulas:

[0073] In the initial stage:

[0074]

[0075] Among them, F s1 k represents the tangential force at the current time step in the initial stage. s,1 Δδ represents the tangential contact stiffness coefficient in the initial stage. s This represents the relative tangential displacement increment;

[0076] During the critical phase:

[0077]

[0078] Among them, F s2 k represents the tangential force at the current time step during the critical phase. s,2 δ is the tangential contact stiffness coefficient in the initial stage. s,TH1 F is the first relative tangential displacement threshold. s,OF2 This represents the tangential force offset during the critical stage.

[0079] During the high-strain phase:

[0080]

[0081] Among them, F s3 k represents the tangential force at the current time step during the high-strain phase. s,3 δ is the tangential contact stiffness coefficient during the high strain stage. s,TH2 F is the second relative tangential displacement threshold. s,OF3 This represents the tangential force offset during the high-strain stage.

[0082] The constitutive model describes the mechanical state between particles at the particle level, and the relative displacement threshold s mentioned is... TH1 and s TH2 Parameter calibration is required based on the actual situation in order to fit the actual data curve.

[0083] In some feasible embodiments, S2 specifically includes:

[0084] The Discrete Element Method (DEM) is used to calculate the force and displacement states between discrete particles. Based on the force-displacement law and Newton's second equation, the DEM first updates the contact forces generated at the particle contact points using the force-displacement law, and then updates the particle position information using Newton's second equation, readjusting the contact relationships between particles. Each time step requires traversing all particles until the model reaches equilibrium or fails.

[0085] For the contact between particles, there exists a relationship:

[0086]

[0087] in, The normal force between the i-th pair of particles, K represents the increment of the tangential force between the i-th pair of particles. n and k s These are the normal contact stiffness coefficient and the tangential contact stiffness coefficient, respectively, n i Let U be the contact normal vector between the i-th pair of particles. n Let ΔU be the relative displacement in the contact normal direction. s This represents the relative displacement in the tangential direction of contact.

[0088] When the particle has not reached equilibrium, its average acceleration and rotational acceleration are as follows:

[0089]

[0090] in, Let be the average acceleration at time t0. Let F be the rotational acceleration at time t0. x and M x These represent the unbalanced force and unbalanced torque in the horizontal direction, respectively. x It is angular momentum.

[0091] Integrating the above equation using the forward difference method, we can calculate... The velocity and rotational speed of the particle in the horizontal direction:

[0092]

[0093] in, Let ω be the velocity of the particle at time t1 along the horizontal axis. x (t1) represents the rotational speed of the particle at time t1 along the horizontal axis.

[0094] In this embodiment, based on the force results of the particle, Newton's second law is used to calculate the particle velocity v:

[0095]

[0096] Where f is the final calculated force on the particle, and d m Δt is the damping coefficient, Δt is the unit time step, and m is the unit mass of the particle. f is obtained from the previously obtained F. n and F s It is obtained by calculating the component forces in the x and y directions and then taking the resultant force. The calculation method is as follows:

[0097]

[0098] Among them, l n F represents the number of bonds formed by the particles. n and F s These are the normal and tangential forces of the bond, θ. l The angle of inclination for bonding.

[0099] f is f x and f y The combined force.

[0100] Based on the particle velocity results, the relative particle displacement Δs at each time step is obtained by multiplying the particle velocity by the unit time step Δt:

[0101] Δs=vΔt

[0102] Simulations were performed using discrete element method (DEM) software. This software utilizes C++ in Visual Studio as its development platform, and the graphical interface is implemented using an OpenGL library based on C++. The software simulation is implemented through multiple modules, including modeling, contact determination, constitutive model application, force-displacement calculation, wall loading, and historical monitoring. Each module is interconnected but can also be edited independently, facilitating adaptation to different model conditions and enabling the solution of various problems.

[0103] To address the mechanical properties of rocks in karst regions, this invention uses the Hoek-Brown failure criterion to determine the failure state between particles.

[0104] The Hoek-Brown failure criterion is an empirical criterion derived by fitting curves from triaxial rock test data. It reflects the nonlinear empirical relationship between the ultimate principal stresses at rock failure. The specific formula for the Hoek-Brown failure criterion is as follows:

[0105]

[0106] Where, σ1 ′ and σ3′ These represent the major and minor principal stresses at failure under triaxial testing conditions; σ uni m is the uniaxial compressive strength of the rock. b , s, and α are all rock strength parameters.

[0107] Rock strength parameter m b , α depends on the rock's geological strength index GSI, disturbed silver D, and intact frictional strength parameter m. i The parameters are defined as follows:

[0108]

[0109] The Geological Intensity Index (GSI) is a parameter used to evaluate the integrity of rocks; a higher GSI value indicates a more intact rock.

[0110] After defining the basic mechanical properties of the constitutive model, the concept of damage value is added to reflect the mechanical state of rock fracturing in karst areas. In the discrete element method, a damage value is assigned to each particle. This value is determined by the bonding state of the current particle and surrounding particles, and its magnitude is between 0 and 1. This value is a dimensionless quantity used to describe the degree of damage to the material in a certain state. d = 0 represents the material being in an undamaged state, and D = 1 represents the material being completely destroyed. The damage value D is calculated as follows:

[0111]

[0112] Where n0 is the number of bonds generated by the particle in the initial state of the model, and n1 is the number of bonds remaining in the particle at the current time step. The larger the damage value, the greater the degree of damage around the particle. Both the force on the particle's bonds and the failure are affected by this. When the damage value of the particle is relatively large, the force on its corresponding bonds will decrease, and the peak strength at failure will also decrease.

[0113] In the constitutive model of rock damage and fragmentation, the magnitude of the force on the particles after damage is calculated according to the different degrees of particle damage, such as... Figure 2 As shown, taking normal force as an example, the normal force generated by damaged particles will be lower than that of undamaged particles to a certain extent. In damage mechanics, it is assumed that the mechanical behavior of a material can be described by an effective undamaged material, and the stress of this effective material is called the "effective stress". The relationship between effective stress and actual stress is usually described by damage variables. In the calculation of the magnitude of the normal force between particles, the interparticle cohesion will be multiplied by a decay factor, which is determined by the damage value. The original normal force is defined as F. n Let the damage value be D, then the normal force F generated by the particle after damage... nmodified It will be calculated by the following formula:

[0114]

[0115] Among them, F nmodified For the corrected normal force, F smodified For the corrected tangential force, F n,pure For the uncorrected normal force, F s,pure This is the uncorrected tangential force.

[0116] At each time step, the damage value and force of the particles are dynamically calculated, and the normal force is adjusted according to the damage value. Stress calculation is performed using the dynamically updated normal and tangential forces. If the damage values ​​of two particles corresponding to the bond are inconsistent, the average of the two values ​​is used to balance the damage effects of the two particles.

[0117]

[0118] D1 and D2 represent the damage values ​​of two different particles, respectively, to reasonably reflect the average damage effect of the interaction between the two particles, making the model more physically meaningful and reasonable.

[0119] When the particle damage value changes, the stress-strain curve will show a sudden change. During the penetration test, if large-scale bond failure occurs in a certain area, the average particle damage value in that area will increase sharply, thereby causing a decrease in the stress of the specimen.

[0120] However, in karst rock penetration tests, samples with the same GSI value may exhibit different peak penetration forces, resulting in significant dispersion in the test results and making simulation difficult. This is because karst rocks contain numerous discontinuities or fractures. The distance between the penetration point and these discontinuities or fractures varies, leading to substantial differences in the peak penetration force. Generally, a smaller peak value corresponds to a smaller distance between the penetration point and the discontinuity, while a larger peak value corresponds to a larger distance.

[0121] To address this phenomenon, the location of discontinuous surfaces needs to be considered during the modeling process. For example... Figure 3 As shown, if a circle is drawn with the penetration point as the center, and the discontinuities are distributed on the circle, then the larger the radius of the circle, the greater the distance between the penetration point and the discontinuities, and the greater the peak penetration force obtained from the penetration test will be. This can effectively reflect the influence of the position of the discontinuities on the test results.

[0122] This invention defines the region on the circle as an area with a high degree of fragmentation, thereby increasing the particle damage value in this region. Considering the location of discontinuities, the calculation methods for the normal and tangential forces of the particles are as follows:

[0123]

[0124] Among them, Fnmodified For the corrected normal force, F smodified For the corrected tangential force, F n,pure For the uncorrected normal force, F s,pure For the uncorrected tangential force, D f The damage coefficient for the discontinuity is between 0 and 1, where R is the distance between the discontinuity and the probe point. b Where d is the particle radius BC The distance between the particle and the probe point.

[0125] In summary, in some feasible embodiments, S2 specifically includes:

[0126] S21. Based on the constitutive model of fracture damage reflecting the location of rock discontinuities in karst areas, initialize the parameters of the numerical model.

[0127] S22. Input material parameters and determine model dimensions;

[0128] Taking a planar problem as an example, it is necessary to determine the length, width, elastic modulus, density, damping coefficient, etc. of the model;

[0129] S23. Determine the coordinates of discrete particles and construct a numerical model;

[0130] S24. Number the discrete particles and input the corresponding data information;

[0131] S25. Apply initial boundary conditions and determine the time step;

[0132] The bottom of the sample is a fixed boundary, while the rest are free boundaries. Too short a time step will result in excessive instantaneous displacement of the particles, while too long a time step will result in low computational efficiency. Therefore, a suitable time step value needs to be selected.

[0133] Determine if the time step has been reached; if the time step has not been reached, proceed to: S26, calculate the particle stress according to the calculation rule of the fracture damage constitutive model reflecting the position of the discontinuity surface of rock in karst areas; if the time step has been reached, proceed directly to: S4, output the calculation results.

[0134] After S26 is completed, S27 is also included to correct the force on the particles, using the formula described above:

[0135]

[0136] In some feasible embodiments, S3 specifically includes:

[0137] After calculating the forces on the particles according to the constitutive model calculation rules, the program determines whether it has traversed all particles. If it has not traversed all particles, the forces on the particles are recalculated according to the constitutive model calculation rules. If it has traversed all particles, the particle velocity and displacement are calculated.

[0138] In some feasible embodiments, S5 specifically includes:

[0139] To verify the simulation accuracy of the developed software and the rationality of the constitutive model, numerical simulation data and experimental data were compared. Based on the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, simulated rock penetration tests were conducted, and the results were compared with those obtained from actual tests to verify the reliability of the constitutive model.

[0140] The tests were conducted using an MTS-E45.305 electronic universal testing machine. The equipment consists of a mechanical frame, servo motor system, testing system, and automatic control system, capable of conducting force, displacement, or strain controlled experiments within a force range of 300 kN, with a displacement measurement accuracy of 0.001 mm. During the tests, the penetration load and depth were measured and recorded in real time by built-in sensors, greatly facilitating test operation and control. All tests were conducted in displacement control mode, with a set piston speed of 0.5 mm / min, until the rock sample was completely destroyed. Throughout the process, important data such as penetration force and tool displacement were recorded in real time.

[0141] The penetration tool used is a ball-shaped tool made of high-strength steel, such as... Figure 5 As shown, the tip of the spherical cutting tool is machined to form a spherical tip with a radius of 5mm, and the side is ground into a 60° cone.

[0142] In this experiment, the penetration depth was taken as the tool displacement measured by the electronic universal testing machine, and the penetration pressure was the ratio of the penetration force to the projected area of ​​the tool penetrating the rock. The penetration pressure P i,exp It can be obtained from the following formula:

[0143]

[0144] Where F is the penetration force, R ind Let d be the tip radius of the spherical cutting tool, and d be the penetration depth.

[0145] The rock samples were obtained from the face of the tunnel during the advanced drilling construction of highway tunnels in karst areas. The rock samples were cut and processed into a cube shape of 100mm×100mm×50mm with a GSI value of 90.

[0146] For the distances from different discontinuities to the penetration point, the experimental and simulation results are as follows: Figures 6-8As shown. In the initial penetration stage, the penetration force exhibits a basically linear increasing trend with the increase of penetration depth until it reaches a maximum value, and then there will be a sudden drop, that is, the rock undergoes macroscopic brittle fracture. In some tests, a local drop in penetration force occurs before reaching the maximum value. Through observation of the test process, we found that the reason for this phenomenon is that the rock sample has local small rock fragments fracturing. However, since the rock sample is not completely destroyed at this time, the test will continue until macroscopic tensile cracks are generated, eventually leading to complete brittle fracture of the sample.

[0147] contrast Figure 6 , Figure 7 , Figure 8 It can be seen that the greater the distance between the discontinuity and the penetration point, the greater the peak value of the penetration pressure, which is consistent with the experimental results and can better reflect the influence of the location of the rock discontinuity and the degree of fracture damage on the sample strength in karst areas.

[0148] A discrete element method (DEM)-based system for simulating rock sample penetration damage in karst areas is provided for executing a DEM-based method for simulating rock sample penetration damage in karst areas. Specifically, the DEM-based system for simulating rock sample penetration damage in karst areas includes:

[0149] The constitutive model construction module is used to construct a fracture damage constitutive model that reflects the location of rock discontinuities in karst areas;

[0150] The first calculation module, based on the discrete element method and the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, takes input parameters and conditions to calculate the force information of particles.

[0151] The second calculation module is used to calculate the particle velocity and displacement based on the particle force information and obtain the calculation results.

[0152] The result output module is used to output the calculation results.

[0153] The content of the above method embodiments is applicable to this system embodiment. The specific functions implemented in this system embodiment are the same as those in the above method embodiments, and the beneficial effects achieved are also the same as those achieved in the above method embodiments.

[0154] A discrete element method (DEM) simulation device for rock penetration damage in karst areas includes:

[0155] At least one processor;

[0156] At least one memory for storing at least one program;

[0157] When the at least one program is executed by the at least one processor, the at least one processor implements a rock penetration damage simulation method based on discrete elements as described above.

[0158] The content of the above method embodiments is applicable to the device embodiments. The specific functions implemented by the device embodiments are the same as those of the above method embodiments, and the beneficial effects achieved are also the same as those achieved by the above method embodiments.

[0159] A storage medium storing processor-executable instructions, which, when executed by a processor, are used to implement a discrete element method for simulating rock penetration damage as described above.

[0160] The content of the above method embodiments is applicable to this storage medium embodiment. The specific functions implemented in this storage medium embodiment are the same as those in the above method embodiments, and the beneficial effects achieved are also the same as those achieved in the above method embodiments.

[0161] The embodiments described above are merely examples of several implementations of the present invention, and while the descriptions are relatively specific and detailed, they should not be construed as limiting the scope of the invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and the present invention also intends to include these modifications and variations.

Claims

1. A method for simulating the damage of rock sample penetration in karst area based on discrete elements, characterized in that: Includes the following steps: S1. Construct a constitutive model of fracture damage that reflects the location of rock discontinuities in karst areas; S2. Based on the discrete element method and the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, input parameters and conditions to calculate particle stress information; wherein, S2 includes the following steps: S21. Initialize the fracture damage constitutive model that reflects the location of rock discontinuities in karst areas; S22, Input material parameters; S23. Determine the coordinates of discrete particles; S24. Number the discrete particles respectively and input the corresponding data information; S25. Apply initial boundary conditions and determine the time step; S26. Calculate the particle stress information according to the calculation rules of the fracture damage constitutive model that reflects the position of rock discontinuities in karst areas; S27. Correct the forces acting on the particles. The formulas for calculating the normal and tangential forces of the particles are as follows: where F nmodified is the corrected normal force, F smodified is the corrected tangential force, F n,pure is the uncorrected normal force, F s,pure is the uncorrected tangential force, D is the damage value of the particle, D f is the damage coefficient of the discontinuity surface, which is between 0 and 1, R is the distance between the discontinuity surface and the probe point, r b is the radius of the particle, d BC is the distance between the particle and the probe point; S3. Calculate the particle velocity and displacement based on the particle force information to obtain the calculation results; S4. Output the calculation results. 2.The karst region rock sample penetration damage simulation method based on discrete elements according to claim 1, characterized in that: The constitutive model of fracture damage reflecting the location of rock discontinuities in karst areas is divided into an initial stage, a critical stage, and a high-strain stage. 3.The karst region rock sample penetration damage simulation method based on discrete elements according to claim 2, characterized in that: In the rock damage and fracture model, the change in normal force within the parallel key is as follows: In the initial stage: where F n1 is the normal force at the current time step in the initial phase; F n is the normal force at the previous time step, k n,1 is the normal contact stiffness coefficient in the initial phase, A is the cross-sectional area, Δδ n is the relative normal displacement increment, is a parallel key; During the critical phase: where F n2 is the normal force at the current time step in the critical phase, k n,2 is the normal contact stiffness coefficient in the critical phase, δ n,TH1 is the first relative normal displacement threshold, F n,OF2 is the normal force offset in the critical phase; During the high-strain phase: where F n3 is the normal force at the current time step in the high strain phase, k n,3 is the normal contact stiffness coefficient in the high strain phase, δ n,TH2 is a second relative normal displacement threshold, F n,OF3 is the normal force offset in the high strain phase.

4. The discrete element-based simulation method of karst region rock sample penetration damage according to claim 2, characterized in that: In the rock damage and fracture model, the variation of the internal tangential force of the parallel bond is as follows: In the initial stage: where F s1 is the tangential force at the current time step in the initial phase, k s,1 is the tangential contact stiffness coefficient in the initial phase, Δδ s is the relative tangential displacement increment; During the critical phase: where F s2 is the tangential force at the current time step in the critical phase, k s,2 is the tangential contact stiffness coefficient in the initial phase, δ s,TH1 is the first relative tangential displacement threshold, F s,OF2 is the tangential force offset in the critical phase; During the high-strain phase: where F s3 is the tangential force at the current time step in the high strain phase, k s,3 is the tangential contact stiffness coefficient in the high strain phase, δ s,TH2 is a second relative tangential displacement threshold, F s,OF3 is the tangential force offset in the high strain phase.

5. The discrete element-based simulation method of karst region rock sample penetration damage according to claim 1, characterized in that: Also includes: S5. In conjunction with the indoor penetration test scheme, the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas is verified.

6. A discrete element method (DEM) simulation system for simulating penetration damage of rock samples in karst areas, characterized in that: Used to perform the discrete element method for simulating rock sample penetration damage in karst areas as described in any one of claims 1-5.

7. The karst area rock sample penetration damage simulation system according to claim 6, characterized in that: include: The constitutive model construction module is used to construct a fracture damage constitutive model that reflects the location of rock discontinuities in karst areas; The first calculation module, based on the discrete element method and the fracture damage constitutive model reflecting the location of rock discontinuities in karst areas, takes input parameters and conditions to calculate the force information of particles. The second calculation module is used to calculate the particle velocity and displacement based on the particle force information and obtain the calculation results. The result output module is used to output the calculation results.

Images