Fdem simulation method and system for hydraulic coupling fracture of tunnel fractured rock mass

Abstract

The application discloses a tunnel fissure rock mass hydraulic coupling fracture FDEM simulation method and system, and relates to the technical field of rock mechanics and rock engineering; the method comprises the following steps: obtaining fissure information of a target region in a target tunnel and mechanical property parameters of a surrounding rock sample; determining a fissure type according to the fissure occurrence and connectivity; constructing an excavation numerical model of the target tunnel by using an FDEM simulation program according to the mechanical property parameters and the fissure type; performing hydraulic coupling fracture simulation based on the excavation numerical model; in the simulation process, calculating normal stress and tangential stress generated by water flow on adjacent triangular units on both sides of each target unit according to the fissure type; the target unit is a quadrilateral unit with a fissure path; converting the normal stress and the tangential stress into corresponding node forces; applying the node forces to the corresponding target unit and updating the fissure opening degree. The application can accurately simulate the instability mechanism of deep tunnel surrounding rock under the combined action of water pressure and ground stress.

Description

Technical Field

[0001] This invention relates to the fields of rock mechanics and rock engineering, and in particular to a method and system for FDEM simulation of hydraulic coupling fracture of fractured rock mass in tunnels. Background Technology

[0002] The surrounding rock of deep underground engineering projects such as transportation tunnels, water conservancy tunnels, and mine roadways generally contains a network of primary fractures. When these fractures connect with karst caves, faults, underground rivers, or surface water bodies, they are often in a state of water filling or even saturation. As the depth of the project increases, the water pressure in the fractures increases significantly. During tunnel excavation, exposure of rock fractures may lead to high-pressure water inrush disasters, seriously threatening construction safety.

[0003] For the surrounding rock of deep, water-bearing tunnels, the stress state is extremely complex. On the one hand, it is subjected to high ground stress, and on the other hand, it is affected by the seepage of high-pressure fracture water. After excavation and unloading, the redistribution of stress in the surrounding rock causes changes in fracture aperture, which in turn alters the seepage path and pressure distribution. The changes in the seepage field, in turn, affect the rock stress field, forming a two-way water-mechanical coupling effect. This coupling process of mutual influence and interaction between the stress field and the seepage field is complex and continuous until a new equilibrium state is reached. Numerical simulation is an important research method to accurately grasp the dynamic changes of the stress field and seepage field of the tunnel surrounding rock throughout the entire space and throughout the excavation process.

[0004] However, neither continuous nor discontinuous numerical methods can effectively simulate the complex two-way coupling process between the stress field and the seepage field. Although the finite element-discrete element coupled method (FDEM) has been proven to be advantageous in simulating the fracturing and swelling process of tunnel surrounding rock and has been applied, it still has significant shortcomings in simulating hydraulic coupling.

[0005] Existing finite element-discrete element (FDEM) coupling methods mainly focus on the fracturing process of surrounding rock under stress fields, and have not yet systematically integrated a seepage field coupling model. In addition, most current FDEM hydraulic coupling simulation studies are limited to simulating the splitting behavior of small-sized rock samples under water pressure, which differs significantly from the engineering scale and mechanism of actual tunnel excavation unloading-progressive fracturing of surrounding rock-seepage dynamic evolution. Summary of the Invention

[0006] To address the problems existing in the prior art, this invention provides a method and system for FDEM simulation of hydraulic coupling fracturing in tunnel fractured rock masses. This method can accurately simulate the instability mechanism of surrounding rock in deep tunnels under the combined action of water pressure and geostress, providing theoretical support for the prevention and control of hydraulic disasters.

[0007] To achieve the above objectives, the present invention provides a tunnel excavation FDEM simulation method, comprising:

[0008] Acquire fracture information in the target area of ​​the target tunnel, as well as the mechanical property parameters of the surrounding rock sample; the fracture information includes: fracture orientation, connectivity, and initial opening.

[0009] The fracture type is determined based on the fracture orientation and connectivity.

