A tunnel fluid-structure damage coupling simulation method based on entity-structure unit dynamic mapping

Abstract

This invention relates to the field of underground engineering and computational geotechnical mechanics numerical simulation technology, and discloses a method for simulating fluid-solid damage coupling in tunnels based on dynamic mapping of solid-structure units. The method first constructs a three-dimensional spatial discrete model, including three-dimensional solid units such as ordinary surrounding rock, consolidated grouting ring or loosened surrounding rock, and lining, as well as internal reinforced structural units, and sets stepped water pressure boundary conditions. It then constructs a non-common-node mapping dictionary based on dual constraints of spatial angle and distance. After performing fluid-solid alternation solutions, it extracts characteristic parameters to calculate the joint damage evolution of the consolidated grouting ring solid unit or the loosened surrounding rock solid unit and the lining solid unit. For the lining solid unit, it performs permeability characteristic reconstruction based on dual control of damage variables and crack opening; for the consolidated grouting ring solid unit or the loosened surrounding rock solid unit, it performs permeability coefficient updates based on effective stress and damage variables, achieving high-precision prediction of the joint disaster evolution of high-pressure water-filled tunnels.

Description

Technical Field

[0001] This invention relates to the field of underground engineering and computational geotechnical mechanics numerical simulation technology, specifically to a tunnel fluid-solid damage coupling simulation method based on dynamic mapping of solid-structure units. Background Technology

[0002] Underground engineering projects such as deep-buried water diversion tunnels and pumped storage power station waterways often bear extremely high internal water pressure during operation. Under the action of high-pressure water, the lining is prone to tensile cracking, which leads to high-pressure water seeping into the surrounding rock, causing joint damage or even instability of the lining and the surrounding rock. Therefore, conducting accurate numerical simulation of tunnel fluid-structure interaction is of great significance for ensuring structural safety. However, existing methods for simulating fluid-structure interaction in tunnels generally suffer from the following technical shortcomings at the level of refined calculation: 1. The challenge of data interaction and mapping between non-co-node meshes of solids and structural elements in multiphysics coupling. In refined 3D numerical simulations, the lining is typically discretized into 3D solid elements, while the internal steel reinforcement network is simulated using 1D structural elements. Due to the complex geometry and reinforcement arrangement of the tunnel, the mesh nodes of the solid elements and structural elements often cannot completely overlap (non-common nodes). Although existing numerical simulation methods can achieve mechanical displacement coordination between non-common nodes through built-in constraint equations or interpolation functions, efficient data mapping and dynamic addressing mechanisms for heterogeneous meshes are generally lacking in multiphysics coupled evolution analysis. This lack of underlying data channels makes it difficult for the system to accurately and automatically extract the real-time axial force of the tensile reinforcement adjacent to a specific lining solid element during fluid-structure interaction solution. Consequently, in fluid-structure interaction calculations, the system cannot establish a dynamic closed-loop evolution model that extracts the actual reinforcement force, calculates the physical crack opening of the lining solid element, and feeds back and drives the seepage field update. 2. The seepage evolution model of the lining and surrounding rock is isolated and simplistic, lacking a closed loop of multi-field interaction feedback between deep rock mass and support structure; Current conventional tunnel fluid-structure interaction simulations often oversimplify and fragment the handling of permeability evolution. On the one hand, for lining solid elements, a constant amplification factor is often applied directly to the permeability coefficient after material yielding, failing to establish a feedback mechanism for the tension of steel reinforcement, calculation of physical crack opening, and nonlinear jumps in the permeability coefficient. On the other hand, for surrounding rock solid elements, existing methods mostly rely on simple empirical formulas for macroscopic porosity, ignoring the macroscopic damage caused by equivalent plastic strain in the surrounding rock under complex stress paths, as well as the exponential jump behavior of the permeability coefficient jointly dominated by this damage and the effective principal stress. This simplification causes existing numerical methods to sever the strong hydraulic coupling between lining tension cracking and leakage and surrounding rock damage and deterioration, failing to realistically reproduce the nonlinear joint catastrophic process of high-pressure internal water breaking through the cracked lining, rapidly infiltrating and deteriorating the surrounding rock. Summary of the Invention

[0003] To address the shortcomings of existing technologies, this invention provides a tunnel fluid-structure interaction damage simulation method based on dynamic mapping of solid-structure units. This method has advantages such as solving the technical problem of lacking efficient dynamic addressing and data interaction mechanisms for heterogeneous meshes in multiphysics calculations in existing numerical simulations, thus solving the aforementioned technical problems.

[0004] To achieve the above objectives, the present invention provides the following technical solution: a method for simulating fluid-structure interaction damage in tunnels based on dynamic mapping of solid-structure units, comprising the following steps: S1: Establish a three-dimensional fluid-structure interaction numerical model of the tunnel and set the stepped water pressure boundary conditions; S2: Construct a non-common node mapping dictionary for entity-structural units based on dual constraints of spatial angle and distance, calculate the polar coordinate angle and radial distance between the lining entity unit and the steel structure unit, and bind the attributes of the optimally matched lining entity unit and steel structure unit within the set angle tolerance and distance addressing upper limit to obtain the non-common node mapping dictionary for entity-structural units. S3: Perform alternating fluid-solid solution and extract characteristic parameters; after alternating solution of mechanical and seepage fields, extract the equivalent tensile strain of the lining solid element and the equivalent plastic strain of the consolidated grouting ring solid element or the surrounding rock loosening ring solid element, and extract the maximum axial tensile stress of the steel reinforcement structural element corresponding to the lining solid element through the solid-structural element non-common node mapping dictionary. S4: Calculate the joint damage evolution of the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock and the solid unit of the lining; for the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock and the solid unit of the lining, respectively, use the exponential damage evolution equation and the four-stage tensile damage model to calculate the damage variables of the corresponding material domains and dynamically reduce the mechanical properties; at the same time, calculate the dynamic physical crack opening on the tension side of the lining solid unit based on the extracted axial tensile stress of the steel reinforcement. S5: Closed-loop update of permeability characteristics driven by multiple mechanical responses; Calculate the matrix damage permeability coefficient and the crack permeability coefficient considering roughness of the lining solid element according to the damage variable and the physical crack opening state, and update the equivalent macroscopic permeability coefficient of the lining solid element using the logarithmic mixing rule; At the same time, update the nonlinear permeability coefficient of the consolidated grouting ring solid element or the surrounding rock loosening ring solid element based on the damage jump mechanism and effective stress, and assign it back to the fluid calculation domain. S6: Output coupled evolution data and tunnel fluid-structure damage coupled analysis model.

