Three-dimensional continuous-discontinuous coupling numerical simulation method for shield TBM excavation jam

Abstract

The application discloses a kind of three-dimensional continuous-discontinuous coupling numerical simulation methods of shield type TBM excavation jamming machine, and relates to the field of rock mechanics numerical method.This method includes the following steps: ①establishing research area geological model, TBM shield model and supporting structure model;②input parameter calibration;③geostress balance;④TBM excavation simulation;⑤TBM shield jamming determination.The application has high calculation precision, and the simulation effect is more consistent with engineering practice, is suitable for TBM excavation surrounding rock large deformation, shield-surrounding rock interaction and shield jamming simulation, and can provide decision-making guidance for TBM tunnel excavation safety construction, jamming early warning and accident handling.

Description

Technical Field

[0001] This invention relates to the field of numerical methods in rock mechanics, and more particularly to a three-dimensional continuous-discontinuous coupled numerical simulation method for a shield-type TBM excavator. Background Technology

[0002] Currently, a number of large and extra-large transportation and water conservancy projects have commenced construction, while coal resources have also entered the deep mining stage. Deeply buried long tunnels are often key control points in the construction of these large and extra-large transportation, water conservancy, and mining projects, and their construction safety constrains the safety and progress of major national infrastructure construction. Compared to the traditional drill-and-blast method, the TBM (Tunnel Boring Machine) method has advantages such as economy, efficiency, safety, and high quality, and has become the preferred choice for the construction of deep-buried long tunnels. However, TBMs are prone to large deformation and jamming disasters when excavating in soft strata with high ground stress, accounting for as much as 37% of the engineering geological disasters encountered during TBM construction. The large deformation process caused by compression in deeply buried soft strata is usually accompanied by the gradual failure and movement of the surrounding rock. The surrounding rock blocks, separated by existing discontinuous structural planes and fractures, undergo mutual compression, slippage, separation, and rotation under high ground stress. As the support effect at the tunnel face disappears, the fracture network further expands. Simultaneously, when the internal stress of the blocks exceeds the rock strength, further fracturing occurs, making the mechanical response of the fractured surrounding rock system extremely complex. This complex mechanical behavior is difficult to quantitatively describe through theoretical analysis and laboratory experiments. The fundamental reason for the frequent occurrence of TBM jamming disasters in deeply buried soft strata is the current lack of scientifically reliable methods for analyzing and predicting TBM jamming during excavation in such strata. Therefore, there is an urgent need to develop three-dimensional numerical simulation methods for the spatiotemporal evolution of surrounding rock fracturing during high-stress TBM excavation in soft surrounding rock and for the interaction between the shield and the surrounding rock.

[0003] Most existing TBM excavation simulations are based on continuum mechanics methods, such as the finite element method (FEM) and the finite difference method (FD). However, traditional continuum mechanics methods struggle to simulate rock fracture and the contact compression behavior of fractured blocks, while discontinuous media mechanics methods generally suffer from computational efficiency bottlenecks and difficulties in parameter determination. The continuous-discontinuous coupled numerical simulation method, represented by the Combined Finite-Discrete Element Method (FDEM), uses solid elements to simulate the continuous deformation of rock and zero-thickness joint elements inserted at the boundaries of solid elements to simulate the damage and fracture process of rock. The FDEM method combines the advantages of the finite element method in handling continuous media and the discrete element method in handling discontinuous media, and can simulate the process of material transitioning from continuous to discontinuous, which is more consistent with engineering reality. It has achieved great success in simulating the progressive failure process of rock in two dimensions. However, research on three-dimensional FDEM is still insufficient, and there are no successful cases of TBM tunnel excavation simulation. Reliable theoretical foundations are still lacking in areas such as stress balance, parameter calibration, and the characterization of TBM shields and supports in three-dimensional FDEM.

[0004] Therefore, to address the aforementioned shortcomings, this invention proposes a three-dimensional continuous-discontinuous coupled numerical simulation method for shield-type TBM excavation machines based on three-dimensional FDEM, including a complete set of simulation technologies such as model establishment, parameter calibration, ground stress balance, TBM continuous excavation simulation, and TBM shield machine determination. Summary of the Invention

