Abstract

The interval section of the Liuyang River flood control levee project of the Changsha Metro Line 6 is used as the engineering background of this study. A three-dimensional finite element numerical model of a tunnel shield containing complex interfaces is established by using the multifield coupling software COMSOL. The paper studied the deformation of flood control levees under shield tunnel excavation working conditions. The results show that when the left line shield machine is excavated below the Liuyang River flood control levee, the deformation of the flood control levee and the surrounding rock of the tunnel is biased towards the built right line tunnel. And it has an impact on the bridge pile near the tunnel. When the shield of the left line crosses the flood control levee, it can easily cause a large deformation and displacement of the levee above the area between the two shields. To ensure the controlled deformation of the levee, the construction should ensure the spacing between the two tunnel palm faces as far as possible and should be far away from other structures on the levee. When the shield crosses the levee, the deformation that occurs at the base of the levee is significantly higher than the deformations that occur at other locations of the levee. Displacement monitoring and secondary reinforcement at the base of the flood control levee appear to be necessary. The numerical simulation results validate the feasibility of using the multifield coupling software COMSOL to study the construction modeling of shield tunnels through rivers and its advantages. This method provides a practical framework for similar tunnel performance engineering or displacement monitoring in future projects.

1. Introduction

China’s latest fourteenth 5-year plan and the 2035 vision clearly stated the following: speed up the metropolitan rail transit network. This policy is closely related to subway construction. Shield construction is a common subway construction method. There have been many studies on the influence of tunnel shield construction on surface volume loss and the common problems of shield construction [13]. Previous studies have found that the seepage field is one of the important factors affecting the engineering properties of tunnels. The coupling of the stress field and seepage field for special geological conditions has gradually become the focus of research. Biot [4] considered that the effective stress is equal to the total pressure of the upper layer minus the equivalent pore pressure. Based on the experimental study of triaxial compression mechanics, a modified effective stress principle was proposed. The greatest value of this theory is to lay a theoretical foundation for the coupling of fluid and solid in porous media. Since then, many scholars have conducted in-depth research on the coupling problem of the stress field and seepage field. Although the coupling of the stress field and seepage field is complex, in the actual shield construction process, the volume loss of the shallow surface is mostly determined based on previous experience and theoretical calculation. Bobet [5] calculated the analytical solution of ground stress and displacement under different boundary conditions based on the elastic plane strain model of shallow buried circular section tunnel lining. Pinto and Whittle [6] proposed the calculation method of tunnel deformation based on the theory of sphere expansion after considering the spatial location of tunnel engineering and calculated the analytical solution of three-dimensional deformation of the stratum caused by tunnel excavation. In the meantime, the research on the influence response of the surrounding structures of the tunnel has not been widely studied, especially in the complex multifield coupling conditions of the geological environment. Vu et al. [7] believed that it is difficult to obtain the parameters affecting the ground displacement and surface volume loss through experience in tunnel construction under complex geological conditions. In order to study the stability of surrounding rock in composite strata under excavation disturbance, Shi et al. [8] carried out the soil mechanics model test. The results showed that the response of surrounding rock in composite strata was significantly different. Ng et al. [9] studied the influence of tunnel construction, confined water mining, soil consolidation deformation, and train cyclic load on tunnel settlement and found that the compression deformation of the stratum caused by confined water mining was the main reason for tunnel settlement.

With the maturity of computer technology and the development of various types of numerical simulation software, the establishment of a numerical model based on finite element analysis to simulate shield tunnel construction under complex conditions has become the mainstream. Finite element analysis can be used to study the parameters of tunnel lining deformation, ground settlement [10], and structural response to the tunnel [11]. Many scholars used FLAC3D, Plaxis, COMSOL, and other numerical simulation software to study the fault zone [12], longitudinal displacement [13], and pressure distribution of the rock mass around the tunnel [14]. Tunnel excavation can induce soil displacement, which affects the displacement of surrounding buildings. At the same time, the weight and stiffness of the building also affect the displacement of the tunnel. Haji et al. [15] proposed a two-stage mixed empirical numerical technique and studied the relationship between tunnel excavation and surrounding buildings through numerical analysis. Lueprasert et al. [16] proposed a tunnel deformation evaluation method, using three-dimensional elastic-plastic numerical analysis to study the influence of the tunnel on adjacent piles. Thongraksa et al. [17] established a numerical model to study the failure initiation behavior of surrounding rock in highly confined pressurized caverns. It was found that the initial position of surrounding rock instability was affected by rock strength and initial stress. According to the theory of a beam standing on an elastic foundation, Zhao et al. [18] established the mechanical model of a double-layer presupport system via numerical simulation and obtained the deformation characteristics of the tunnel under different factors. Sagaseta [19] studied the ground settlement and displacement trend of the tunnel face under different supporting forces via numerical simulation. The study shows that the stratum loss caused by shield tunnel construction is the main influencing factor of tunnel face instability. Salehnia et al. [20] obtained the shear strain region through numerical simulation and, thus, obtained the deformation and failure region of soil caused by tunnel excavation.