[0010] Based on the mechanical performance parameters and fracture type, an excavation numerical model of the target tunnel is constructed using an FDEM simulation program; the excavation numerical model includes triangular elements and quadrilateral elements; the FDEM simulation program embeds calculation formulas for water pressure, water velocity, and water flow rate.

[0011] Hydraulic coupling rupture simulation is performed based on the excavation numerical model; during the simulation, the following operations are performed iteratively until the termination condition is met:

[0012] Based on the fracture type, calculate the normal and tangential stresses generated by the water flow on the adjacent triangular elements on both sides of each target element; the target element is a quadrilateral element with fracture paths.

[0013] The normal stress and tangential stress are converted into corresponding nodal forces;

[0014] The nodal force is applied to the corresponding target element, and the crack opening is updated.

[0015] Optionally, the mechanical property parameters include elastic modulus, Poisson's ratio, compressive strength, tensile strength, cohesion, and internal friction angle; the mechanical property parameters are obtained by conducting uniaxial compression, Brazilian splitting, and triaxial compression tests on the surrounding rock sample.

[0016] Optionally, the formula for calculating the normal contact stiffness of the triangular element is:

[0017]

[0018] In the formula, P n P is the normal stiffness of the triangular element contact. b For the basic stiffness, and P b =0.1448P f P f The penalty value is for quadrilateral elements; assuming triangular elements i and j are adjacent triangular elements on both sides of the target element, Δh i Let Δh be the height of triangular element i in the contact surface direction. j Let Δl be the height of triangular element j in the contact surface direction. i Let be the length of the contact surface between triangular element i and triangular element j.

[0019] Optionally, the fracture types include Class I, Class II, Class III, Class IV, and Class V fractures; Class I fractures contain a living water source with a constant water flow; Class II fractures are located on the tunnel floor and remain statically filled with water after pressure relief; Class III fractures are connected to Class I fractures and can retain water; Class IV fractures are connected to either Class I or Class II fractures and only partially retain water; Class V fractures are connected to Class I fractures and are completely drained after tunnel excavation, thus unable to retain water.

[0020] Optionally, during the simulation, only the hydrostatic pressure is calculated for the Class II, III, and IV fractures; the hydrostatic pressure and hydrodynamic pressure are both 0 for the Class V fracture.

[0021] Optionally, the water pressure calculation formula is: ΔP = ρ w gH; The formula for calculating the water flow velocity is: The formula for calculating water flow rate is: Q = vu0W;

[0022] In the formula, ΔP is the water pressure difference around the tunnel, and ρ w Where is the density of water, g is the acceleration due to gravity, H is the head height of the water outlet at the fissure around the tunnel, v is the flow velocity of the water at the water outlet at the fissure around the tunnel, Q is the flow rate, u0 is the initial opening of the fissure, and W is the width of the fissure.

[0023] Optionally, the formula for calculating the normal stress is:

[0024]

[0025] The formula for calculating the tangential stress is:

[0026]

[0027] In the formula, σ ni H is the normal stress of crack i. i Let τ be the water head height of fracture i. i Let μ be the tangential stress of crack i. i Let ρ be the dynamic viscosity of water in fracture i. w Let g be the density of water, g be the acceleration due to gravity, and v be the acceleration due to gravity. i Let u be the water flow velocity in fissure i. i Let be the crack opening of crack i.

[0028] Optionally, the formula for calculating the nodal force is as follows:

[0029]

[0030] In the formula, f xi0 f xi1 f xi2 and f xi3The numbers represent the nodal forces in the x-direction at nodes i0, i1, i2, and i3, respectively; f yi0 f yi1 f yi2 and f yi3 The nodes i0, i1, i2 and i3 are represented in the y direction, respectively. i0, i1, i2 and i3 are four nodes of a quadrilateral element arranged counterclockwise. Δx is the coordinate difference of the quadrilateral element nodes in the x direction and Δy is the coordinate difference of the quadrilateral element nodes in the y direction.

[0031] This invention also provides a hydraulic coupling fracturing FDEM simulation system for tunnel fractured rock mass, comprising:

[0032] The data acquisition unit acquires fracture information of the target area in the target tunnel, as well as the mechanical property parameters of the surrounding rock sample; the fracture information includes: fracture orientation, connectivity, and initial opening.

[0033] A type determination unit is used to determine the fracture type based on the fracture orientation and connectivity.

[0034] The model building unit is used to construct an excavation numerical model of the target tunnel using an FDEM simulation program based on the mechanical performance parameters and crack type; the excavation numerical model includes triangular elements and quadrilateral elements; the FDEM simulation program embeds calculation formulas for water pressure, water velocity, and water flow rate.

[0035] The simulation unit is used to perform hydraulic coupling rupture simulation based on the excavation numerical model; during the simulation process, the following operations are performed iteratively until the termination condition is met:

[0036] Based on the fracture type, calculate the normal and tangential stresses generated by the water flow on the adjacent triangular elements on both sides of each target element; the target element is a quadrilateral element with fracture paths.

[0037] The normal stress and tangential stress are converted into corresponding nodal forces;

[0038] The nodal force is applied to the corresponding target element, and the crack opening is updated.

[0039] According to specific embodiments provided by the present invention, the present invention discloses the following technical effects:

[0040] The FDEM simulation method for hydraulically coupled fracture of fractured rock mass in tunnels provided by this invention embeds fracture type classification and a hydraulic-coupled calculation iteration mechanism into the FDEM framework. For the first time, it systematically realizes the dynamic bidirectional coupling of stress field and seepage field during the unloading process of deep tunnel excavation in numerical simulation. This invention performs refined classification of fractures based on fracture orientation and connectivity, and applies differentiated calculation formulas for water pressure, velocity, and flow rate for different fracture types. This solves the problem of insufficient seepage field integration in existing FDEM methods and overcomes the limitations of simulating a single stress field.

[0041] By iteratively executing a closed-loop process of "stress calculation - nodal force conversion - fracture opening update," the normal and tangential stresses generated by water flow on the fractured rock mass (adjacent triangular elements) are dynamically converted into nodal loads of the target quadrilateral element, and the fracture opening is updated in real time and fed back to the hydraulic calculation. This mechanism accurately reproduces the interaction between the progressive fracturing of the surrounding rock and the evolution of the seepage path in actual engineering, overcoming the shortcomings of existing methods that can only simulate water pressure fracturing of small-sized samples and are out of touch with engineering scale. This invention achieves high-precision simulation of fracture propagation, dynamic reconstruction of seepage, and evolution of sudden water inrush risk in the fracture network under water-mechanical coupling throughout the entire process of deep tunnel excavation, providing reliable data support for the spatiotemporal prediction and prevention and control strategy formulation of high-pressure sudden water inrush disasters. Attached Figure Description

[0042] The above and other objects, features and advantages of the present invention will become more apparent from the more detailed description of exemplary embodiments of the invention in conjunction with the accompanying drawings, wherein the same reference numerals generally represent the same parts.

[0043] Figure 1 This is a schematic diagram of the method flow for FDEM simulation of hydraulic coupling fracture of fractured rock mass in tunnels, as shown in an embodiment of the present invention.

[0044] Figure 2 This is a schematic diagram illustrating the types of cracks in an embodiment of the present invention;

[0045] Figure 3 This is a schematic diagram of the FDEM triangular element normal contact stiffness calculation model shown in an embodiment of the present invention;

[0046] Figure 4 This is a schematic diagram of the normal stress of a quadrilateral element shown in an embodiment of the present invention;

[0047] Figure 5 This is a schematic diagram of the tangential stress of a quadrilateral element shown in an embodiment of the present invention;

[0048] Figure 6 This is a schematic diagram illustrating the arrangement of quadrilateral unit nodes according to an embodiment of the present invention;

[0049] Figure 7 This is a schematic diagram of the crack opening degree shown in an embodiment of the present invention;

[0050] Figure 8 This is a schematic diagram illustrating the failure mode of fractured rock mass in a tunnel under hydraulic coupling, as shown in an embodiment of the present invention.