[0005] The purpose of this invention is to overcome the problems existing in the prior art and provide a three-dimensional continuous-discontinuous coupled numerical simulation method for shield-type TBM excavation trucks.

[0006] The objective of this invention is achieved as follows:

[0007] Specifically, this method includes the following steps:

[0008] ①Establish geological models, TBM shield models, and support structure models for the study area.

[0009] Modeling a three-dimensional continuous-discontinuous coupled numerical model involves two steps:

[0010] a. Establish a tetrahedral element finite element model and generate it using common finite element modeling software;

[0011] b. Insert cohesive elements with no initial thickness between tetrahedral elements, implemented through a program, for damage and fracture simulation;

[0012] The geological model is large enough to offset boundary effects. The geological structure and stratigraphic structure in the model are consistent with reality. It is divided into tetrahedral element meshes and cohesive elements are inserted between the tetrahedral elements.

[0013] The TBM shield model is generated using triangular element meshes. Since the shield material is high-strength steel, only elastic deformation is considered in the simulation, while the interaction force between the surrounding rock and the shield is recorded.

[0014] The support structure model described is consistent with the actual site, including gravel filling grouting and segments, all of which are divided into tetrahedral element meshes, with cohesion elements inserted between the tetrahedral elements;

[0015] ② Input parameter calibration

[0016] The input parameters mainly include elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, internal friction angle, contact stiffness, normal and tangential deformation stiffness of cohesion unit, tensile fracture energy, and shear fracture energy.

[0017] Complete rock samples from different strata in the field study area were collected and subjected to indoor triaxial compression and direct tensile tests to obtain the elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, and internal friction angle of different strata in the geological model.

[0018] Samples with different structural planes were taken from the field study area, and indoor structural plane shear tests were carried out to obtain the cohesion and internal friction angle of different structural planes in the geological model.

[0019] The contact stiffness of rocks in different strata in the geological model was set to 10 times the rock elastic modulus, and the normal and tangential deformation stiffness of the cohesive unit was set to 500 times the rock elastic modulus. Triaxial compression and direct tension numerical simulations were carried out to calibrate the tensile fracture energy and shear fracture energy of rocks in different strata.

[0020] ③Ground stress balance

[0021] The process of achieving geostress equilibrium is as follows:

[0022] a. Conduct on-site in-situ stress tests to obtain in-situ stress values ​​for the study area;

[0023] b. Based on the finite element geological model without the insertion of initial cohesive elements of zero thickness, surface forces opposite to the actual geostress state are applied to the boundary of the finite element geological model. That is: if the geostress at the actual boundary is tensile stress, then a compressive surface force of the same magnitude is applied; if the geostress at the actual boundary is compressive stress, then a tensile surface force of the same magnitude is applied; if there is shear stress at the actual boundary, then a shear force equal in magnitude but opposite in direction to the actual shear stress is applied; if gravity is considered, then gravity is reversed in the finite element geological model.

[0024] c. Apply vertical constraints to the bottom surface of the finite element geological model, and perform static equilibrium analysis using the Neo-Hookean hyperelastic constitutive model, as shown in the following expression:

[0025]

[0026] In the formula, μ and λ are Lamé constants, determined by the elastic modulus and Poisson's ratio; η is viscous damping, determined by density and elastic modulus; I is the unit tensor, D is the deformation velocity tensor, determined by the element node velocity; J is the deformation gradient determinant, determined by the element node coordinates; T is the Cauchy stress tensor, and B is the Green strain tensor; static equilibrium analysis can be solved using common finite element software or programs.

[0027] d. Insert cohesive elements and extract the coordinates of the deformed finite element geological model nodes from the static equilibrium analysis results as the initial coordinates of the nodes in the three-dimensional continuous-discontinuous coupled numerical model, while the original coordinates of the undeformed nodes of the finite element geological model are used as the current coordinates of the nodes in the three-dimensional continuous-discontinuous coupled numerical model.