In the process of tunnel construction, due to the influence of the high hydraulic gradient, the tunnel lining is prone to seepage accidents, and the outburst of the soil-water slurry occurs. Huang et al. [21] studied the damage mechanism and structural stability of the tunnel lining via numerical simulation. Lee et al. [22] studied the horizontal component of the seepage force of the tunnel face during shield tunnel construction via numerical simulation. The results show that groundwater has a great influence on the horizontal seepage flow of the tunnel face. Chen et al. [23] proposed a numerical simulation method for tunnel-induced slurry fracturing. The numerical model fully considers the coupling between the stress distribution of the rock mass and fluid seepage. In order to study the stability of surrounding rock when the shield tunnel is close to a cave with confined water, Li et al. [24] established a three-dimensional geomechanical model to test the process of water in the surrounding rock when a compound earth pressure balance shield machine is tunneling. Kargar et al. [25] adopted the COMSOL finite element software to simulate the concrete lining and surrounding rock that are considered to be isotropic mean homogenous materials and found the time discretization of the interaction between the tunnel advancing rate and rock support. Zhang et al. [26] used COMSOL multiphysics and ABAQUS software to study potential hazards to a tunnel caused by adjacent reservoirs. The seepage path, seepage field, and stress field of the shallow tunnel section at different construction and operation stages in the study area were studied via numerical simulation to verify the semianalytical solution of groundwater entering the lining tunnel in a semi-infinite aquifer based on the conformal mapping technology. Ying et al. [27] used COMSOL to establish a numerical model and carried out parameter analysis on seabed thickness, permeability coefficient, tunnel buried depth, and lining thickness involved in the construction of subsea tunnels. COMSOL multiphysics simulation results showed that the error caused by seabed thickness could not be ignored in tunnel construction.

Dynamic numerical simulation of shield tunneling through existing buildings is a major difficulty. The complex geological conditions and the influence of the fluid-solid coupling problem increase the difficulty of research. The multifield coupling software COMSOL provides feasibility for the study of such problems. This paper is aimed at predicting the influence of the river-crossing shield tunnel in the Liuyang River section of the Changsha Rail Transit Line 6 on the flood control levees and surrounding structures. Based on fluid mechanics, solid mechanics, and fluid-solid coupling theory, the three-dimensional finite element numerical model of the tunnel shield under a complex geological environment is established using COMSOL. Combined with multifield coupling software, a nonlinear iterative solution method for the rock mass problem under the fluid-solid coupling condition is proposed. The influence of shield tunnel excavation on the stability of the complex rock mass is studied. The numerical simulation of shield tunnel engineering is carried out to study the influence of shield tunnel excavation construction on the stability of the complex rock and soil mass. It solves the problem of displacement prediction of the surrounding rock of the flood control embankment and achieves the purpose of providing early warning for tunnel shield construction.

2. Project Overview

2.1. Engineering Situations

This project is the construction project of the Liuyang River section of the Changsha Rail Transit Line 6. The aerial view of the subway line and the surrounding environment of the tunnel construction site is shown in Figure 1. The subway line goes east along the south side of Renmin Road and crosses the Liuyang River after crossing the flood control levee on the west side. The part of the tunnel crossing the Liuyang River is parallel to the Renmin East Road Liuyang River Bridge and continues to cross the east flood control levee after crossing the Liuyang River.

The shield construction method is adopted in this project. The total length of the left tunnel is 1012.807 m, and the total length of the right tunnel is 1007.285 m. The horizontal spacing between tunnels is about 13.2~15.2 m, and the depth from the surface is about 9.2 m~26.9 m. In actual construction, the construction of the two tunnels is carried out by two shield machines with a horizontal spacing of about 80 m between the palm surfaces of the tunnel.