[0005] As a preferred technical solution of the present invention, step S1 includes the following steps: S1.1: A three-dimensional numerical model of the tunnel and surrounding rock structure was established using HyperMesh and the mesh was discretized. Then, it was imported into FLAC3D 7.0 software. The three-dimensional numerical calculation model of the tunnel and surrounding rock structure includes solid units of the loosened surrounding rock zone or solid grouting zone, solid units of the lining, and solid units of the tunnel. S1.2: After applying the initial geostress field and displacement boundary conditions, the tunnel solid element is set to empty, and mechanical equilibrium calculation of excavation unloading is performed. After excavation equilibrium, the displacement field of the whole model and the plastic strain field of the solid element of the consolidated grouting ring or the solid element of the surrounding rock loosening ring are forcibly cleared to determine the zero point reference state. S1.3: Based on the spatial polar coordinate geometric equation, one-dimensional steel reinforcement structural units arranged in a ring are automatically generated inside the lining solid unit at a set interval to simulate the steel reinforcement network, and the corresponding steel reinforcement cross-sectional area and elastic modulus are assigned. A contact surface with normal stiffness, tangential stiffness and tensile strength is established between the lining solid unit and the consolidated grouting ring solid unit or the surrounding rock loosening ring solid unit, and the fluid-structure interaction parameters are configured. S1.4: Set stepped water pressure boundary conditions, establish the groundwater static water level, and set the initial conditions of the tunnel inner wall surface as a drainage boundary with zero pore water pressure. During the water filling simulation phase, a stepped pressurization method is used to increase the target design internal water pressure. Divided into The loading step, in the _th loading step, in the _th In each loading step, the current internal water pressure .

[0006] As a preferred technical solution of the present invention, step S2 specifically includes the following steps: S2.1: Obtain the geometric feature parameters of the unit space and extract the coordinates of the tunnel center axis. ; Traverse the lining solid elements and steel reinforcement structural elements in the discretized model to establish the dynamic mapping relationship between the lining solid elements and the steel reinforcement structural elements; S2.2: Polar coordinate transformation and axial pre-screening. For the lining solid element, calculate its polar coordinate angle relative to the tunnel center. Set the upper limit for axial distance search With angle search tolerance Perform Z-axis pre-screening: determine the absolute axial distance between the reinforced concrete structural unit and the lining solid unit. Does it meet the requirements? If the condition is met, proceed to S2.3; otherwise, discard the structural unit. S2.3: For the first batch of axially pre-screened... Calculate the polar coordinate angle of each reinforced concrete structural element. Calculate the absolute value of the polar coordinate angle difference between solid elements and structural elements. ,when At that time, a geometric correction is performed on the angle difference. If the angle difference satisfies... Then calculate the three-dimensional Euclidean distance between the lining solid element and the reinforced concrete element. Among all structural elements that satisfy the angle constraints, select Minimum and The optimal matching steel reinforcement structural element; S2.4: Bind the global identification code of the best matching steel reinforcement structural unit obtained in S2.3 to the attribute list of the current lining entity unit; traverse all lining entity units and repeat the above steps to finally generate the entity-structural unit non-common node mapping dictionary.

[0007] As a preferred technical solution of the present invention, in step S2.2, the polar coordinate angle of the lining solid element relative to the tunnel center is calculated. The specific expression is as follows: in, Indicates the coordinates of the tunnel's center axis. , The centroid spatial coordinates of the lining solid unit, Represents the arctangent trigonometric function; No. Polar coordinate angles of each reinforced concrete structural unit The specific expression is as follows: in, , For the first The centroid spatial coordinates of a steel reinforcement structural unit.

[0008] As a preferred technical solution of the present invention, S3 includes: The equivalent plastic strain of the solid element of the consolidated grouting ring was calculated using the Von Mises yield criterion. The specific expression is as follows: in, Indicates the principal strain of lateral expansion. Indicates the intermediate principal strain. Indicates the principal strain under axial compression; in, Represents the cumulative shear plastic strain. Indicates tensile plastic strain. Indicates the shear expansion angle. Represents the sine function; Equivalent tensile strain synthesis of lining solid elements: in, This represents the equivalent tensile strain of the lining solid element. This function takes a positive value; it outputs the corresponding value when the value inside the function is greater than 0, and outputs 0 if the value is less than or equal to 0. , and Represents the total principal strain tensor; S3 also includes calling the entity-structural element non-common node mapping dictionary generated in step S2; for any lining entity element currently being calculated, the optimal matching steel reinforcement structural element bound to it in space is directly located through the entity-structural element non-common node mapping dictionary, and the transient axial tensile force of that structural element is extracted. and combined with the cross-sectional area of ​​the reinforcing bars The actual axial tensile stress of the reinforcing steel at the corresponding location of the solid element is calculated. If multiple adjacent reinforcing bars are mapped, the maximum value of the axial tensile stress is taken as the maximum axial tensile stress. .

