A method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model
By constructing a variable contact stiffness model that is related to shear stiffness and normal stress, the simulation error of the variation law of rock tensile strength with confining pressure was solved, and the accuracy of rock failure simulation and the stability of confining pressure application were improved.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-08-19
- Publication Date
- 2026-04-03
AI Technical Summary
Existing methods for simulating rock failure based on the discrete element method fail to accurately simulate the non-monotonic variation of rock tensile strength with confining pressure, which increases first and then decreases. This is mainly because the Mohr-Coulomb contact model, which is independent of shear stiffness and normal stress, cannot consider the correlation between shear stiffness and normal stress.
A variable contact stiffness model was adopted to establish a linear relationship between shear stiffness and normal stress. By calibrating and adjusting parameters such as cohesion and internal friction angle through micro-parameter calibration, the shear stiffness under confining pressure exhibits a pattern of first increasing and then decreasing tensile strength.
The simulation results of rock tensile strength are more consistent with experimental results, the stability of confining pressure is improved, the micro-parameter calibration method is effective, and the simulation accuracy is improved.
Smart Images

Figure CN120927443B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of rock failure simulation technology, specifically relating to a method for simulating multiaxial tensile failure of rocks based on a variable contact stiffness model. Background Technology
[0002] Tensile failure of rock is a key failure mode in various geotechnical engineering processes such as hydraulic fracturing and tunnel excavation. Due to the gravity of the overlying rock mass, tensile failure usually occurs under confining pressure, meaning that one principal stress direction is tensile stress, while the other two directions are compressive stress. Numerous experimental studies have shown that the strength envelope of rock exhibits a parabolic shape, showing that tensile strength initially increases and then decreases with increasing confining pressure.
[0003] Currently, the Mohr-Coulomb contact model is commonly used in rock failure simulation based on the discrete element method (DEM), which assumes that the shear stiffness is a constant independent of the normal stress. However, previous studies have found that the shear stiffness of rock is positively correlated with the normal stress. Because existing models do not consider this correlation between shear stiffness and normal stress, it is difficult to accurately simulate the non-monotonic change in rock tensile strength with confining pressure, which initially increases and then decreases. Summary of the Invention
[0004] To address the aforementioned shortcomings in existing technologies, the rock multiaxial tensile failure simulation method based on a variable contact stiffness model provided by this invention solves the problem that existing methods cannot simulate the initial increase and subsequent decrease of the tensile strength of rock under confining pressure.
[0005] To achieve the aforementioned objectives, the present invention employs the following technical solution: a method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model, comprising the following steps:
[0006] S1. Construct the geometric model corresponding to the rock sample and discretize it into several Thiessen polygons;
[0007] S2. Set the boundary conditions for the direct tensile test with confining pressure;
[0008] S3. Construct a variable contact stiffness model in which shear stiffness and normal stress are linearly related, and calibrate the micro-parameters.
[0009] S4. Based on the set boundary conditions and calibration parameters, a variable contact stiffness model is used to perform multiaxial tensile failure numerical simulation on rock samples discretized as Thiessen polygons.
[0010] Furthermore, in step S2, setting boundary conditions includes:
[0011] Apply confining pressure at a constant rate until the target value is reached and then maintain it constant.
[0012] The bottom edge is fixed, and a constant velocity is applied to the top edge to simulate quasi-static loading conditions.
[0013] Furthermore, in step S3, the constructed variable contact stiffness model is expressed as:
[0014]
[0015] In the formula, Indicates shear stiffness. Represents the normal stress on the contact surface Initial shear stiffness at time, This represents the enhancement coefficient.
[0016] Further, in step S3, regarding the shear stiffness An upper limit constraint is set. and lower bound constraints .
[0017] Furthermore, in step S3, for the variable contact stiffness model, the contact stiffness of the contact surface is determined by the normal stiffness. and shear stiffness Control, contact strength is based on tensile strength Cohesion and internal friction angle Determined, and the shear strength of the contact surface. .
[0018] Furthermore, in step S3, for the variable contact stiffness model, the calibrated micro parameters include: normal stiffness. Shear stiffness Tensile strength Cohesion internal friction angle and shear stiffness Upper limit constraint and lower bound constraints ;
[0019] The method for calibrating microscopic parameters is as follows:
[0020] pass The initial normal stiffness was calculated, and then the normal stiffness and initial shear stiffness were adjusted to make the simulated Young's modulus and Poisson's ratio consistent with the experimental results; among them, and These represent the bulk modulus and shear modulus of the rock, respectively. This represents the minimum width of adjacent elements along the normal direction of the contact surface. ;
[0021] Adjustment The value of is chosen to ensure that the tensile strength without confining pressure is consistent with the experimental results;
[0022] The reinforcement factor is determined based on the slope of tensile strength versus confining pressure. , to cohesion internal friction angle Set the lower limit constraint to a larger value that will not cause shear failure. Set to 0, upper limit constraint Set to a larger value, and adjust the enhancement coefficient accordingly. This ensures that the slope of the tensile strength as the confining pressure increases is consistent with the experimental results;
[0023] Adjusting upper limit constraints To ensure that the maximum tensile strength matches the experimental results, the lower limit constraint is set. Calibrated as initial shear stiffness 10%;
[0024] During the process of decreasing tensile strength, the cohesion is adjusted according to the slope and intercept of the tensile strength decrease segment. and internal friction angle This ensures that the results match the experimental results.
[0025] Further, step S4 specifically includes:
[0026] According to the experimental scheme, set the corresponding confining pressure value, carry out direct tensile numerical test with confining pressure, monitor axial stress, axial strain and transverse strain in real time, and record the peak stress of axial stress. When the axial stress drops to 70% of the peak stress, the rock sample is judged to be completely destroyed, and the numerical simulation of rock tensile failure is completed.
[0027] Specifically, the axial strain is calculated based on the initial height of the specimen by removing the top edge, and the average value of the virtual transverse extensometer measurement results is used as the transverse strain.
[0028] The beneficial effects of this invention are as follows:
[0029] (1) The present invention can simulate the law that the tensile strength of rock increases first and then decreases under confining pressure by designing a variable contact stiffness model, which is more consistent with the actual experimental results.
[0030] (2) The confining pressure application method provided by the present invention can achieve uniform and stable application of confining pressure, effectively solving the problem of confining pressure instability caused by traditional methods.
[0031] (3) The calibration method for micro parameters in the variable contact stiffness model proposed in this invention can realize the rapid calibration of relevant micro parameters. Attached Figure Description
[0032] Figure 1The flowchart of the rock multiaxial tensile failure simulation method based on the variable contact stiffness model provided by the present invention is shown.
[0033] Figure 2 A schematic diagram of the numerical model and boundary condition settings for the direct tensile test with confining pressure provided by the present invention.
[0034] Figure 3 Provided by the present invention and Relationship diagram.
[0035] Figure 4 The envelope of tensile strength under confining pressure (granite) obtained from experiments and numerical simulations provided for this invention.
[0036] Figure 5 The envelope of tensile strength under confining pressure (sandstone) obtained from experiments and numerical simulations provided for this invention. Detailed Implementation
[0037] The specific embodiments of the present invention are described below to enable those skilled in the art to understand the present invention. However, it should be understood that the present invention is not limited to the scope of the specific embodiments. For those skilled in the art, various changes are obvious as long as they are within the spirit and scope of the present invention as defined and determined by the appended claims. All inventions utilizing the concept of the present invention are protected.
[0038] This invention provides a method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model, such as... Figure 1 As shown, it includes the following steps:
[0039] S1. Construct the geometric model corresponding to the rock sample and discretize it into several Thiessen polygons;
[0040] S2. Set the boundary conditions for the direct tensile test with confining pressure;
[0041] S3. Construct a variable contact stiffness model in which shear stiffness and normal stress are linearly related, and calibrate the micro-parameters.
[0042] S4. Based on the set boundary conditions and calibration parameters, a variable contact stiffness model is used to perform multiaxial tensile failure numerical simulation on rock samples discretized as Thiessen polygons.
[0043] In step S1 of this embodiment of the invention, a corresponding geometric model is established in UDEC (Universal Distinct Element Code) or other discrete element software according to the geometry of the rock sample. For example, for a standard cylindrical sample (100 mm high and 50 mm in diameter), a corresponding two-dimensional rectangular sample is established.
[0044] Based on the constructed geometric model, the entire geometric model is discretized into multiple Thiessen polygons using the Voronoi method. The Voronoi method can naturally account for the shape of rock grains and the resulting heterogeneity, thus more accurately reflecting the interaction between rock grains and the failure mechanism of the rock. The model width direction should contain at least 10 Thiessen polygons; therefore, it is recommended to set the equivalent diameter of the Thiessen polygons to 3 mm.
[0045] In step S2 of this embodiment of the invention, if stress conditions are directly applied to the left and right sides of the specimen in the UDEC software, the confining pressure will not be stable during axial loading. To solve this problem, as follows: Figure 2 As shown, the boundary conditions set in this invention include:
[0046] Apply confining pressure at a constant rate until the target value is reached and then maintain it constant.
[0047] The bottom edge is fixed, and a constant velocity is applied to the top edge to simulate quasi-static loading conditions.
[0048] In this embodiment, the specific method for applying confining pressure at a constant rate until the target value is reached and maintained is as follows:
[0049] The coordinates of the nodes on both sides of the sample are extracted, and a constant speed (0.1 mm / s) is applied to each node to simulate the uniform loading process of the confining pressure. The reaction forces of the nodes on both sides are extracted in real time and summed to calculate the current confining pressure. When the confining pressure value approaches 10% of the target value, the system automatically adjusts the node speed according to the remaining difference. If the confining pressure exceeds the target value, the speed is adjusted to a negative value to appropriately reduce the pressure and maintain the stability of the confining pressure.
[0050] In this embodiment, a constant speed of 0.01 mm / s is applied to the top edge.
[0051] In step S3 of this embodiment of the invention, the constructed variable contact stiffness model is expressed as follows:
[0052]
[0053] In the formula, Indicates shear stiffness. Represents the normal stress on the contact surface Initial shear stiffness at time, This represents the enhancement coefficient.
[0054] For the variable contact stiffness model constructed above, from a practical physical perspective, the shear stiffness... It is impossible to increase or decrease indefinitely; based on this, such as Figure 3 As shown, in this embodiment, the shear stiffness... An upper limit constraint is set. and lower bound constraints .
[0055] In this embodiment, for the variable contact stiffness model constructed above, the contact stiffness of the contact surface is determined by the normal stiffness. and shear stiffness Control, contact strength is based on tensile strength Cohesion and internal friction angle Determined, and the shear strength of the contact surface. .
[0056] In step S3 of this embodiment of the invention, for numerical simulation based on the discrete element method, the microscopic parameters of the contact are crucial to the accuracy and rationality of the simulation results. For the variable contact stiffness model proposed in this invention, the calibrated microscopic parameters include: normal stiffness. Shear stiffness Tensile strength Cohesion internal friction angle and shear stiffness Upper limit constraint and lower bound constraints ;
[0057] The method for calibrating microscopic parameters is as follows:
[0058] Through formula The initial normal stiffness was calculated, and then the initial normal stiffness and initial shear stiffness were adjusted to make the simulated Young's modulus and Poisson's ratio consistent with the experimental results; among which, and These represent the bulk modulus and shear modulus of the rock, respectively. This represents the minimum width of adjacent elements along the normal direction of the contact surface. ;
[0059] Adjustment The value of is chosen to ensure that the tensile strength without confining pressure is consistent with the experimental results;
[0060] When the confining pressure is small, the tensile strength gradually increases with the increase of the confining pressure. The reinforcement factor is determined according to the slope of the tensile strength versus the confining pressure. , to cohesion internal friction angle Set the lower limit constraint to a larger value that will not cause shear failure. Set to 0, upper limit constraint Set to a larger value, and adjust the enhancement coefficient accordingly. This ensures that the slope of the tensile strength as the confining pressure increases is consistent with the experimental results;
[0061] Adjusting upper limit constraints To ensure that the maximum tensile strength matches the experimental results, the lower limit constraint is set. Calibrated as initial shear stiffness 10%;
[0062] When the confining pressure increases to a certain level, the tensile strength begins to decrease, indicating the onset of shear failure. During this decrease in tensile strength, the cohesion is adjusted based on the slope and intercept of the decreasing tensile strength segment. and internal friction angle This ensures that the results match the experimental results.
[0063] Step S4 in this embodiment of the invention is specifically as follows:
[0064] According to the experimental scheme, set the corresponding confining pressure value, carry out direct tensile numerical test with confining pressure, monitor axial stress, axial strain and transverse strain in real time, and record the peak stress of axial stress. When the axial stress drops to 70% of the peak stress, the rock sample is judged to be completely destroyed, and the numerical simulation of rock tensile failure is completed.
[0065] Specifically, the axial strain is calculated based on the initial height of the specimen by removing the top edge, and the average value of the virtual transverse strain gauge measurements is used as the transverse strain; specifically, refer to Figure 2 Three positions at different heights were selected in the middle of the specimen. The relative displacement between the leftmost and rightmost nodes was measured and divided by the initial width of the specimen to obtain the transverse strain at the corresponding position. The average of the three was taken as the overall transverse strain. The axial stress was calculated by extracting the sum of the reaction forces of each node at the bottom edge of the specimen and dividing it by the length of the bottom edge.
[0066] In one specific embodiment of the present invention, Figures 4-5 The experimental results obtained from numerical simulations of rock failure using the variable contact stiffness model of this invention and the traditional Mohr-Coulomb model are presented. As can be seen from the figures, when using the variable contact stiffness model of this invention, the tensile strength exhibits a trend of first increasing and then decreasing, which is highly consistent with the experimental results. In contrast, the results obtained using the Mohr-Coulomb model differ significantly from the experimental results and cannot simulate the phenomenon of tensile strength gradually increasing with increasing confining pressure. Specifically, when the confining pressure of granite is 6 MPa, the tensile strength reaches its maximum value, with the numerical error of the improved model being 1.67%, while the error of the Mohr-Coulomb model is as high as 12.68%. Similarly, when the confining pressure of sandstone is 50 MPa, the error of the improved model is only 0.14%, while the error of the Mohr-Coulomb model is 18.90%. Clearly, the improved model performs better in simulating tensile failure under confining pressure.
[0067] Specific embodiments have been used to illustrate the principles and implementation methods of this invention. The descriptions of the embodiments above are only for the purpose of helping to understand the method and core ideas of this invention. At the same time, for those skilled in the art, there will be changes in the specific implementation methods and application scope based on the ideas of this invention. Therefore, the content of this specification should not be construed as a limitation of this invention.
[0068] Those skilled in the art will recognize that the embodiments described herein are intended to help the reader understand the principles of the invention, and should be understood that the scope of protection of the invention is not limited to such specific statements and embodiments. Those skilled in the art can make various other specific modifications and combinations based on the technical teachings disclosed in this invention without departing from the spirit of the invention, and these modifications and combinations are still within the scope of protection of this invention.
Claims
1. A method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model, characterized in that, Includes the following steps: S1. Construct the geometric model corresponding to the rock sample and discretize it into several Thiessen polygons; S2. Set the boundary conditions for the direct tensile test with confining pressure; S3. Construct a variable contact stiffness model in which shear stiffness and normal stress are linearly related, and calibrate the micro-parameters. S4. Based on the set boundary conditions and calibration parameters, a variable contact stiffness model is used to perform multiaxial tensile failure numerical simulation on rock samples discretized as Thiessen polygons. In step S3, the constructed variable contact stiffness model is expressed as follows: In the formula, Indicates shear stiffness. Represents the normal stress on the contact surface Initial shear stiffness at time, Indicates the enhancement coefficient; For the variable contact stiffness model, the calibrated micro parameters include: normal stiffness. Shear stiffness Tensile strength Cohesion internal friction angle and shear stiffness Upper limit constraint and lower bound constraints ; The method for calibrating microscopic parameters is as follows: pass The initial normal stiffness was calculated, and then the normal stiffness and initial shear stiffness were adjusted to make the simulated Young's modulus and Poisson's ratio consistent with the experimental results; among them, and These represent the bulk modulus and shear modulus of the rock, respectively. This represents the minimum width of adjacent elements along the normal direction of the contact surface. ; Adjustment The value of is chosen to ensure that the tensile strength without confining pressure is consistent with the experimental results; The reinforcement factor is determined based on the slope of tensile strength versus confining pressure. , to cohesion internal friction angle Set the lower limit constraint to a larger value that will not cause shear failure. Set to 0, by adjusting the enhancement coefficient. This ensures that the slope of the tensile strength as the confining pressure increases is consistent with the experimental results; Adjusting upper limit constraints To ensure that the maximum tensile strength matches the experimental results, the lower limit constraint is set. Calibrated as initial shear stiffness 10%; During the process of decreasing tensile strength, the cohesion is adjusted according to the slope and intercept of the tensile strength decrease segment. and internal friction angle This ensures that the results match the experimental results.
2. The method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model according to claim 1, characterized in that, In step S2, setting boundary conditions includes: Apply confining pressure at a constant rate until the target value is reached and then maintain it constant. The bottom edge is fixed, and a constant velocity is applied to the top edge to simulate quasi-static loading conditions.
3. The method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model according to claim 1, characterized in that, In step S3, regarding the shear stiffness An upper limit constraint is set. and lower bound constraints .
4. The method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model according to claim 1, characterized in that, In step S3, for the variable contact stiffness model, the contact stiffness of the contact surface is determined by the normal stiffness. and shear stiffness The strength of the contact surface is controlled, including tensile strength and shear strength, wherein the shear strength is determined based on the tensile strength. Cohesion and internal friction angle Determined, and shear strength .
5. The method for simulating multiaxial tensile failure of rock based on a variable contact stiffness model according to claim 1, characterized in that, Step S4 specifically involves: According to the experimental scheme, set the corresponding confining pressure value, carry out direct tensile numerical test with confining pressure, monitor axial stress, axial strain and transverse strain in real time, and record the peak stress of axial stress. When the axial stress drops to 70% of the peak stress, the rock sample is judged to be completely destroyed, and the numerical simulation of rock tensile failure is completed. Specifically, the axial strain is calculated based on the initial height of the specimen by removing the top edge, and the average value of the virtual transverse extensometer measurement results is used as the transverse strain.
Citation Information
Patent Citations
Deep high-temperature and high-pressure environment rock stretching and tension-compression cyclic mechanics experiment device
CN111307606A
Simulation method and system for water-induced rock strength degradation based on discrete element method
CN115964901A