[0051] Figure 9 This is a schematic diagram of the module structure of the FDEM simulation system for hydraulic coupling fracture of tunnel fractured rock mass, as shown in an embodiment of the present invention. Detailed Implementation

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

[0053] See Figure 1 , Figure 1 This is a schematic diagram of the FDEM simulation method for hydraulically coupled fracturing of fractured rock mass in tunnels. The FDEM simulation method for hydraulically coupled fracturing of fractured rock mass in tunnels includes the following steps:

[0054] S101: Obtain fracture information of the target area in the target tunnel, as well as the mechanical property parameters of the surrounding rock sample.

[0055] The fracture information includes fracture attitude, connectivity, and initial opening. Fracture attitude refers to the spatial orientation and geometric characteristics of the fracture, such as strike and dip angle, used to describe the spatial distribution of fractures. Connectivity refers to the hydraulic connection or linkage between fractures, such as whether a fracture is connected to a water source or other fractures, which directly affects the formation of seepage paths. Initial opening refers to the initial width of the fracture, reflecting its physical state when undisturbed. This information can be obtained through on-site engineering surveys, such as geological exploration, borehole sampling, and 3D scanning, ensuring that the data accurately reflects the actual geological conditions of the tunnel surrounding rock. Accurate acquisition of fracture information provides crucial input for subsequent fracture classification and seepage field simulation, avoiding simulation distortion due to data bias.

[0056] The mechanical properties of surrounding rock samples can be obtained through laboratory tests. These parameters include core indices such as elastic modulus, Poisson's ratio, compressive strength, tensile strength, cohesion, and angle of internal friction. Specifically, uniaxial compression tests can be conducted to determine the elastic modulus, Poisson's ratio, and uniaxial compressive strength of the surrounding rock samples; Brazilian splitting tests can be performed to determine tensile strength; and triaxial compression tests can be conducted to determine cohesion and angle of internal friction. As core inputs for constructing numerical models, these mechanical properties directly affect the accuracy of rock mass fracture behavior. Therefore, it is understandable that the testing process must strictly adhere to rock mechanics testing standards, ensuring that samples are taken from the in-situ surrounding rock of the target tunnel area to guarantee the accuracy of the parameters.

[0057] Fracture information describes the rock mass structure from a macro-geological perspective, while mechanical property parameters quantify the rock mass's mechanical response from a micro-material perspective. By combining the two, the hydraulic coupling characteristics of the tunnel surrounding rock can be comprehensively characterized. For example, the initial opening degree and mechanical parameters jointly determine the deformation behavior of fractures under water pressure; parameters such as connectivity and cohesion synergistically affect the bidirectional coupling effect of seepage and stress.

[0058] S102: Determine the fracture type based on fracture orientation and connectivity.

[0059] In applications, fracture type classification can be achieved by combining geological mapping and borehole exploration data, and analyzing the connectivity paths and water supply relationships of the fracture network. The fracture orientation determines the spatial distribution of fractures and the direction of potential hydraulic gradients, while connectivity directly determines the water supply, runoff, and discharge paths in the fracture network. Both factors jointly influence the hydraulic behavior of fractures.

[0060] See Figure 2 , Figure 2 This diagram illustrates the different types of fractures. By comprehensively analyzing the fracture orientation and connectivity of the rock mass in the target area, fractures can be classified into five types with different hydraulic characteristics: Type I, Type II, Type III, Type IV, and Type V fractures. Specifically, Type I fractures refer to fractures with a constant internal water source (such as those connected to underground rivers), whose water flow remains constant during the simulation and serve as the main water supply channel. Type II fractures specifically refer to fractures located at the tunnel floor that maintain a static water-filled state after excavation and depressurization, typically connected to closed water bodies or locally saturated areas. Type III fractures refer to fractures connected to Type I fractures (with an active water source) and effectively retain water, maintaining a certain water pressure and flow rate under excavation disturbance. Type IV fractures refer to fractures that, although connected to Type I or Type II fractures, can only partially retain water, with relatively limited conductivity and storage capacity. Class V fissures refer to fissures that, although connected to Class I fissures, are completely drained of water and unable to retain water after tunnel excavation due to direct exposure or connection to drainage channels.