2.1.1. Engineering Geology and Engineering Hydrology

According to the geological exploration report of geotechnical layers in this project interval, the composition, distribution position, and average thickness of each rock layer from top to bottom of the stratum within this project interval are summarized as shown in Table 1.

The general geological distribution map is drawn according to the geological survey of the interval section of this project, as shown in Figure 2.

The width of the Liuyang River within the construction interval section is about 220 m, the average annual water level elevation is 28.68 m, the average annual water depth is 3.64 m, and the water depth range is 0.90~5.50 m. The water level of the Liuyang River varies from 2 m to 4 m throughout the year. The pore water in the soft soil layer of this project section is connected with the rock aquifer, and no obvious impermeable layer was found during the geological exploration.

Most of the metro shield tunneling in the interval section of this project was carried out in strongly weathered muddy silt stone. A few localized areas of the shield tunnel traversed more complex geological rock formations. Some of these rocks will soften to varying degrees when exposed to water, which will lead to a significant reduction in strength when immersed in water for too long, and will lead to damage from water loss and dry cracking when exposed to air for too long.

3. Process-Oriented Finite Element Simulation of Shield Tunneling

3.1. Model Theory

The project underneath the shield construction section of the Liuyang River is affected by gravity and river seepage force in the process of shield excavation. The whole engineering problem can be regarded as a fluid-solid coupling problem. Based on the fluid-solid coupling theory, the finite element numerical model of the construction process is established by using the multifield coupling software COMSOL.

3.1.1. Hydrodynamic Theory

The groundwater flow module of COMSOL contains multiple relevant models and interfaces based on the law of fluid mechanics. These models can be used to describe various porous media seepage problems such as saturated, unsaturated, free flow, and fracture flow. These models can be used to simulate the seepage of rock and soil layer when the shield tunnel passes through the Liuyang River.

The Darcy law model built into the software describes the seepage of fluid in saturated porous media affected by water head gradient and gravity potential energy. Based on the liquid continuity equation, the Darcy law is optimized, and the seepage partial differential equation suitable for simulating the fluid-solid coupling of the rock and soil is established. The steady-state equation expression can be written as where represents the source-sink term, represents the density of the fluid, represents the permeability, represents the viscosity of the fluid and represents the pressure .

In porous media fluid flow, higher flow rates cause an increase in viscous shear stress. So Brinkman introduced a modification of Darcy’s law by introducing the shear stress and proposed the Brinkman equation [28]. The Brinkman equation built into the software describes the rapid flow of a fluid in a porous medium due to the fast flow rate and the energy dissipation caused by external forces being ignored. The modified Brinkman equation can be expressed as where represents the percolation flow rate at any spatial point and represents the effective dynamic viscosity of the fluid in the porous medium.

Richards [29] proposed the Richards equation based on the law of mass conservation and Darcy’s law under unsaturated conditions. The software’s built-in Richards equation describes the pressure, head gradient, and flow of a fluid in an unsaturated porous medium or a variably saturated porous medium. The equation can be expressed as follows: where represents the soil volumetric water content; is equal to 1 when the porous medium is saturated and 0 when it is unsaturated; represents the elastic water release coefficient; represents time; represents the unsaturated hydraulic conductivity; represents the total head, where the total head can be decomposed into the soil water pressure head and the location head ; and represents the source-sink term.

According to the geological exploration of the actual project, combined with the built-in model of the software to simulate the fluid-solid coupling problem of geotechnical engineering, different forms of seepage equations are used for fluid-solid coupling calculations for geotechnical bodies in different environments. The equations used for the seepage characteristics in each geotechnical layer are shown in Figure 3.

3.1.2. Solid Mechanics Theory

The strength criterion is used to characterize the relationship between the material strength parameters and stress characteristics of the geotechnical body at the ultimate stress state. The Mohr-Coulomb strength criterion is one of the most widely used strength criteria. The equation is expressed as follows: where represents the shear stress on the shear surface, represents the cohesive force, represents the positive stress in the vertical direction, and represents the angle of internal friction.