[0009] As a preferred technical solution of the present invention, step S4 includes the following steps: S4.1: Calculate the joint damage evolution of the consolidated grouting ring solid element and the lining solid element, specifically including: Calculate the damage variables of solid elements in the consolidated grouting ring. : when hour, ; when hour, ; when hour, ; in, The initial damage strain threshold, For the ultimate destructive strain, Damage index decay coefficient Represents equivalent plastic strain; The equivalent tensile strain of the synthesized lining solid element in step S3 is calculated. Substitute the values ​​into the four-stage evolution model to calculate the damage variables of the lining solid element; Elastic phase : ; Linear microcrack stage : ; Macroscopic cracking acceleration stage : ; Complete destruction stage : ; in, The initial damage threshold, The strain at the end of the linear segment, For ultimate destructive strain; This serves as a reference damage value for the transition point. The exponential decay coefficient is... This represents the equivalent tensile strain of the lining solid element; S4.2: Obtain the updated damage variables of the lining solid element Damage variables of solid elements in the consolidation grouting ring Subsequently, the elastic modulus E, cohesion c, and internal friction angle of the lining solid unit and the consolidated grouting ring solid unit were determined. Dynamic reduction is performed, and the specific expression is as follows: in, , , These are the initial elastic modulus, initial cohesion, and initial internal friction angle, respectively. , These represent residual cohesion and residual internal friction angle, respectively. This represents damage variables, including damage variables of lining solid elements. Damage variables of solid elements in the consolidation grouting ring ; S4.3 Calculation of maximum physical crack opening based on actual stress: For the tensile cracking zone of the lining unit, read the maximum axial tensile stress of the reinforcing steel extracted in step S3. ,when When the set minimum threshold is exceeded, the coefficient of non-uniformity of steel strain is calculated sequentially. Maximum crack opening and crack spacing : in, This refers to the standard value of the tensile strength of concrete. The elastic modulus of the steel reinforcement. To achieve an effective reinforcement ratio, The diameter of the reinforcing bar. The coefficients for calculating the stress characteristics of the lining solid element are as follows: The coefficients for calculating the stress characteristics of the lining solid element are as follows: This refers to the surface shape factor of the reinforcing steel. like If no macroscopic crack has occurred, then it is determined that no macroscopic crack has been generated. .

[0010] As a preferred technical solution of the present invention, step S5 includes the following steps: S5.1: Update the equivalent macroscopic permeability coefficient of the lining solid element using the logarithmic mixed law, specifically: based on the maximum physical crack aperture output in step S4. Extracting the local average crack aperture Calculate the matrix damage permeability coefficient The specific expression is as follows: in, The initial permeability coefficient of the lining unit. and This is the empirical coefficient for matrix permeability evolution. Represents the damage variable of the lining solid element; When the maximum physical crack opening When the minimum opening threshold is less than or equal to the threshold value, calculate the equivalent macroscopic permeability coefficient of the lining unit. The expression is as follows: in, Represents a logarithm; When the maximum physical crack opening When the fracture permeability coefficient is greater than the minimum aperture threshold, calculate the fracture permeability coefficient. And based on the matrix damage permeability coefficient Comprehensive calculation of the equivalent macroscopic permeability coefficient of the lining solid unit Specifically, it is expressed as follows: in, The density of water, Let be the dynamic viscosity coefficient of water. The roughness of the crack. Indicates the local average crack aperture; S5.2: Based on the damage jump mechanism and effective stress, update the nonlinear permeability coefficient of the consolidated grouting ring solid element and assign it back to the fluid computation domain, according to the damage variables of the consolidated grouting ring solid element. Calculate the permeability coefficient and jump coefficient : when hour, ; when hour, ; when hour, in, This represents the set maximum jump coefficient, used to extract the maximum effective principal stress of the solid element in the consolidated grouting ring. Combined with the initial permeability coefficient With coupling coefficient Calculate the nonlinear permeability coefficient of the updated consolidated grouting ring solid element: in, Represents the natural constant.

[0011] Compared with existing technologies, this invention provides a method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units, which has the following advantages: 1. This invention constructs a solid-structure unit mapping dictionary based on dual constraints of polar coordinate angle and spatial distance. This method cleverly utilizes the geometric features of tunnel structures and successfully solves the data interaction problem caused by the non-common nodes between three-dimensional solid units and one-dimensional structural units in fluid-structure interaction dynamic calculations without increasing the computational overhead. This enables real-time, accurate, and targeted extraction of the internal structural forces.

[0012] 2. This invention abandons the crude approach of simply amplifying the permeability coefficient in traditional models. For the lining unit, it calculates the physical crack opening based on the actual stress of the structure, and then uses the crack opening and damage variables to drive the nonlinear jump of the matrix-crack dual-mode logarithmic hybrid permeability coefficient. For the consolidated grouting ring or the loosened surrounding rock ring, it establishes a stress-damage coupled permeability mutation model controlled by effective stress and damage variables. This mechanism is more in line with the actual physical disaster evolution law of deep rock engineering, and highly integrates the processes of tunnel step filling, reinforcement stress extraction, damage evolution, permeability jump and interface separation early warning. The system can not only capture the cracking and leakage of the lining, but also simultaneously reproduce the accelerated deterioration process of the consolidated grouting ring or the loosened surrounding rock ring after high-pressure water breakthrough, providing a scientific and quantitative prediction tool for pressure-limited operation and support optimization of large-scale projects such as water diversion tunnels. Attached Figure Description

[0013] Figure 1 This is the overall flowchart of this method; Figure 2 This is a schematic diagram illustrating the principle of non-common node mapping between solid and structural units; Figure 3 The logic diagram for jointly updating the permeability coefficients of the lining solid unit and the consolidated grouting ring solid unit; Figure 4 Damage variable cloud map of solid element of the consolidation grouting ring; Figure 5 For damage variable cloud diagrams of lining solid elements; Figure 6 The graph shows the variation of the maximum damage variable for the solid element of the consolidated grouting ring. Figure 7 This is a graph showing the variation of the maximum damage variable for the lining solid element. Figure 8 This is a diagram showing the maximum stress variation in the reinforcing steel. Figure 9 This is a diagram of relative radial displacement. Detailed Implementation

[0014] 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.