[0061] By classifying fractures, the dynamic response differences of different fractures in the hydraulic coupling process can be distinguished. Type I fractures are the source of the dominant flow; Types II, III, and IV fractures involve water storage and seepage, but with different water pressure states (Type II is static, Type III maintains pressure, and Type IV is partially depressurized); Type V fractures lose their hydraulic function due to evacuation. The classification results directly determine the method and computational intensity of applying hydraulic loads (static and dynamic pressure) to each type of fracture in subsequent simulations. For example, Types II, III, and IV only require calculation of static pressure, Type V ignores hydraulic loads, while Type I requires complete hydrodynamic calculations.

[0062] S103: Based on the mechanical performance parameters and crack type, construct the excavation numerical model of the target tunnel using the FDEM simulation program.

[0063] In applications, the calculation formulas for water pressure, water velocity, and water flow can be embedded into the FDEM simulation program. Based on this, according to the actual geological profile and engineering dimensions of the target tunnel, a corresponding three-dimensional excavation numerical model is established using FDEM. The surrounding rock elements in the excavation numerical model are also known as mesh elements, including triangular elements and quadrilateral elements. Triangular elements represent continuous rock mass, while quadrilateral elements connecting these triangular elements represent potential fracture interfaces.

[0064] When calibrating the Type II fracture energy corresponding to different cohesion and the Type I fracture energy corresponding to different tensile strengths in the excavation numerical model, the calibration method adopted can be referred to the standard procedure described in the invention patent "A Precise Calibration Method for Fracture Energy Eliminating Loading Rate Effect" (Patent No.: ZL202211201119.6). Viscous damping, joint penalty value, and tangential contact stiffness can be calculated as follows:

[0065]

[0066] In the formula: μ is viscous damping, h is the mesh size, i.e., the side length of the triangular element, E is the elastic modulus, ρ is the rock density, and P is the rock density. f To reduce the penalty value, P t This refers to the tangential contact stiffness.

[0067] See Figure 3 , Figure 3 This is a schematic diagram of the FDEM triangular element normal contact stiffness calculation model. Assuming triangular element i and triangular element j are adjacent triangular elements on either side of the target element, the normal contact stiffness of triangular element i or j can be calculated using the following formula:

[0068]

[0069] In the formula, P n P is the normal stiffness of the triangular element contact.b For the basic stiffness, and P b =0.1448P f P f The penalty value for the quadrilateral element; Δh i Let Δh be the height of triangular element i in the contact surface direction. j Let Δl be the height of triangular element j in the contact surface direction. i Let be the length of the contact surface between triangular element i and triangular element j.

[0070] S104: Perform hydraulic coupling rupture simulation based on excavation numerical model.

[0071] During the simulation, the following operations are performed iteratively until the termination condition is met:

[0072] Based on the fracture type, calculate the normal and tangential stresses generated by the water flow on the adjacent triangular elements on both sides of each target element; the target element is a quadrilateral element with fracture paths.

[0073] Convert normal stress and tangential stress into corresponding nodal forces;

[0074] Apply nodal forces to the corresponding target elements and update the crack opening.

[0075] The termination condition can be set as follows: during the simulation, no new cracks are generated within N consecutive calculation time steps, and the rate of change of the stress (normal stress and tangential stress) generated by the water flow on the triangular element is less than or equal to the first preset value, and the rate of change of the seepage field parameters (water pressure, flow velocity, flow rate) is less than or equal to the second preset value.

[0076] The values ​​of the first preset value, the second preset value, and N can be flexibly set according to the actual situation. For example, the first and second preset values ​​can both be set to 0, and N can be set to 3. It can be understood that if the rate of change of stress is less than the first preset value, it means that the rate of change of both normal stress and tangential stress must be less than the first preset value; the same applies to the seepage field parameters.

[0077] During the simulation, only hydrostatic pressure was calculated for Class II, III, and IV fractures; both hydrostatic and hydrodynamic pressures were 0 for Class V fractures.