It is well known that the Mohr-Coulomb criterion is a hexagonal pyramid in the π-plane. It is difficult to deal with the constitutive equation from the perspective of numerical calculation; that is, in some cases, the numerical calculation is not convergent. Meanwhile, the Mohr-Coulomb criterion does not consider the effect of principal stress on the mechanical properties of the geotechnical soil. For this reason, Drucker and Prager proposed the famous Drucker-Prager criterion [30] . This criterion is an improvement from the von Mises yield criterion, and the categories considered by the Mohr-Coulomb criterion are covered. The criterion can be represented as a smooth cone in the π-plane. The criterion introduces the effect of intermediate principal stresses; at the same time, it considers the action of hydrostatic pressure on the geotechnical body. It is widely used in the field of geotechnical engineering research, especially in the field of geotechnical engineering numerical simulation research. The Drucker-Prager criterion can be expressed as where and represent the constant parameters related to material cohesion and internal friction angle, respectively; represents the first invariant; and represents the second invariant of the stress bias.

In geotechnical projects with complex geological and hydrological conditions, it is difficult to determine the advantages and disadvantages of the Mohr-Coulomb criterion and Drucker-Prager criterion. The original purpose of this numerical simulation model is to study whether displacement occurs in the tunnel envelope and the magnitude at which it occurs. In the model-solving setup, linear elasticity can be understood as a special case of plasticity, where most structural problems can be viewed as linear elastic problems under very small forces or very small deformations. Therefore, the model setup can be based on the linear elastic model with the addition of the soil plasticity property setup to analyze the plastic deformation of the model. The numerical model is set to be calculated based on the Drucker-Prager criterion, matching the Mohr-Coulomb criterion. In this case, the internal friction angle and cohesion of the Mohr-Coulomb criterion replace, respectively, the experimental constants and related to the angle of internal friction and bonding of rocks in the Drucker-Prager criterion. Based on the above approach, the convergence requirement of the model calculation results can be satisfied, and a certain calculation accuracy can be achieved.

3.1.3. Coupled Fluid-Solid Theorem

The advantages of COMSOL in solving multifield coupled problems make it possible to perform direct coupled computational analysis of multiple physical fields. The software has a built-in PDE module, which enables fluid-solid coupling through various forms of partial differential equation-solving modules. In this paper, the following two types of partial differential equations are chosen for calculations related to fluid-solid coupling.

For general simulation, the general form of partial differential equations can be used, and the equation is as follows: where represents the mass coefficient term, represents the damping coefficient term, represents the field variable, represents the gradient operator, represents the conserved flux term, and represents the source term.

The general form of the partial differential equation is solved by using the scatter of the flux vector and the conservation laws that control multiple physical fields. In practical applications, flux vectors are usually fluxes of conserved quantities such as momentum, mass, and heat, which are typically associated with dependent variables through material laws.

For the numerical simulation of some special problems, the coefficient partial differential equation can be used, and the equation is as follows: where represents the gradient operator, represents the conserved flux diffusion coefficient term, represents the conserved flux convection coefficient term, represents the conserved flux source term, represents the convection coefficient term, and represents the absorption coefficient term.

Coefficient-based partial differential equations are used to specify and solve some partial differential equations in the form of coefficients, and the software defines such partial differential equations by specifying different coefficients. It is easy to see via the form of expression of the partial differential equation: the general partial differential equation is obtained by integrating the diffusion coefficient term, the conservation flux convection coefficient term, the conservation flux source term, the convection coefficient term, and the absorption coefficient term of the coefficient-type partial differential equation into an open-sourced conservation flux term. This leaves the general form partial differential equation without the constraint of fixed terms, making its coverage area much broader. Considering that this experimental study mainly uses partial differential equations in rock mechanics and seepage mechanics, combined with the needs of its own model, this paper selects general form partial differential equations and coefficient-type partial differential equations for the numerical simulation solution of fluid-solid coupling.

3.2. Finite Element Model
3.2.1. Basic Assumptions

(1)Considering the geotechnical body as a porous medium and that the permeability coefficient is isotropic, the fluid seepage in the geotechnical body conforms to the seepage law set during the modeling theory analysis(2)The actual project uses the earth pressure balance shield method, assuming that the pressure at the tunnel palm surface is homogeneous, without considering the seepage balance time during the shield excavation(3)Considering that the shield tunnel is an impermeable lining, the lining perimeters of the shield tunnel and the palm surface are set as impermeable interfaces(4)Since the shield tunnel alignment is almost perpendicular to the flood control levees on the east and west banks of the Liuyang River, the effect of the deviation of the angle between the shield tunnel and the flood control levees on the stability of the levees is ignored

3.2.2. Establishment of Geometric Model