[0028] e. Apply normal displacement constraints on the boundary surface corresponding to the boundary surface where surface forces are applied in the finite element geological model in the three-dimensional continuous-discontinuous coupled numerical model, and apply vertical displacement constraints on the bottom surface of the three-dimensional continuous-discontinuous coupled numerical model; if gravity is considered, apply actual gravity.

[0029] f. Start the main program of the three-dimensional continuous-discontinuous coupled numerical simulation method and calculate until the geostress equilibrium is reached; as the calculation time step increases, the total kinetic energy of the system decreases rapidly, and the geostress equilibrium condition is that the total kinetic energy (e) of the system is less than a given threshold (e). cr ):

[0030]

[0031] In the formula, m is the total mass of the model, and v cr This is a speed threshold, generally not greater than 0.001 m / s;

[0032] ④ TBM excavation simulation

[0033] Based on the stress balance in step ③, the tunnel is first excavated at a distance three times the tunnel diameter in one go. At the same time, gravel filling grouting and segmentation are activated within the excavation range. The calculation is performed until equilibrium is reached to eliminate the excavation boundary effect. Then, the tunnel face and TBM shield are advanced at a constant excavation rate. The tetrahedral elements in contact with the tunnel face at the front of the shield are deleted. After the TBM has excavated one ring of segments, gravel filling grouting and segmentation are activated behind the shield to simulate the continuous tunneling mode consistent with the actual site.

[0034] ⑤ TBM Shield Delay Judgment

[0035] In step ④, during the TBM excavation simulation, the interaction force between the shield and the surrounding rock is recorded in real time. When the sum of the frictional force components along the tunnel axis is greater than the TBM's escape thrust, it indicates that the TBM shield is stuck. The criterion for shield sticking is:

[0036]

[0037] In the formula, F0 is the TBM's extrication thrust, F′ is the frictional force on the shield surface along the tunnel axis, and f s This refers to the frictional force component along the tunnel axis, representing the interaction force between the shield surface and the surrounding rock at their contact points.

[0038] This invention has the following advantages and positive effects:

[0039] 1. The geostress equilibrium method used in this invention has high calculation accuracy and efficiency. The size of the model after geostress equilibrium is consistent with the size of the study area. It is applicable to different geological models and can be extended to two-dimensional cases.

[0040] 2. This invention proposes a method for representing TBM shields and support structures in a three-dimensional excavation model, realizing integrated simulation of shield-type TBM excavation, support, and jamming detection;

[0041] 3. Compared with existing TBM tunnel excavation simulation, this invention realizes the simulation of the entire process of tunnel surrounding rock excavation, unloading, damage, deformation, and fracture, as well as the interaction between broken blocks and the shield. It can accurately obtain the evolution law of surrounding rock stress field, displacement field, fracture field, and shield-surrounding rock interaction force during the excavation process.

[0042] In summary, this invention has high calculation accuracy and the simulation effect is more consistent with the actual engineering situation. It is suitable for simulating large deformation of surrounding rock, shield-surrounding rock interaction, and shield jamming in TBM excavation. It can provide decision-making guidance for safe construction, jamming early warning, and accident handling in TBM tunnel excavation. Attached Figure Description

[0043] Figure 1 This is a flowchart illustrating the steps of this method;

[0044] Figure 2-1 This is a geological model of the study area in this embodiment of the invention;

[0045] Figure 2-2 This is a TBM shield model in an embodiment of the present invention;

[0046] Figure 2-3 This is a segment model in an embodiment of the present invention;

[0047] Figure 2-4 This is a gravel-filled grouting model in an embodiment of the present invention;

[0048] Figure 3 The stress balance results of the geological model of the study area in this embodiment of the invention;

[0049] Figure 4-1 This is a simulation result of the TBM excavation shield machine in an embodiment of the present invention;

[0050] Figure 4-2 This is the distribution of frictional force on the shield surface in an embodiment of the present invention;

[0051] English to Chinese

[0052] 1. TBM: Tunnel Boring Machine. It is a new type of advanced tunnel construction machinery that uses rotating cutting tools to excavate, break up the surrounding rock, and excavate to form the entire tunnel cross-section. Compared to commonly used methods, TBM integrates drilling, tunneling, and support, using advanced technologies such as electronics, information, telemetry, and remote control to guide and monitor all operations, ensuring the tunneling process is always in optimal condition. Internationally, it has been widely used in water conservancy and hydropower, mining, transportation, municipal, and defense projects.