[0078] For Class I fractures, the calculation process for water pressure, water velocity, and water flow rate is as follows.

[0079] After the tunnel excavation and unloading, the water pressure around the tunnel drops to 0, and the water pressure difference is:

[0080] ΔP=ρ w gH;

[0081] In the formula, ΔP is the water pressure difference around the tunnel, and ρ wLet ρ be the density of water, g be the acceleration due to gravity, and H be the head height of the water outlet at the fissure around the tunnel.

[0082] Furthermore, based on the water pressure difference, the water flow velocity at the outlet of the fissure around the tunnel can be calculated as follows:

[0083]

[0084] In the formula, v is the water flow velocity at the outlet of the fissure around the tunnel.

[0085] The water flow rate can be calculated from the flow velocity:

[0086] Q = vu0W;

[0087] In the formula, Q is the flow rate, u0 is the initial crack opening, and W is the crack width.

[0088] In the simulation, it is assumed that the flow rate in the same fracture is equal and constant. Thus, the flow velocity at any location in the same fracture can be obtained as follows:

[0089]

[0090] In the formula, v i Let u be the water flow velocity in fissure i. i W represents the crack opening of crack i. i Let be the width of crack i.

[0091] See Figure 4 and Figure 5 The formulas for calculating the normal and tangential stresses of the water flow on the two triangular elements are as follows:

[0092]

[0093] In the formula, σ ni H is the normal stress of crack i. i Let τ be the water head height of fracture i. i Let μ be the tangential stress of crack i. i This is the dynamic viscosity of water.

[0094] Furthermore, the normal and tangential stresses are converted into nodal forces of the quadrilateral elements to calculate the influence of water flow on the fracture opening and velocity. The fracture opening and velocity are then updated again in the next calculation step to achieve the simulation of hydraulic coupling.

[0095] See Figure 6 , Figure 6 This is a schematic diagram showing the node arrangement of a quadrilateral element. The formula for calculating the nodal force is:

[0096]

[0097] In the formula, fxi0 f xi1 f xi2 and f xi3 The numbers represent the nodal forces in the x-direction at nodes i0, i1, i2, and i3, respectively; f yi0 f yi1 f yi2 and f yi3 The nodes i0, i1, i2 and i3 are represented in the y direction, respectively. i0, i1, i2 and i3 are four nodes of a quadrilateral element arranged counterclockwise. Δx is the coordinate difference of the quadrilateral element nodes in the x direction and Δy is the coordinate difference of the quadrilateral element nodes in the y direction.

[0098] The formulas for calculating Δx and Δy are as follows:

[0099]

[0100] like Figure 6 As shown in the formula, x i0 x i1 x i2 and x i3 These represent the coordinates of nodes i0, i1, i2, and i3 in the x-direction, respectively, and y-direction. i0 y i1 y i2 and y i3 These correspond to the coordinates of nodes i0, i1, i2, and i3 in the y-direction, respectively.

[0101] See Figure 7 , Figure 7 This is a schematic diagram of the fracture opening. The total fracture opening is the sum of the new opening caused by water flow after tunnel excavation and the initial opening. The new opening is calculated by the FDEM numerical program, while the initial opening is obtained based on on-site exploration, as shown in the following formula:

[0102] u = u0 + Δu;

[0103] In the formula, u is the total crack opening, u0 is the initial crack opening, and Δu is the new crack opening.

[0104] like Figure 8 As shown, Figure 8 This is a schematic diagram of the failure morphology of tunnel fractured rock mass under hydraulic coupling obtained using the method of the present invention.

[0105] Corresponding to the aforementioned application function implementation method embodiments, the present invention also provides a tunnel fracture rock mass hydraulic coupling FDEM simulation system and corresponding embodiments.

[0106] Please see Figure 9 , Figure 9This is a schematic diagram of the module structure of an FDEM simulation system for hydraulically coupled fracturing of rock mass in tunnels. The system includes:

[0107] Data acquisition unit 91 acquires fracture information of the target area in the target tunnel, as well as the mechanical property parameters of the surrounding rock sample; fracture information includes: fracture orientation, connectivity and initial opening.

[0108] Type determination unit 92 is used to determine the fracture type based on fracture orientation and connectivity;

[0109] Model building unit 93 is used to construct the excavation numerical model of the target tunnel using the FDEM simulation program based on mechanical performance parameters and fracture type; the excavation numerical model includes triangular elements and quadrilateral elements; the FDEM simulation program embeds calculation formulas for water pressure, water velocity, and water flow rate;

[0110] Simulation unit 94 is used to perform hydraulic coupling rupture simulation based on the excavation numerical model; during the simulation, the following operations are performed iteratively until the termination condition is met:

[0111] Based on the fracture type, calculate the normal and tangential stresses generated by the water flow on the adjacent triangular elements on both sides of each target element; the target element is a quadrilateral element with fracture paths.

[0112] Convert normal stress and tangential stress into corresponding nodal forces;

[0113] Apply nodal forces to the corresponding target elements and update the crack opening.

[0114] Regarding the system in the above embodiments, the specific manner in which each unit module performs operations has been described in detail in the embodiments related to the method, and will not be elaborated further here.

[0115] The various embodiments of the present invention have been described above. These descriptions are exemplary and not exhaustive, nor are they limited to the disclosed embodiments. Many modifications and variations will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments. The terminology used herein is chosen to best explain the principles, practical application, or improvement of the technology in the market, or to enable others skilled in the art to understand the embodiments disclosed herein.

Claims

1. A method for FDEM simulation of hydraulically coupled fracturing in fractured rock mass of tunnels, characterized in that, include: Obtain fracture information in the target area of ​​the target tunnel, as well as the mechanical property parameters of the surrounding rock samples; The fracture information includes: fracture orientation, connectivity, and initial opening. Based on the fracture orientation and connectivity, the fracture type is determined; the fracture types include Class I, Class II, Class III, Class IV, and Class V fractures; Class I fractures contain a living water source with a constant water flow; Class II fractures are located on the tunnel floor and remain statically filled with water after pressure relief; Class III fractures are connected to Class I fractures and can retain water; Class IV fractures are connected to either Class I or Class II fractures and only partially retain water; Class V fractures are connected to Class I fractures and are completely drained after tunnel excavation, thus unable to retain water. Based on the mechanical performance parameters and fracture type, an excavation numerical model of the target tunnel is constructed using an FDEM simulation program; the excavation numerical model includes triangular elements and quadrilateral elements; the FDEM simulation program embeds calculation formulas for water pressure, water velocity, and water flow rate. Hydraulic coupling rupture simulation is performed based on the excavation numerical model; during the simulation, the following operations are performed iteratively until the termination condition is met: Based on the fracture type, calculate the normal and tangential stresses generated by the water flow on the adjacent triangular elements on both sides of each target element; the target element is a quadrilateral element with fracture paths. The normal stress and tangential stress are converted into corresponding nodal forces; The nodal force is applied to the corresponding target element, and the crack opening is updated; During the simulation, only the hydrostatic pressure was calculated for Class II, III, and IV fractures; the hydrostatic and hydrodynamic pressures for Class V fractures were both 0. The formula for calculating the normal stress is: ； The formula for calculating the tangential stress is: ; In the formula, σ ni For cracks i normal stress, H i For cracks i water head height, τ i For cracks i tangential stress, μ i For cracks i The dynamic viscosity of greywater ρ w For the density of water, g It is the acceleration due to gravity. v i For cracks i Water flow velocity, For cracks i The degree of crack opening.

2. The FDEM simulation method for hydraulically coupled fracturing of fractured rock mass in tunnels according to claim 1, characterized in that, The mechanical properties include elastic modulus, Poisson's ratio, compressive strength, tensile strength, cohesion, and internal friction angle; these mechanical properties are obtained by conducting uniaxial compression, Brazilian splitting, and triaxial compression tests on the surrounding rock sample.