The shield tunnel project interval ZDK40+620 to ZDK41+120, located below the riverbed of the Liuyang River, was selected for the numerical simulation study. The cross-sectional study area of the model is based on the middle line between the two tunnels, extending 130 m to the left and 70 m to the right. The longitudinal section study area of the model that is based on the horizontal center line of the tunnel cross-section shall prevail, extending 30 m downward, and upward to the actual height of each geotechnical layer and viaduct. The distance between the two tunnels is about 13.2~15.2 m, and the burial depth of the interval tunnel is 9.2 m~25.9 m. The inner diameter of the shield tunnel is 5.5 m, the outer diameter is 6.2 m, the tunnel tube piece is manufactured via high-strength reinforced concrete, the thickness of the tube piece is 0.35 m, and the concrete strength grade of the tube piece is C50. The viaduct bridge piles are C45 reinforced concrete pile foundations with pile lengths of 18.0 to 25.0 m and pile diameters of 2.5 m and 3 m. The geometric model of the final numerical simulation is shown in Figure 4.

3.2.3. Set Material Parameters

The software has built-in most types of materials and their property parameters in various fields of study and can also customize the material parameters according to actual needs. After the material parameters are selected, they are simply assigned to the corresponding geometric space domain of the numerical model under the corresponding physical field, which is the material with the specific properties. The values of the material parameters of each domain in the numerical model are set according to the geological survey data of the actual project as shown in Table 2.

3.2.4. Set Boundary Conditions

When numerical simulation is performed using the finite element method, a finite portion of the infinite space needs to be intercepted to build the computational model. For the different physical fields selected, COMSOL gives the assumed equations for each boundary condition. The model boundary conditions are set in conjunction with the actual construction conditions.

The tunnel is directly below the shield construction, which results in infinity which is negligible, so the fixed constraint without displacement is set at the bottom of the model. The displacement of the bottom of the model in all directions is 0. The effect of the seepage field on the -axis at infinity is negligible, but the influence of the gravitational field is always present. The model is oriented along the -axis in the flow direction of the Liuyang River, which always subjects the model to the coupling effect of seepage and gravity fields along the -axis. In order to ensure the convergence of the model, the YZ plane perpendicular to the -axis is set as the roller support. The XZ plane perpendicular to the -axis is set to the specified displacement with a displacement value of 0 along the -axis direction. The actual project is subject to a gravity field that generates ground stress and sets the gravity properties for the full domain of the model. The acceleration of gravity is set in the direction of the coordinate axis , with the direction vertically downward, and the equation can be written as .

The Liuyang River in the model simulates river seepage by adding additional head pressure to the surface of the geometric domain. The stress boundary is set at the interface of the simulated Liuyang Riverbed in the model, the riverbed is subjected to the force of the Liuyang River water, and the riverbed boundary is set as the boundary load. According to the need of the numerical model seepage field boundary conditions, the standard representation is used to set the constrained boundary conditions as the displacement boundary conditions along the specified direction. Constraint boundary conditions can be expressed as . There are flow-free boundary areas in the model, such as the surface of the flood control levee and the exposed ground surface in the real situation, which are set as flow-free boundaries to ensure the convergence of the model. The boundary conditions for no flow are shown as .

3.2.5. Dividing the Grid

The accuracy of the calculation results of the model built based on FEA is closely related to the delineated mesh. Delineating the grid is a key step in optimizing the model results. Based on the comparison results of model convergence metrics, the convergence of a certain grid result can be judged. When the difference in results becomes small enough, the model is considered to have converged. Finally, three levels of grids are used to divide the whole model, and the grid is encrypted for the tunnel lining as well as the area with large changes in the geological layer of the flood control embankment, according to COMSOL mesh division criteria for finite element calculations. The finite element meshing of the model is shown in Figure 5.

3.2.6. Specify the Solution Steps

The original purpose of the model was to study the magnitude of the displacement occurring in the surrounding rock during the construction of the tunnel shield and the area where the displacement occurs. When analyzing the effect of external loads on the overall stability of the structure, the effect of ground stresses needs to be balanced. Therefore, it is necessary to set up the ground stress of the equilibrium gravity field in the model calculation step first and, after that, set up the model subjected to external loads.

This model uses the addition of prestressing force to balance the ground stress via a step-by-step solution. To simplify the model setup as well as to optimize the model calculation speed, the computational study is completed in two steps. The first step of the computational study is only for the solid mechanics model, and the second step of the computational study is for the solid mechanics and hydrodynamic model after equilibrium over the ground stress.