[0053] 2. FDEM: Combined Finite-Discrete Element Method. Proposed by Professor A. Munjiza, it is currently a hot topic in international numerical methods for rock mechanics. This method combines the advantages of traditional continuum mechanics and discontinuous media mechanics methods, and can simulate the continuous to discontinuous large deformation and failure process of materials, including deformation, damage, fracture, and slippage, collision, and overturning between broken blocks. The simulation results are more consistent with reality. Detailed Implementation

[0054] The following detailed description is provided in conjunction with the embodiments and accompanying drawings.

[0055] like Figure 1 The specific implementation steps of this method are as follows:

[0056] ①Establish geological models, TBM shield models, and support structure models for the study area - 101

[0057] Modeling a three-dimensional continuous-discontinuous coupled numerical model involves two steps:

[0058] a. Establish a tetrahedral element finite element model and generate it using common finite element modeling software;

[0059] b. Insert cohesive elements with no initial thickness between tetrahedral elements, implemented through a program, for damage and fracture simulation;

[0060] The geological model is large enough to offset boundary effects. The geological structure and stratigraphic structure in the model are consistent with reality. It is divided into tetrahedral element meshes and cohesive elements are inserted between the tetrahedral elements.

[0061] The TBM shield model is generated using triangular element meshes. Since the shield material is high-strength steel, only elastic deformation is considered in the simulation, while the interaction force between the surrounding rock and the shield is recorded.

[0062] The support structure model described is consistent with the actual site, including gravel filling grouting and segments, all of which are divided into tetrahedral element meshes, with cohesion elements inserted between the tetrahedral elements;

[0063] In this embodiment, the study area is a TBM tunnel excavation section, such as... Figure 2-1 The tunnel excavation diameter is 4.0m. Due to symmetry, half of the model is used for analysis. Abaqus finite element software is used to establish the finite element analysis model. The model mesh uses three-dimensional tetrahedral elements. To improve simulation accuracy, the elements within a range of 2.5 times the tunnel diameter around the tunnel are refined. Cohesive elements are inserted between the tetrahedral elements using the three-dimensional FDEM cohesive element insertion module. The number of tetrahedral elements and cohesive elements are 1,611,052 and 3,161,650, respectively. To eliminate boundary effects, the model size is set to length × width × height = 30m × 22m × 44m. Vertical displacement constraints are applied to the bottom of the model. The strata in the excavation section of the model are homogeneous and not penetrated by faults.

[0064] The TBM shield is 3.8m in diameter and 5m in length, and is divided using triangular shell units, such as... Figure 2-2 ;

[0065] The outer diameter of the tunnel segment is 3.8m, and the thickness of the tunnel segment is 0.25m. Figure 2-3 ;

[0066] The tunnel arch widening gap is 0.15m, and the bottom widening gap is 0.05m. The widening gap is simulated by grouting with pea gravel filling. Figure 2-4 ;

[0067] ② Input parameter calibration -102

[0068] The input parameters mainly include elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, internal friction angle, contact stiffness, normal and tangential deformation stiffness of cohesion unit, tensile fracture energy, and shear fracture energy.

[0069] Complete rock samples from different strata in the field study area were collected and subjected to indoor triaxial compression and direct tensile tests to obtain the elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, and internal friction angle of different strata in the geological model.

[0070] Samples with different structural planes were taken from the field study area, and indoor structural plane shear tests were carried out to obtain the cohesion and internal friction angle of different structural planes in the geological model.

[0071] The contact stiffness of rocks in different strata in the geological model was set to 10 times the rock elastic modulus, and the normal and tangential deformation stiffness of the cohesive unit was set to 500 times the rock elastic modulus. Triaxial compression and direct tension numerical simulations were carried out to calibrate the tensile fracture energy and shear fracture energy of rocks in different strata.

[0072] In this embodiment, the input parameters calibrated through indoor experiments and numerical simulations are as follows:

[0073] a. The rock's elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, internal friction angle, contact stiffness, cohesion element normal and tangential deformation stiffness, tensile fracture energy, and shear fracture energy are 11.2 GPa, 2380 kg / m³, and 11.2 GPa, 2380 kg / m³, respectively. 3 , 0.32, 1.1MPa, 3.2MPa, 27°, 112GPa, 5600GPa, 5600Gpa, 20J / m 2 h and 70J / m 2 .

[0074] b. The elastic modulus, density, and Poisson's ratio of the shield are 200 GPa, 7500, and 0.2, respectively;

[0075] c. The elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, internal friction angle, contact stiffness, normal and tangential deformation stiffness of the cohesion unit, tensile fracture energy, and shear fracture energy of the segment are 30 GPa and 3000 kg / m, respectively. 3 , 0.25, 3MPa, 15MPa, 35°, 300GPa, 15000GPa, 15000Gpa, 100J / m 2 h and 300J / m 2 .

[0076] ③Geostress Balance - 103

[0077] The process of achieving geostress equilibrium is as follows:

[0078] a. Conduct on-site in-situ stress tests to obtain in-situ stress values ​​for the study area;

[0079] b. Based on the finite element geological model without the insertion of initial cohesive elements of zero thickness, surface forces opposite to the actual geostress state are applied to the boundary of the finite element geological model. That is: if the geostress at the actual boundary is tensile stress, then a compressive surface force of the same magnitude is applied; if the geostress at the actual boundary is compressive stress, then a tensile surface force of the same magnitude is applied; if there is shear stress at the actual boundary, then a shear force equal in magnitude but opposite in direction to the actual shear stress is applied; if gravity is considered, then gravity is reversed in the finite element geological model.

[0080] c. Apply vertical constraints to the bottom surface of the finite element geological model, and perform static equilibrium analysis using the Neo-Hookean hyperelastic constitutive model, as shown in the following expression:

[0081]

[0082] In the formula, μ and λ are Lamé constants, determined by the elastic modulus and Poisson's ratio; η is viscous damping, determined by density and elastic modulus; I is the unit tensor, D is the deformation velocity tensor, determined by the element node velocity; J is the deformation gradient determinant, determined by the element node coordinates; T is the Cauchy stress tensor, and B is the Green strain tensor; static equilibrium analysis can be solved using common finite element software or programs.

[0083] d. Insert cohesive elements and extract the coordinates of the deformed finite element geological model nodes from the static equilibrium analysis results as the initial coordinates of the nodes in the three-dimensional continuous-discontinuous coupled numerical model, while the original coordinates of the undeformed nodes of the finite element geological model are used as the current coordinates of the nodes in the three-dimensional continuous-discontinuous coupled numerical model.

[0084] e. Apply normal displacement constraints on the boundary surface corresponding to the boundary surface where surface forces are applied in the finite element geological model in the three-dimensional continuous-discontinuous coupled numerical model, and apply vertical displacement constraints on the bottom surface of the three-dimensional continuous-discontinuous coupled numerical model; if gravity is considered, apply actual gravity.

[0085] f. Start the main program of the three-dimensional continuous-discontinuous coupled numerical simulation method and calculate until the geostress equilibrium is reached; as the calculation time step increases, the total kinetic energy of the system decreases rapidly, and the geostress equilibrium condition is that the total kinetic energy (e) of the system is less than a given threshold (e). cr ):

[0086]

[0087] In the formula, m is the total mass of the model, and v cr This is the speed threshold, generally not greater than 0.001 m / s.

[0088] In this embodiment, the boundary surface normal of the finite element geological model coincides with the principal stress axis of the actual geostress. There is no shear stress on the boundary surface of the finite element geological model. The geostress (all compressive stresses) in the x, y, and z directions of the actual region are 26.40 MPa, 21.05 MPa (vertical stress on the upper surface of the model), and 35.60 MPa, respectively, taking gravity into account. Therefore, a tensile surface force of 26.40 MPa is applied to the two boundary surfaces in the x direction of the finite element geological model. Due to the vertical displacement constraint at the bottom in the y direction, a tensile surface force of 21.05 MPa is applied only to the top of the finite element geological model. A tensile surface force of 35.60 MPa is applied to the two boundary surfaces in the z direction of the finite element geological model. Gravity is reversed, that is, the direction of gravity is set vertically upward.