[0015] This embodiment takes a deep-buried water diversion tunnel subjected to high internal water pressure as an example. The tunnel diameter is 6.6m, the lining thickness is 0.6m, the tunnel depth is 75m-110m, the lining uses C30 concrete, the reinforcing steel is arranged inside the lining concrete, and the concrete protective layer thickness is 60mm. The initial reinforcement scheme is 8C32, the thickness of the consolidation grouting ring is taken as 2m, and the hydraulic radius is taken as 100 times the outer radius of the lining. The tunnel center is about 85m from the groundwater level. The surrounding rock is mainly Class II with a small amount of Class III. The tunnel is designed to have an internal water pressure of 1.5MPa. In this embodiment, taking the solid unit of the consolidated grouting ring as an example, the calculation process of the solid unit of the loosened surrounding rock ring is completely consistent with that of the solid unit of the consolidated grouting ring. Please see Figure 1 A method for simulating fluid-structure interaction damage in tunnels based on dynamic mapping of solid-structure units includes the following steps: S1: Construct a three-dimensional fluid-structure interaction numerical model of the tunnel and set the stepped water pressure boundary conditions; S1.1: Establishing a three-dimensional spatial discrete model and assigning initial values: Based on the geological survey and design data of the water diversion tunnel, a three-dimensional numerical calculation model of the tunnel and surrounding rock structure was established using HyperMesh and the mesh was discretized. This model was then imported into FLAC3D 7.0 software. The solid elements of this model specifically include four parts: solid elements of the loosened surrounding rock zone or the consolidated grouting zone, solid elements of the lining, and solid elements of the tunnel. In the initial mechanical constitutive setting, an ideal elastoplastic model (such as the Mohr-Coulomb model) was assigned to the ordinary surrounding rock solid elements. To accurately simulate the failure evolution of the rock mass and lining, a strain-softening constitutive model was assigned to the consolidated grouting zone solid elements and the lining solid elements, along with initial elastic modulus, Poisson's ratio, cohesion, and internal friction angle. S1.2 Simulated excavation unloading and calculation field variable reset: After applying the initial geostress field and displacement boundary conditions, the tunnel solid element is set to empty, and mechanical equilibrium calculations for excavation unloading are performed. In order to eliminate the interference of the stress path during construction on the high-pressure water filling simulation during operation, after excavation equilibrium, the displacement field of the entire model (including nodal displacement and velocity in all directions) and the plastic strain field (such as shear plastic strain and tensile plastic strain) of the solid element of the consolidated grouting ring or the solid element of the surrounding rock loosening ring are forcibly cleared to zero, and the zero-point reference state of the fluid-structure interaction simulation during operation is established. S1.3 Construction of the composite support system and initialization of fluid-solid parameters: Based on the spatial polar coordinate geometric equations, one-dimensional steel reinforcement structural units arranged in a ring at a set interval are automatically generated inside the lining solid unit to simulate the steel reinforcement network, and the corresponding steel reinforcement cross-sectional area and elastic modulus are assigned. At the same time, a contact surface with normal stiffness, tangential stiffness and tensile strength is established between the lining solid unit and the consolidated grouting ring solid unit to simulate the possible separation behavior of the lining solid unit and the consolidated grouting ring solid unit under high pressure water. In terms of fluid-structure interaction parameter configuration, the isotropic fluid calculation mode is enabled, and the initial porosity and initial matrix permeability coefficient are set for the consolidated grouting ring solid unit and the lining solid unit, respectively, and the Biot coefficient reflecting the effective stress principle is set. S1.4, Apply a stepped internal water pressure dynamic boundary: A static groundwater level was established, and the initial conditions of the tunnel inner wall surface were set as a drainage boundary with zero pore water pressure; during the water filling simulation stage, a stepped pressurization method was used to increase the target design internal water pressure ( ) divided into The loading step, in the _th loading step, in the _th In each loading step, the current internal water pressure In each loading step, a value equal to [value missing] is synchronously applied to the inner surface of the tunnel. The pore water pressure boundary (driving the seepage field) and the normal mechanical pushing stress (driving mechanical deformation) are used to realistically reproduce the complex load process of the gradual filling of the water diversion tunnel. Step S2: Construct a non-common node mapping dictionary for entity-structural elements based on dual constraints of spatial angle and distance; Due to the non-common node problem between lining entity elements and reinforced structural elements, this step establishes a dynamic mapping dictionary through a polar coordinate geometric addressing algorithm to achieve efficient interaction of multiphysics data. This includes the following sub-steps: S2.1 Obtain the geometric feature parameters of the unit space: Extracting the coordinates of the tunnel's center axis ; Traverse the lining solid elements in the discretized model to obtain the first Centroid space coordinates of each entity unit Similarly, traverse the internal steel reinforcement structural elements to obtain the first... Centroid space coordinates of each structural unit ; S2.2 Polar coordinate transformation and axial pre-screening: Regarding the first For each solid element, calculate its polar coordinate angle relative to the tunnel center. Set the upper limit for axial distance search. With angle search tolerance To improve the computational efficiency of cross-mesh addressing, a Z-axis pre-screening is first performed: determining the absolute axial distance between structural elements and solid elements. Does it meet the requirements? If the conditions are met, proceed to the next level of angle matching; otherwise, discard the structural unit. S2.3 Circumferential Angle Constraints and Spatial Distance Addressing: For the first axial pre-screening Calculate the polar coordinate angle of each structural unit. ; Calculate the absolute value of the polar coordinate angle difference between solid elements and structural elements. To eliminate the phase abrupt change of the arctangent function when it crosses the coordinate system quadrant (i.e. or (flipping problem), when At that time, geometric correction is performed on the angle difference: Apply the dual constraint judgment condition: if the angle difference satisfies Then, the three-dimensional Euclidean distance between the two is further calculated: Among all structural elements that satisfy the angular constraints, the spatial Euclidean distance is selected. Minimum (and The optimal matching structural unit; S2.4 Dynamic mapping dictionary generation: The global identification code (ID) of the best-matched steel structure unit is bound to the attribute list of the current lining entity unit; the above addressing process is repeated for all lining entity units, and finally a globally unique entity-structure unit non-common node mapping dictionary is generated to provide a data channel for subsequent stress extraction. S3: Perform alternating fluid-solid solution and extract characteristic parameters; after alternating solution of mechanical and seepage fields, extract the equivalent tensile strain of the lining solid element and the equivalent plastic strain of the consolidated grouting ring solid element or the surrounding rock loosening ring solid element, and extract the maximum axial tensile stress of the steel reinforcement structural element corresponding to the lining solid element through the solid-structural element non-common node mapping dictionary. The selection of the equivalent plastic strain for extracting the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock is up to the operator. The process described below is for the solid unit of the consolidated grouting ring. The extraction and calculation process of the equivalent plastic strain for the solid unit of the loosened surrounding rock is completely the same as that for the solid unit of the consolidated grouting ring. S3 specifically includes the following sub-steps: S3.1 Solution for alternating fluid-solid freezing: To ensure the convergence and stability of multiphysics calculations, an alternating solution strategy is adopted: First, the mechanical calculation module is frozen, and the seepage calculation module is activated to perform a specified number of fluid calculations based on the current pore water pressure boundary; then, the seepage calculation module is frozen, and the mechanical calculation module is activated to drive the solid skeleton to perform a specified number of stress redistribution and deformation calculations using the updated pore water pressure field; the above seepage field and stress field are calculated independently and alternately throughout the gradual filling process, and the dynamic coupling of the seepage field and stress field is achieved through cross-module data mapping. S3.2 Equivalent plastic strain extraction of solid elements of consolidated grouting zone or solid elements of loosened surrounding rock: For the solid element of the consolidated grouting ring, obtain its cumulative shear plastic strain. Tensile plastic strain and shear dilatation angle ; Reconstructing principal plastic strain components based on dilatation properties: Calculating the principal strain of lateral expansion. and axial compressive principal strain Intermediate principal strain ; The equivalent plastic strain of the solid element of the consolidated grouting ring was calculated using the Von Mises yield criterion: S3.3 Equivalent tensile strain synthesis of lining solid elements: For the lining solid element, its principal stress tensor, elastic modulus, and Poisson's ratio are obtained, and the elastic principal strain in three directions is calculated. Combined with the plastic strain components of the element, its total principal strain tensor is synthesized. , , ; To accurately characterize the tensile damage state of the lining solid element, the equivalent tensile strain of the lining solid element is calculated using the Mazars equivalent tensile strain formula: ; In the formula, Macaulay brackets This represents a function that takes a positive value, i.e., when... hour ,when hour ; S3.4 Dictionary-based dynamic extraction of true axial force in reinforcing bars Call the entity-structural element non-common node mapping dictionary generated in step S2; for any lining entity element currently being calculated, directly locate the optimal matching steel reinforcement structural element that is spatially bound to it through the mapping dictionary; Extract the transient axial tensile force of this structural unit and combined with the cross-sectional area of ​​the reinforcing bars The actual axial tensile stress of the reinforcing steel at the corresponding location of the solid element is calculated. If the dictionary maps multiple adjacent reinforcing bars, then the maximum value of the axial tensile stress is taken as the representative tensile stress of that local area. ; S4: Calculate the joint damage evolution of the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock ring and the solid unit of the lining; This step is based on the extracted stress and strain parameters, and performs damage calculation and macroscopic mechanical parameter reduction for different material domains, and solves the local physical crack opening; For calculating the joint damage evolution of the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock and the solid unit of the lining, the calculation methods and processes are completely consistent. The following steps take the solid unit of the consolidated grouting ring as an example. S4 specifically includes the following sub-steps: S4.1 Calculation of nonlinear damage variables for solid elements of the consolidated grouting ring and solid elements of the lining: To accurately reflect the degradation characteristics of different materials, two independent damage evolution models were established; the damage variable D adopted an irreversible update strategy, that is, the damage variable in the current step... ,in This is the damage value from the previous calculation step. This represents the damage value at the current step. (1) For solid units of consolidated grouting rings or solid units of loosened surrounding rock: the equivalent plastic strain extracted in step S3 is used to... Substitute the threshold-based exponential damage evolution equation; The calculation method for the solid element of the loosened surrounding rock is completely consistent with that for the solid element of the consolidated grouting ring. In this embodiment, the damage variable of the solid element of the consolidated grouting ring is calculated. ) when hour, ; when hour, ; when hour, ; In the formula, The initial damage strain threshold, Ultimately, destructive strain The damage index decay coefficient; (2) For the lining solid element, the equivalent tensile strain synthesized in step S3 is... Substitute the values ​​into the four-stage evolution model to calculate the damage variables of the lining solid element. ): Elastic phase : ; Linear microcrack stage : ; Macroscopic cracking acceleration stage : ; Complete destruction stage : ; In the formula, Initial damage threshold, For the end strain of the linear segment and For ultimate destructive strain; This serves as a reference damage value for the transition point. The exponential decay coefficient; S4.2 Dynamic reduction of macroscopic mechanical properties: Obtain the updated damage variables of the lining solid element Damage variables of solid elements in the consolidation grouting ring Subsequently, the elastic modulus E, cohesion c, and internal friction angle of the lining solid unit and the consolidated grouting ring solid unit were determined. Dynamic reduction: In the formula, , , These are the initial elastic modulus, initial cohesion, and initial internal friction angle, respectively. , These are the residual cohesion and the residual internal friction angle, respectively. This represents damage variables, which in this embodiment include damage variables of the lining solid elements. Damage variables of solid elements in the consolidation grouting ring Each of these is calculated separately during the calculation process; S4.3 Calculation of the maximum physical crack opening based on actual force-driven forces: For the tensile cracking zone of the lining solid unit, the maximum axial tensile stress of the reinforcing steel extracted by the mapping dictionary in step S3 is introduced. ;when When the set minimum threshold is exceeded, the coefficient of non-uniformity of steel strain is calculated sequentially. Maximum crack opening and crack spacing ; In the formula, This refers to the standard value of the tensile strength of concrete. The elastic modulus of the steel reinforcement. To achieve an effective reinforcement ratio, The diameter of the reinforcing bar; Calculation coefficients for the stress characteristics of the lining solid unit. Calculation coefficients for the stress characteristics of the lining solid unit. This is the surface shape factor of the reinforcing steel; if trial calculation If no macroscopic crack has occurred, then it is determined that no macroscopic crack has been generated. ; S5: Closed-loop update of permeability characteristics driven by multiple mechanical responses; Calculate the matrix damage permeability coefficient and the crack permeability coefficient considering roughness of the lining solid element according to the damage variable and the physical crack opening state, and update the equivalent macroscopic permeability coefficient of the lining solid element using the logarithmic mixing rule; At the same time, update the nonlinear permeability coefficient of the consolidated grouting ring solid element or the surrounding rock loosening ring solid element based on the damage jump mechanism and effective stress, and assign it back to the fluid calculation domain. The calculation process for the solid element of the loosened surrounding rock is exactly the same as that for the solid element of the consolidated grouting ring. The specific choice between the solid element of the consolidated grouting ring and the solid element of the loosened surrounding rock in this scheme shall be determined by the operator. S5 specifically includes the following sub-steps: S5.1 Update of the dual-mode logarithmic hybrid permeability coefficient for lining solid elements: Based on the maximum physical crack opening of the lining unit in step S4 Extracting the local average crack aperture (In this embodiment, we take) To address the seepage characteristics of the lining unit at different stress stages, a dual-mode seepage mechanism of "matrix-crack" is adopted: First, calculate the permeability coefficient due to matrix damage caused by micromaterial degradation. In the formula, The initial permeability coefficient of the lining unit. and This is the empirical coefficient for matrix permeability evolution; Subsequently, a minimum opening threshold (e.g.) was introduced. ) when When the value is less than or equal to this threshold, it is determined that no macroscopic permeable cracks have been generated in the lining unit, and crack flow is not considered. The logarithm of its permeability coefficient is used. ; when When the value exceeds this threshold, macroscopic open cracks are determined to have occurred in the lining unit; a modified cubic law considering the surface roughness of the crack is introduced to calculate the crack permeability coefficient. : In the formula, The density of water, Let be the dynamic viscosity coefficient of water. The roughness of the crack; Based on this, the logarithmic mixture rule was used to integrate the permeability contributions of the matrix and cracks, and the equivalent macroscopic permeability coefficient of the lining solid element was calculated. : S5.2 Update the nonlinear permeability coefficient of the consolidated grouting zone solid element or the surrounding rock loosening zone solid element based on the damage jump mechanism and effective stress, and assign it back to the fluid computation domain: For solid units of consolidated grouting rings, the evolution of their permeability characteristics is controlled by the combined effects of damage-driven microcrack propagation and effective stress-driven microcrack closure.

[0016] First, based on the damage variables of the solid element of the consolidated grouting ring. Calculate the permeability coefficient and jump coefficient : when hour, ; when hour, ; when hour, ; In the formula, The maximum jump coefficient is set (e.g., 1000 in this embodiment). Subsequently, the maximum effective principal stress of the solid element of the consolidated grouting ring was extracted. Combined with the initial permeability coefficient With coupling coefficient Calculate the nonlinear permeability coefficient of the updated consolidated grouting ring solid element. : S5.3 Cross-physics parameter feedback and evolution iteration: The equivalent macroscopic permeability coefficient of the lining solid element and the nonlinear permeability coefficient of the consolidated grouting ring solid element obtained from the calculation update are assigned to the corresponding fluid calculation unit; this completes the closed-loop iteration of "mechanical stress extraction → damage crack calculation permeability attribute mutation". Subsequently, based on the updated fluid and mechanical parameter field boundary, the system will automatically advance to the next level of internal water pressure loading step, realizing progressive multi-field coupled calculation; S6: Output coupled evolution data; execute iteratively until the target water pressure is reached and perform long-term coupled calculations to extract parameters such as the state of the tunnel contact surface, the steel reinforcement stress at key locations, the maximum steel reinforcement stress, the permeability coefficient of the lining solid element and the consolidated grouting ring solid element, and the damage variable cloud map. In this embodiment, after each stage of internal water pressure loading and the final long-term fluid-structure interaction calculation, the system automatically extracts and processes macro- and micro-evolution data to achieve a comprehensive assessment of the safety of the tunnel lining and consolidation grouting ring system. This includes the following sub-steps: S6.1 Dynamic monitoring of multidimensional structural response parameters: Macroscopic radial displacement data of key parts of the tunnel (including the arch crown, left and right arch waists, and arch bottom) are extracted; at the same time, the average axial tensile stress and peak tensile stress of the steel structural units corresponding to the above key parts are extracted simultaneously through the entity-structural unit mapping dictionary. S6.2 Determination of the contact surface state of the solid unit of the lining solid unit - the solid unit of the consolidated grouting ring: To assess the risk of structural separation caused by the outward pushing of high-pressure water, the radial displacement of the outer edge solid element of the lining and the inner edge solid element of the consolidated grouting ring were extracted along the circumferential direction at a designated axial monitoring section of the tunnel; the relative radial displacement difference between the two was then calculated. And a contact tolerance threshold is introduced (set to 0.5mm in this embodiment): when When the local lining unit and the solidified grouting ring unit are in contact and under stress, it is determined that they are in a state of "contact" and coordinated stress. when When the area is in a detached state, the angle range of the detached zone is output to provide an early warning of possible splitting damage caused by high-pressure water. S6.3 Statistical analysis of spatial evolution characteristics of damage and seepage field: By traversing the entire computational domain, the system extracts the damage variables of the lining solid elements. Damage variables of solid elements of the consolidated grouting ring The calculation results include: maximum, minimum, and average values ​​of the lining solid element; permeability coefficient K (extreme and average values); steel reinforcement stress in key parts of the lining solid element (arch crown, arch waist, and arch bottom); and the contact state between the lining solid element and the consolidated grouting ring solid element. The results show that the maximum damage values ​​of the lining solid element and the consolidated grouting ring solid element are 0.2732 and 0.002347, respectively; and their maximum permeability coefficients are 1.49038 × 10⁻⁶. -9 and 7.09622×10 -10The maximum steel stress within the lining unit reached 144.364 MPa. After the calculation of each internal water pressure loading level was completed, the system called the graphics rendering module to extract the damage variable D, and set the parallel projection and fixed top view. It automatically exported the damage variable cloud map of the lining unit and the damage variable cloud map of the consolidated grouting ring unit, which included quantitative legends, in order to intuitively map the crack propagation distribution and damage evolution range under complex stress paths.

[0017] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A method for simulating fluid-structure interaction damage in tunnels based on dynamic mapping of solid-structure units, characterized in that: Includes the following steps: S1: Establish a three-dimensional fluid-structure interaction numerical model of the tunnel and set the stepped water pressure boundary conditions; S2: Construct a non-common node mapping dictionary for entity-structural units based on dual constraints of spatial angle and distance, calculate the polar coordinate angle and radial distance between the lining entity unit and the steel structure unit, and bind the attributes of the optimally matched lining entity unit and steel structure unit within the set angle tolerance and distance addressing upper limit to obtain the non-common node mapping dictionary for entity-structural units. S3: Perform alternating fluid-solid solution and extract characteristic parameters; after alternating solution of mechanical and seepage fields, extract the equivalent tensile strain of the lining solid element and the equivalent plastic strain of the consolidated grouting ring solid element or the surrounding rock loosening ring solid element, and extract the maximum axial tensile stress of the steel reinforcement structural element corresponding to the lining solid element through the solid-structural element non-common node mapping dictionary. S4: Calculate the joint damage evolution of the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock and the solid unit of the lining; for the solid unit of the consolidated grouting ring or the solid unit of the loosened surrounding rock and the solid unit of the lining, respectively, use the exponential damage evolution equation and the four-stage tensile damage model to calculate the damage variables of the corresponding material domains and dynamically reduce the mechanical properties; at the same time, calculate the dynamic physical crack opening on the tension side of the lining solid unit based on the extracted axial tensile stress of the steel reinforcement. S5: Closed-loop update of permeation characteristics driven by multiple mechanical responses; Based on the damage variables and physical crack opening status, the matrix damage permeability coefficient and the crack permeability coefficient considering roughness of the lining solid element are calculated. The equivalent macroscopic permeability coefficient of the lining solid element is updated using the logarithmic mixture rule. At the same time, the nonlinear permeability coefficient of the consolidated grouting ring solid element or the surrounding rock loosening ring solid element is updated based on the damage jump mechanism and effective stress, and then assigned back to the fluid calculation domain. S6: Output coupled evolution data and tunnel fluid-structure damage coupled analysis model.

2. The method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units according to claim 1, characterized in that: S1 includes the following steps: S1.1: A three-dimensional numerical model of the tunnel and surrounding rock structure was established using HyperMesh and the mesh was discretized. Then, it was imported into FLAC3D 7.0 software. The three-dimensional numerical calculation model of the tunnel and surrounding rock structure includes solid units of the loosened surrounding rock zone or solid grouting zone, solid units of the lining, and solid units of the tunnel. S1.2: After applying the initial geostress field and displacement boundary conditions, the tunnel solid element is set to empty, and mechanical equilibrium calculation of excavation unloading is performed. After excavation equilibrium, the displacement field of the whole model and the plastic strain field of the solid element of the consolidated grouting ring or the solid element of the surrounding rock loosening ring are forcibly cleared to determine the zero point reference state. S1.3: Based on the spatial polar coordinate geometric equation, one-dimensional steel reinforcement structural units arranged in a ring are automatically generated inside the lining solid unit at a set interval to simulate the steel reinforcement network, and the corresponding steel reinforcement cross-sectional area and elastic modulus are assigned. A contact surface with normal stiffness, tangential stiffness and tensile strength is established between the lining solid unit and the consolidated grouting ring solid unit or the surrounding rock loosening ring solid unit, and the fluid-structure interaction parameters are configured. S1.4: Set stepped water pressure boundary conditions, establish the groundwater static water level, and set the initial conditions of the tunnel inner wall surface as a drainage boundary with zero pore water pressure. During the water filling simulation phase, a stepped pressurization method is used to increase the target design internal water pressure. Divided into The loading step, in the _th loading step, in the _th In each loading step, the current internal water pressure .

3. The method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units according to claim 2, characterized in that: S2 specifically includes the following steps: S2.1: Obtain the geometric feature parameters of the unit space and extract the coordinates of the tunnel center axis. ; Traverse the lining solid elements and steel reinforcement structural elements in the discretized model to establish the dynamic mapping relationship between the lining solid elements and the steel reinforcement structural elements; S2.2: Polar coordinate transformation and axial pre-screening. For the lining solid element, calculate its polar coordinate angle relative to the tunnel center. Set the upper limit for axial distance search With angle search tolerance Perform Z-axis pre-screening: determine the absolute axial distance between the reinforced concrete structural unit and the lining solid unit. Does it meet the requirements? If the condition is met, proceed to S2.3; otherwise, discard the structural unit. S2.3: For the first batch of axially pre-screened... Calculate the polar coordinate angle of each reinforced concrete structural element. Calculate the absolute value of the polar coordinate angle difference between solid elements and structural elements. ,when At that time, a geometric correction is performed on the angle difference. If the angle difference satisfies... Then calculate the three-dimensional Euclidean distance between the lining solid element and the reinforced concrete element. Among all structural elements that satisfy the angle constraints, select Minimum and The optimal matching steel reinforcement structural element; S2.4: Bind the global identification code of the best matching steel reinforcement structural unit obtained in S2.3 to the attribute list of the current lining entity unit; traverse all lining entity units and repeat the above steps to finally generate the entity-structural unit non-common node mapping dictionary.

4. The method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units according to claim 3, characterized in that: In S2.2, the polar coordinate angle of the lining solid element relative to the tunnel center is calculated. The specific expression is as follows: in, Indicates the coordinates of the tunnel's center axis. , The centroid spatial coordinates of the lining solid unit, Represents the arctangent trigonometric function; No. Polar coordinate angles of each reinforced concrete structural unit The specific expression is as follows: in, , For the first The centroid spatial coordinates of a steel reinforcement structural unit.

5. The method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units according to claim 4, characterized in that: S3 includes: The equivalent plastic strain of the solid element of the consolidated grouting ring was calculated using the Von Mises yield criterion. The specific expression is as follows: in, Indicates the principal strain of lateral expansion. Indicates the intermediate principal strain. Indicates the principal strain under axial compression; in, Represents the cumulative shear plastic strain. Indicates tensile plastic strain. Indicates the shear expansion angle. Represents the sine function; Equivalent tensile strain synthesis of lining solid elements: in, This represents the equivalent tensile strain of the lining solid element. This function takes a positive value; it outputs the corresponding value when the value inside the function is greater than 0, and outputs 0 if the value is less than or equal to 0. , and Represents the total principal strain tensor; S3 also includes calling the entity-structural element non-common node mapping dictionary generated in step S2; for any lining entity element currently being calculated, the optimal matching steel reinforcement structural element bound to it in space is directly located through the entity-structural element non-common node mapping dictionary, and the transient axial tensile force of that structural element is extracted. and combined with the cross-sectional area of ​​the reinforcing bars The actual axial tensile stress of the reinforcing steel at the corresponding location of the solid element is calculated. If multiple adjacent reinforcing bars are mapped, the maximum value of the axial tensile stress is taken as the maximum axial tensile stress. .

6. The method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units according to claim 5, characterized in that: S4 includes the following steps: S4.1: Calculate the joint damage evolution of the consolidated grouting ring solid element and the lining solid element, specifically including: Calculate the damage variables of solid elements in the consolidated grouting ring. : when hour, ; when hour, ; when hour, ; in, The initial damage strain threshold, For the ultimate destructive strain, Damage index decay coefficient Represents equivalent plastic strain; The equivalent tensile strain of the synthesized lining solid element in step S3 is calculated. Substitute the values ​​into the four-stage evolution model to calculate the damage variables of the lining solid element; Elastic phase : ; Linear microcrack stage : ; Macroscopic cracking acceleration stage : ; Complete destruction stage : ; in, The initial damage threshold, The strain at the end of the linear segment, For ultimate destructive strain; This serves as a reference damage value for the transition point. The exponential decay coefficient is... This represents the equivalent tensile strain of the lining solid element; S4.2: Obtain the updated damage variables of the lining solid element Damage variables of solid elements in the consolidation grouting ring Subsequently, the elastic modulus E, cohesion c, and internal friction angle of the lining solid unit and the consolidated grouting ring solid unit were determined. Dynamic reduction is performed, and the specific expression is as follows: in, , , These are the initial elastic modulus, initial cohesion, and initial internal friction angle, respectively. , These represent residual cohesion and residual internal friction angle, respectively. This represents damage variables, including damage variables of lining solid elements. Damage variables of solid elements in the consolidation grouting ring ; S4.3 Calculation of maximum physical crack opening based on actual stress: For the tensile cracking zone of the lining unit, read the maximum axial tensile stress of the reinforcing steel extracted in step S3. ,when When the set minimum threshold is exceeded, the coefficient of non-uniformity of steel strain is calculated sequentially. Maximum crack opening and crack spacing : in, This refers to the standard value of the tensile strength of concrete. The elastic modulus of the steel reinforcement. To achieve an effective reinforcement ratio, The diameter of the reinforcing bar. The coefficients for calculating the stress characteristics of the lining solid element are as follows: The coefficients for calculating the stress characteristics of the lining solid element are as follows: This refers to the surface shape factor of the reinforcing steel. like If no macroscopic crack has occurred, then it is determined that no macroscopic crack has been generated. .

7. The method for simulating tunnel fluid-structure interaction damage based on dynamic mapping of solid-structure units according to claim 6, characterized in that: S5 includes the following steps: S5.1: Update the equivalent macroscopic permeability coefficient of the lining solid element using the logarithmic mixed law, specifically: based on the maximum physical crack aperture output in step S4. Extracting the local average crack aperture Calculate the matrix damage permeability coefficient The specific expression is as follows: in, The initial permeability coefficient of the lining unit. and This is the empirical coefficient for matrix permeability evolution. Represents the damage variable of the lining solid element; When the maximum physical crack opening When the minimum opening threshold is less than or equal to the threshold value, calculate the equivalent macroscopic permeability coefficient of the lining unit. The expression is as follows: in, Represents a logarithm; When the maximum physical crack opening When the fracture permeability coefficient is greater than the minimum aperture threshold, calculate the fracture permeability coefficient. And based on the matrix damage permeability coefficient Comprehensive calculation of the equivalent macroscopic permeability coefficient of the lining solid unit Specifically, it is expressed as follows: in, The density of water, Let be the dynamic viscosity coefficient of water. The roughness of the crack. Indicates the local average crack aperture; S5.2: Based on the damage jump mechanism and effective stress, update the nonlinear permeability coefficient of the consolidated grouting ring solid element and assign it back to the fluid computation domain, according to the damage variables of the consolidated grouting ring solid element. Calculate the permeability coefficient and jump coefficient : when hour, ; when hour, ; when hour, in, This represents the set maximum jump coefficient, used to extract the maximum effective principal stress of the solid element in the consolidated grouting ring. Combined with the initial permeability coefficient With coupling coefficient Calculate the nonlinear permeability coefficient of the updated consolidated grouting ring solid element: in, Represents the natural constant.

Images