The stress-strain results calculated in the first step have an impact on the properties of the fluid in the fluid model in the second step. In the first step of the calculation, the geotechnical stress state in the unexcavated state of the shield tunnel is calculated. The distribution of the stress state in the unexcavated geotechnical body is determined. This step of the operation can be understood as balancing the ground stress by loading the previous stresses in the form of a prestressing force into the next calculation study. The coupling between the two can be easily achieved via the second Piola-Kirchhoff formula built into the model.

3.2.7. Modeling of the Excavation Process

In order to simulate the actual process of shield excavation, a numerical model of shield excavation is established based on the geometric model of the tunnel and the relevant parameters of the shield machine and the construction process. As shown in Figure 6, the shield consists of three main parts: the cutout ring, the shield machine shell, and the high-strength concrete lining at the end of the shield. During the excavation of the geotechnical body, the shield machine shell and the high-strength concrete lining form a support for the surrounding rock. The geotechnical body is excavated via the cutout ring of the shield machine.

In the numerical model of shield excavation, the excavated part of the geotechnical body is supported by the shield machine shell or high-strength concrete lining. Fixed displacement is set at the palm face to simulate earth pressure balance excavation construction. The excavation simulation of the geotechnical body to be excavated is performed by the software activation operation of the geotechnical body.

When Shi et al. [31] studied the influence of shield tunnel construction in saturated soil on the spatial and temporal distribution of excess pore water pressure, a three-dimensional finite element model was used to simulate the construction process of shield excavation, tunnel propulsion, and segment installation. When using a continuous shell with no joints, aligned with a longitudinal joint configuration with normal ring joints or staggered longitudinal joint configuration with normal ring joints to simulate shield tunnel lining, the numerical simulation of tunnel deformation and stress has similar numerical simulation results [32]. Considering the convenience and formality of modeling, the aligned longitudinal joint configuration with normal ring joints was used to build the tunnel lining.

The detailed simulation of the shield excavation process is shown in Figure 7. When , the shield machine has been installed to the th shield tunnel segment. At this point, the excavated surrounding rock is supported by the shield machine shell and the subsequent high-strength concrete segment. The geotechnical body to be excavated is at the front of the shield machine shell, and the length of a single excavation is the same as the length of one shield tunnel segment. The activate operation in the material model is used to simulate the excavation process of the shield machine. The activate operation in COMSOL software is similar to the delete operation. The activated geometric domain will lose its carrying capacity. The simplified steps of the shield tunnel simulation are as follows: (1) One ring of the shield machine is advanced to activate one ring of the geotechnical body. (2) The shield shell is moved forward one ring distance. (3) One ring of concrete lining is added behind the shell. After the completion of the th shield excavation, the th tunnel shield excavation is performed. The simulation is carried out in the same way as the previous version. The high-strength concrete lining, the shield machine shell, and the area to be excavated are moved forward by one shield tunnel length. The above steps are repeated so that the excavation simulation of the geotechnical body can be realized.

4. Results and Discussions

4.1. Model Analysis in the Initial State before Shield Excavation

Before the study of simulating shield machine excavation of rock and soil, the stress distribution of the model under the initial state is analyzed. The stress distribution of the model under the action of the gravity field is shown in Figure 8. With the increase of depth, the ground stress is increasing. Due to the difference of mechanical properties of rock and soil layers, the stress distribution increases nonlinearly, especially in the Liuyang River flood control levee body and pier part, as shown in Figure 8(a), which details part of the amplification diagram. There is obvious stress redistribution at the built tunnel. The stress increase in the inner wall of the shield tunnel is most pronounced. The stress at the arch waist of the tunnel also increases significantly. The stresses at the top of the tunnel vault and the elevation arch are significantly reduced, as shown in Figure 8(c).

From the provisions of the solution steps, it is clear that the ground stresses need to be balanced before proceeding to simulate the actual shield excavation. After ground stress equilibrium is performed, the model will be in equilibrium without stress and strain if shield excavation is not performed. In this way, in the simulated shield excavation model, stress and strain changes will only occur in the tunnel lining and the surrounding geotechnical body where displacements are likely to occur, and no changes will occur in other areas.

4.2. Analysis of Model Calculation

The analysis of the model calculation results is based on the example of the west side flood control levee. The analysis of the simulation results for the east side flood control levee is similar and will not be repeated.