[0089] In this embodiment, the static equilibrium analysis is solved using Abaqus finite element software. During the calculation process, the software can automatically determine each parameter in equation (1) based on the input material parameters and element node information.

[0090] In this embodiment, v cr Take 0.001 m / s, e cr =34.56 J, the geostress equilibrium time step is set to 20000, and the calculation time should not exceed 10 minutes. The geostress equilibrium results are as follows: Figure 3 As can be seen from the figure, the magnitude of the equilibrium geostress is consistent with the actual geostress in the region.

[0091] ④TBM Continuous Excavation Simulation - 104

[0092] Based on the stress balance in step ③, the tunnel is first excavated at a distance three times the tunnel diameter in one go. At the same time, gravel filling grouting and segmentation are activated within the excavation range. The calculation is performed until equilibrium is reached to eliminate the excavation boundary effect. Then, the tunnel face and TBM shield are advanced at a constant excavation rate. The tetrahedral elements in contact with the tunnel face at the front of the shield are deleted. After the TBM has excavated one ring of segments, gravel filling grouting and segmentation are activated behind the shield to simulate the continuous tunneling mode consistent with the actual site.

[0093] ⑤ TBM shield jamming judgment -105

[0094] In step ④, during the TBM excavation simulation, the interaction force between the shield and the surrounding rock is recorded in real time. When the sum of the frictional force components along the tunnel axis is greater than the TBM's escape thrust, it indicates that the TBM shield is stuck. The criterion for shield sticking is:

[0095]

[0096] In the formula, F0 is the TBM's extrication thrust, F′ is the frictional force on the shield surface along the tunnel axis, and f sThis refers to the frictional force component along the tunnel axis, representing the interaction force between the shield surface and the surrounding rock at their contact points.

[0097] In this embodiment, the maximum free-traffic thrust of the TBM is 60MN. Figure 4-1 It can be seen that after the TBM has tunneled a certain distance, the surrounding rock undergoes severe deformation and fracturing, with the deformation exceeding the gap between the shield and the surrounding rock, resulting in a strong interaction between the surrounding rock and the shield; the frictional force distribution acting on the shield surface is as follows. Figure 4-2 The calculated total friction force is 72MN, exceeding the maximum freeing thrust that the TBM can provide, thus confirming that the shield is stuck. Therefore, a three-dimensional continuous-discontinuous coupled numerical simulation method for shield-type TBM excavation stuck can achieve integrated simulation of shield-type TBM excavation, support, and stuck detection, and can accurately obtain the evolution process of surrounding rock stress field, displacement field, fracture field, and shield-surrounding rock interaction force during excavation.

Claims