3. The FDEM simulation method for hydraulically coupled fracturing of fractured rock mass in tunnels according to claim 1, characterized in that, The formula for calculating the normal contact stiffness of the triangular element is: ； In the formula, P n The normal stiffness of the triangular element contact. P b For the basic stiffness, and P b =0.1448 P f , P f The penalty value is for quadrilateral elements; assuming a triangular element... and triangular unit For the triangular elements adjacent to the target element on both sides, Δ h i Triangular unit The height Δ in the direction of the contact surface h j Triangular unit The height Δ in the direction of the contact surface l i Triangular unit and triangular unit The length of the contact surface between them.

4. The FDEM simulation method for hydraulically coupled fracturing of fractured rock mass in tunnels according to claim 1, characterized in that, The formula for calculating water pressure is: The formula for calculating the water flow velocity is: ; The formula for calculating water flow rate is: ; In the formula, Δ P The water pressure difference around the tunnel, ρ w For the density of water, g It is the acceleration due to gravity. H The water head height at the outlet of the fissure around the tunnel; v The water flow velocity at the outlet of the fissure around the tunnel; Q For water flow rate, The initial opening of the crack. W The width is the crack width.

5. The FDEM simulation method for hydraulically coupled fracturing of fractured rock mass in tunnels according to claim 1, characterized in that, The formula for calculating the nodal force is: ； ； In the formula, f xi0 , f xi1 , f xi2 and f xi3 Representing nodes in sequence i 0、 i 1. i 2 and i 3 in x Nodal forces in the direction; f yi0 , f yi1 , f yi2 and f yi3 Representing nodes in sequence i 0、 i 1. i 2 and i 3 in y Nodal forces in the direction; i 0、 i 1. i 2 and i 3 represents the four nodes of a quadrilateral element arranged counterclockwise, Δ x For quadrilateral element nodes in x The difference in coordinates of direction, Δ y For quadrilateral element nodes in y The difference in coordinates of direction.

6. A hydraulically coupled fracture FDEM simulation system for tunnel fractured rock mass, characterized in that, include: The data acquisition unit acquires fracture information of the target area in the target tunnel, as well as the mechanical property parameters of the surrounding rock sample; The fracture information includes: fracture orientation, connectivity, and initial opening. A type determination unit is used to determine the fracture type based on the fracture orientation and connectivity. The fracture types include Class I, Class II, Class III, Class IV, and Class V fractures. Class I fractures contain a living water source with a constant water flow. Class II fractures are located on the tunnel floor and remain statically filled with water after pressure relief. Class III fractures are connected to Class I fractures and can retain water. Class IV fractures are connected to either Class I or Class II fractures and only partially retain water. Class V fractures are connected to Class I fractures and are completely drained after tunnel excavation, thus unable to retain water. The model building unit is used to construct an excavation numerical model of the target tunnel using an FDEM simulation program based on the mechanical performance parameters and crack type; the excavation numerical model includes triangular elements and quadrilateral elements; the FDEM simulation program embeds calculation formulas for water pressure, water velocity, and water flow rate. The simulation unit is used to perform hydraulic coupling rupture simulation based on the excavation numerical model; during the simulation process, the following operations are performed iteratively until the termination condition is met: Based on the fracture type, calculate the normal and tangential stresses generated by the water flow on the adjacent triangular elements on both sides of each target element; the target element is a quadrilateral element with fracture paths. The normal stress and tangential stress are converted into corresponding nodal forces; The nodal force is applied to the corresponding target element, and the crack opening is updated; During the simulation, only the hydrostatic pressure was calculated for Class II, III, and IV fractures; the hydrostatic and hydrodynamic pressures for Class V fractures were both 0. The formula for calculating the normal stress is: ； The formula for calculating the tangential stress is: ； In the formula, σ ni For cracks i normal stress, H i For cracks i water head height, τ i For cracks i tangential stress, μ i For cracks i The dynamic viscosity of greywater ρ w For the density of water, g It is the acceleration due to gravity. v i For cracks i Water flow velocity, For cracks i The degree of crack opening.

Images