According to the actual engineering site survey and engineering experience analysis, the flood control levee is most likely to produce large displacements under three working conditions: the right tunnel crossing the flood control levee, the left tunnel crossing the flood control levee before, and both tunnels have crossed the flood control levee. The numerical simulation results of the longitudinal section, cross-section, and top view of the tunnel under the above three conditions are interpreted and described. The position of the screenshot in the model is shown in Figure 9.

4.2.1. Right Line Tunnel through the Flood Control Levee

When the right line shield crosses the Liuyang River flood control levee at this time, the stability of the flood control levee is mainly influenced by the shield in the right line tunnel. The simulation results are shown in Figure 10.

The longitudinal section of the right tunnel is shown in Figure 10(a). The shield housing and the tunnel-surrounding rock in front of it deformed the most, with a maximum displacement of about 6 mm. Since the Liuyang River flood control levee is mainly composed of miscellaneous fill, plain fill, and silty clay, the strength of the soil is generally not high. Under the influence of shield excavation, the soil appears to be displaced in favor of the front of the shield housing, and the further down the soil is displaced, the greater the displacement. The maximum displacement settlement of the levee surface loss is about 4 mm.

The longitudinal section of the left tunnel is shown in Figure 10(b). The shield in the left line has not yet been dug into the flood control levee, and at this time, the displacement of the surrounding rock of the shield shell part is similar to the surrounding rock of the tunnel in the right line. Since the soil above the tunnel in the left line is dominated by round gravel and plain fill, and the thickness of the soil layer is lower, the maximum displacement settlement lost in the surface above is smaller than that in the right line, which is about 3 mm.

The displacement convergence diagram of the right tunnel palm face is shown in Figure 10(c). From the tunnel vault to the surface, the displacement of the surrounding rock is gradually decreasing, and the overall appearance is elliptical. The displacement at the top of the tunnel arch is the largest, about 6 mm, and the displacement of the elevation arch is relatively small.

The displacement convergence diagram of the left line tunnel palm face is shown in Figure 10(d). Slightly different from the right line palm face, the surrounding rocks on the right side of the palm face of the left tunnel, especially the looser surrounding rocks on the upper right side, appeared to be significantly larger than the displacement of the surrounding rocks at the relative position on the left side of the palm face. The reason for the above phenomenon is that the right shield tunnel was constructed before the left shield tunnel, and the excavation of the right shield tunnel and the high-strength concrete support changed the gravity field and seepage field of the original surrounding rock. During the construction of the left line shield tunnel, attention should be paid to the displacement and deformation of the surrounding rock on the right side of the tunnel.

The top view of the tunnel model displacement is shown in Figure 10(e). When the shield machine in the right tunnel crossed directly below the Liuyang River flood control levee, the displacement plastic deformation caused by the shield excavation was not large in scope, and the ground loss settlement was within a relatively small range.

4.2.2. Before the Left Line Tunnel Crosses the Flood Control Levee

When the shield machine in the left line is bored directly below the flood control levee, the stability of the levee is affected by both the shield machine boring construction in the left line tunnel and the structure of the completed tunnel in the right line. The simulation results are shown in Figure 11.

The longitudinal section of the right tunnel is shown in Figure 11(a). Influenced by the seepage field, the displacement of the surrounding rock in the upper part of the shield shell was significantly higher than that when it did not cross the Liuyang River previously. The plastic displacement range of the tunnel-surrounding rock is significantly expanded, with a maximum displacement of about 8 mm. At the same time, the uplift deformation of the surrounding rock in the lower part of the shell has expanded. The displacement deformation of the flood control levee on the right tunnel profile is around 3.5 mm due to the shield tunneling in the left line.

The longitudinal section of the left tunnel is shown in Figure 11(b). The form of displacement of the surrounding rock near the shield housing is similar to the previous displacement. The maximum settlement displacement for surface loss of the flood control levee is about 4 mm. Unlike before, the range of plastic displacement is significantly expanded. These ranges are mainly the surrounding rocks near the two shield housings and between the two shields, with a displacement size of about 2 mm, as shown in the red dashed box in the figure. The surrounding rock between the two shields is just below the flood control levee. During the actual construction, the site is very noteworthy, and the surrounding rock between the two shields tended to show a wide range of overall displacement.