1. A numerical simulation method of three-dimensional continuous-discontinuous coupling for shielded TBM excavation jumbo, characterized in that It comprises the following steps: ①establishing a geological model of the study area, a TBM shield model and a support structure model The three-dimensional continuous-discontinuous coupled numerical model modeling comprises two steps: a. establishing a tetrahedral element finite element model, which is generated by a common finite element modeling software; b. inserting initial non-thickness cohesive elements between the tetrahedral elements, which is realized by a program and used for damage and fracture simulation; The geological model is large enough to offset the boundary effect, the geological structure and stratum structure in the model are consistent with the actual situation, the tetrahedral element grid is adopted for partitioning, and cohesive elements are inserted between the tetrahedral elements; The TBM shield model is generated by a triangular element grid, since the shield material is high-strength steel, only elastic deformation is considered in the simulation, and the interaction force between the surrounding rock and the shield is recorded; The support structure model is consistent with the actual site, including pea gravel filling grouting and pipe segments, both of which are partitioned by a tetrahedral element grid, and cohesive elements are inserted between the tetrahedral elements; ②input parameter calibration The input parameters mainly include elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion, internal friction angle, contact stiffness, normal and tangential deformation stiffness of the cohesive element, tensile fracture energy and shear fracture energy; Complete rock samples of different strata in the study area are taken to carry out indoor triaxial compression experiments and direct tension experiments, so as to obtain the elastic modulus, density, Poisson's ratio, uniaxial tensile strength, cohesion and internal friction angle of different strata in the geological model; Different structural surface samples are taken from the study area to carry out indoor structural surface shear experiments, so as to obtain the cohesion and internal friction angle of different structural surfaces in the geological model; The contact stiffness of the rock in different strata in the geological model is set to 10 times the elastic modulus of the rock, and the normal and tangential deformation stiffness of the cohesive element is 500 times the elastic modulus of the rock; triaxial compression and direct tension numerical simulation are carried out to calibrate the tensile fracture energy and shear fracture energy of different strata; ③geostress balance The geostress balance adopts the following steps: a. carrying out field geostress testing to obtain the geostress value of the study area; b. based on the finite element geological model without inserting initial non-thickness cohesive elements, a surface force opposite to the actual geostress state is applied on the boundary of the finite element geological model, that is, if the geostress on the boundary of the actual area is tensile stress, a compressive surface force of the same size is applied; If the geostress on the boundary of the actual area is compressive stress, a tensile surface force of the same size is applied; If there is shear stress on the boundary of the actual area, a shear force of the same size but opposite direction is applied; if gravity is considered, the gravity is reversed in the finite element geological model; c. vertical constraints are applied on the bottom surface of the finite element geological model, and Neo-Hookean hyperelasticity constitutive is used for static equilibrium analysis, and the expression is as follows: In the formula, μ and λ are Lame constants, which are determined by elastic modulus and Poisson's ratio; η is viscous damping, which is determined by density and elastic modulus; I is unit tensor, D is deformation velocity tensor, which is determined by element node velocity; J is deformation gradient determinant, which is determined by element node coordinate; T is Cauchy stress tensor, B is Green strain tensor; static equilibrium analysis can be solved by common finite element software or program; d. Insert cohesive element, extract the deformed finite element geological model node coordinates in the static equilibrium analysis results as the initial coordinates of the three-dimensional continuous-discontinuous coupled numerical model nodes, and the undeformed node original coordinates of the finite element geological model as the current coordinates of the three-dimensional continuous-discontinuous coupled numerical model nodes; e. Apply normal displacement constraints on the boundaries corresponding to the boundary surfaces on which surface forces are applied in the finite element geological model in the three-dimensional continuous-discontinuous coupled numerical model, and apply vertical displacement constraints on the bottom surface of the three-dimensional continuous-discontinuous coupled numerical model; If gravity is considered, true gravity is applied; f. Start the main program of the three-dimensional continuous-discontinuous coupled numerical simulation method, and calculate to the geostress balance; As the computation time step increases, the total kinetic energy of the system rapidly decreases, and the ground stress equilibrium condition is that the total kinetic energy e of the system is less than a given threshold value e cr : where m is the total mass of the model, v cr is the velocity threshold, not greater than 0.001 m / s; ④TBM excavation simulation On the basis of step ③ geostress balance, first excavate a distance of 3 times the tunnel diameter at one time, and activate the gravel filling grouting and segment in the one-time excavation range at the same time, and calculate to balance to eliminate the excavation boundary effect; Then push the tunnel face and the TBM shield forward at a constant excavation rate, and delete the tetrahedral elements in contact with the tunnel face at the front end of the shield; When the TBM advances a segment distance, activate the gravel filling grouting and segment behind the shield to realize consistent continuous excavation mode simulation with the field; ⑤TBM shield jamming determination In the process of step ④ TBM excavation simulation, the shield- surrounding rock interaction force is recorded in real time, and when the total sum of the friction force component along the axial direction of the tunnel is greater than the TBM escape thrust, it is indicated that TBM shield jamming occurs, and the shield jamming criterion is: In the formula, F0 is the TBM escape thrust, F' is the friction force of the shield surface along the tunnel axis direction, f s is the friction force component of the interaction force between the shield surface and the surrounding rock at the contact position along the tunnel axis direction.

Images