The displacement convergence diagram of the right line tunnel palm face is shown in Figure 11(c). Compared with before crossing the Liuyang River, the displacement of the surrounding rocks above the palm face became significantly larger, with a maximum displacement of about 8 mm. The form of the surrounding rock deformation is also no longer elliptical but funnel-shaped expanding outward.

The displacement convergence diagram of the left line tunnel palm face is shown in Figure 11(d). Influenced by the shield construction of the right line tunnel, the high-strength concrete lining support, and the overlying soil layer, the displacement trend of the palm surface of the left line tunnel appeared to be significantly to the right. The bridge piles not far away also showed significant plastic deformation displacement.

The top view of the tunnel model displacement is shown in Figure 11(e). After the right line shield is dug below the Liuyang Riverbed, the value and range of plastic displacement are expanded. The excavation of the shield in the left line has caused a displacement of the plain fill layer of the flood control levee in favor of the middle of the two tunnels, and part of the plastic deformation zone caused by the construction excavation of the two shields is continuously penetrating. The through area is within the plain fill layer just below the levee slope, as shown in the red dashed boxes in Figures 11(b) and 11(e).

4.2.3. Both Tunnels Have Crossed the Flood Control Levee

When the shield of the left line passes through the flood control levee, the stability of the levee is less affected by the construction of the shield. The shield machines in both tunnels have been driven under the riverbed of the Liuyang River. The simulation results are shown in Figure 12.

The longitudinal section of the right tunnel is shown in Figure 12(a). The displacement of the surrounding rock at the shield shell of the right line is basically similar to the situation when the shield of the right line first enters below the riverbed of the Liuyang River.

The shield of the left line is excavated to the junction between the base of the flood control levee and the riverbed, as shown in Figure 12(b). The slope of the flood control levee shows obvious plastic deformation toward the slope, and the displacement deformation of the plain fill layer is about 3 mm. The displacement of the upper part of the shield housing is the largest, with a maximum displacement of about 7 mm. The displacement deformation at the base of the flood control levee is about 5 mm.

The displacement convergence diagrams of the right and left tunnel palm faces are shown in Figures 12(c) and 12(d), respectively. The displacement deformation form of the surrounding rock at the palm face of both tunnels shows a funnel shape expanding outward. As the shield machine has passed through a more complex geological area, the surrounding rock of the left tunnel no longer has a significant plastic displacement deflected to one side, but the shield excavation of the left tunnel still has an impact on the lining structure of the right tunnel. Due to the seepage field and the nature of the geotechnical body, the displacement value and range of the right tunnel palm face, which is already under the riverbed, are significantly larger than those of the left tunnel palm face.

The top view of the tunnel model displacement is shown in Figure 12(e). Due to the change in the geological stratum in which it is located, the geotechnical body between the two shields no longer has a penetrating plastic deformation zone. When the two shields passed through the levee, the slope of the levee directly above the tunnel showed a large range of displacement deformation, as shown in the red dashed box in Figure 12(e).

Based on the results of numerical simulation, it can be inferred that during the shield construction, the levee base and slope of the flood control levee on the west side of the Liuyang River are highly susceptible to large geotechnical displacement deformation when the left line shield machine crosses the flood control levee.

5. Conclusions

The deformation of the flood control embankment under the condition of shield tunneling across the river is studied. Taking the shield construction of Liuyang River as the engineering background, the multifield coupling software COMSOL is used to establish the shield construction model. The main conclusions are as follows: (1)When the shield machine of the left line is excavated below the flood control levee, the displacement of the surrounding rock at the palm face of the tunnel is obviously biased towards the completed tunnel on the right and has an impact on the bridge piles not far away(2)When the shield machine of the left line is dug into the levee foundation of the flood control levee directly below, the levee foundation is disturbed by two excavations, and a large range of deformation occurs. The deformation value and range of the geotechnical body are higher than in other areas(3)The tunneling construction speed of the shield machine should be reasonably arranged and the construction distance between the palm faces of the shield machine lengthened, especially when it enters the range of the flood control levee. In the process of right line tunnel construction, a good job should be done of the grouting support for the left curved wall of the right line tunnel(4)In case of special conditions, the secondary lining support of the tunnel will be done in time. Particular attention should be paid to the two points in time when the left line shield crosses the flood control levee and when the shield crosses the base of the flood control levee

Data Availability

Data are available on request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant numbers 11772358 and 52208333) and the Water Conservancy Science and Technology Major Project of Hunan Province, China (Grant No. XSKJ2019081-10).