Study on mechanical response and crack development law of tunnel

： Considering the linear elastic fracture mechanics and complex function theory, considering the non hydrostatic pressure field, the stress solutions of lining and surrounding rock of deep buried circular underground cavern with lining under the action of external internal pressure are studied. The fracture mechanics model of underground cavern supporting structure is established, and the numerical verification is carried out. The analytical solution of circular pressure cavern with lining shows that when the internal pressure is 0, it can degenerate into the existing classical solution. The results of finite element calculation and analytical calculation show that the tensile stress of lining is very large at the vault and arch bottom under the action of strong external supporting force. The calculation shows that the cracks of lining will expand and open under the action of tensile stress, and the lining is in the shape of "flat duck egg", However, if there is no through rupture, there will be great stress concentration at the top and bottom of the arch, and the tensile stress value is much larger than that calculated by elastic mechanics. Therefore, it is a good method to configure a certain amount of reinforcement to enhance the stiffness and tensile strength of the lining.


Introduction
The theoretical research on underground circular tunnel excavation and support started from the famous Fenner formula [1] and Kastner formula [2] , and these classical solutions have still been extensively applied so far. Ren et al. [3][4] pointed out that in the deeply buried hydraulic tunnel, the Fenner formula should have different expressions under different internal pressure loads, and then they derived the elastic-plastic solution to a lined circular tunnel on this basis. Furthermore, for a deeply buried underground circular tunnel with a high internal water pressure, the analytical elastic-plastic solution to the deeply buried unlined circular cavern was derived in consideration of the influence of seepage field [5] . Wang and Li [6] derived a new elastic-plastic solution to the lined circular pressure tunnel under the hydrostatic pressure field, and this solution, which considered the release of ground stress and time sequence of lining installation, could be reduced to classical solutions like Fenner and Kirsch. Taking a deeply buried hydraulic pressure tunnel as the study object, Wang et al. [7] derived the elastic-plastic solution to the circular tunnel with internal pressure based on the D-P criterion, and analyzed the influence of change in hydraulic pressure on the lining and plastic zone of surrounding rock. Based on the twin-shear unified strength theory, Fan W and Yu [8] derived the analytical elastic-plastic solution to a pressure tunnel in consideration of material shear dilatancy and softening. Considering the nonhydrostatic pressure field, Chu et al. [9] derived the elastic stress solution and viscoelastic stress solution of surrounding rock and lining in a deeply buried circular tunnel, but the influence of internal pressure was neglected. After considering the heat exchange of compressed air with lining and surrounding rock, Zhou and Zhou et al. [10][11] derived the formulas for the stress field and temperature field distributions of the surrounding rock and lining in an underground compressed air energy storage cavern under typical operating conditions. The lined underground circular cavern tunnel is prone to cracking, thus destructing the integrity of lining structure and greatly shortening the service life of tunnel. The generation and propagation of lining cracks are a very complex nonlinear process. Based on a statistical analysis of lining crack distribution characteristics in the existing highway tunnels, Huang et al. [12] used the extended finite element method to study the influences of factors like eccentric compression, cavity behind lining and relaxation ground pressure on the crack propagation and distribution characteristics. Zhang et al. [13] simulated the crack distribution in reinforced concrete using a discrete crack model, and analyzed the influence of reinforcement amount on the crack width. Zhu et al. [14] analyzed the characteristics and mechanisms of construction cracks, eccentric compression cracks and differential settlement cracks of a double-arched tunnel, and proposed the corresponding prevention and control measures.
Most of the above studies have assumed that the underground circular cavern is under a uniform hydrostatic pressure field, but actually, most deeply buried circular caverns are under the nonhydrostatic pressure field. The Kirsch solution is a classical solution to unlined circular caverns in the nonhydrostatic pressure field, but this solution does not consider the effect of lining. As for the theoretical derivation in this study, it was assumed that the underground cavern was under the nonhydrostatic pressure field, the effect of high pressure inside the cavern was taken into account, the stress solutions to the lining and surrounding rock in the lined underground circular cavern were derived, followed by the finite element verification. On this basis, a numerical model of lining cracking and crack propagation under the action of internal pressure was established through the extended finite element method.

Mechanical Response Analysis of Surrounding Rock and Lining Structure
A mechanical planar strain model was constructed for the deeply buried circular cavern under the joint action of nonhydrostatic pressure field and high internal pressure, as shown in Figure 1. As shown in Figure 1, assume that the lining is sealed and keep a close contact with the surrounding rock, the frictional force on the contact surface between lining and surrounding rock is neglected, the outer diameter r 0 is equivalent to the tunnel excavation radius, and the air storage cavern is located in a nonhydrostatic pressure field with the lateral pressure coefficient of λ. The boundary conditions of the model are as follows: The analytical function of surrounding rock, the boundary conditions of which are combined, is given through the following equation [15] : (2) The stress solution to the surrounding rock is as follows: (3) With the boundary conditions combined, the analytical function of lining is given by the following equation: The stress solution to the lining is as below: (5) γ, β, δ, H and A1-A6 are expressed by Equation (6) and Equation (7): where σcθ and σθ are the circumferential normal stresses of lining structure and surrounding rock, respectively; τcθ and τθ are the shear stresses of lining structure and surrounding rock, respectively; Gc and G represent the shear moduli of lining structure and surrounding rock, respectively; κc and κ are the Poisson's ratios of lining structure and surrounding rock, respectively.

Example Verification and Finite Element Analysis
In order to verify the correctness of theoretical solutions given in this study, the theoretically derived results will be compared with the finite element calculation results.

Analytic calculation results
In the calculation, the parameters were selected as follows: The outer radius and inner radius of lining were r 0 =2.4 m and r 1 =2.05 m, respectively, the internal pressure of cavern was p i =10 MPa, the vertical pressure at cavern top was P=6.3 MPa, the elasticity moduli of surrounding rock and lining structure were E=5 GPa and E r =30 GPa, respectively, the Poisson's ratios of surrounding rock and lining structure were κ=0.32 and κ c =0.2, respectively, and the lateral pressure coefficient of surrounding rock was λ=0.47. Figure 2 presents the radial and circumferential stress distributions within the lining range at the angle of θ=90°. It could be known that the maximum circumferential tensile stress of lining was 29.33 MPa, appearing at inner side of lining top. With the increase in the radius r, the circumferential tensile stress of lining was gradually reduced, and that at the outer side of lining was 16.84 MPa. Under the radial compression of lining, the radial pressure inside the lining was equivalent to internal pressure (10 MPa) of cavern. As the radius was increased, the radial pressure gradually declined, and that at the outer side of lining was 5.14 MPa. The radial stress and circumferential stress distribution laws of lining at different angles are presented in Figure 3- Figure 4. It could be observed from Figure 3 that the radial stress of lining changed little with the angle, the internal surface pressure of lining was closely related to the change of angle. On the internal surface of lining, the circumferential tensile stress was the maximum (about 29.33 MPa) at the apex of arch (θ=0°). On the external surface of lining, the circumferential tensile stress reached the maximum value of about 23.98 MPa at the arch foot (θ=180°).
The radial stress and circumferential stress distributions laws of surrounding rock at the angle of θ=90° are shown in Figure 5. It could be seen that the circumferential tensile stress of surrounding rock was 2.49 MPa under r=2.4 m. As the radius was increased, the circumferential tensile stress of surrounding rock was gradually reduced, declined to 0 nearby r=2.8 m, and then was turned into compressive stress. When r tended to be infinite, the circumferential tensile stress of surrounding rock tended to be λp=2.96MPa, according with the assumption of boundary conditions. When the surrounding rock was under radial compression, the radial stress of surrounding rock was -5.14 MPa under r=2.4 m, being equivalent to the radial stress at the outer side of lining, this conformed to the assumption of lining-surrounding rock contact. With the increase in the radius, the radial pressure of surrounding rock was firstly reduced and then increased. As r tended to be infinitely great, the radial stress of surrounding rock was approximate to p=-6.3 MPa. Figure 6-Figure 7 give the radial stress and circumferential stress distribution laws of surrounding rock at different angles. It could be observed that the radial stress of surrounding rock was always compressive stress, which was firstly increased and then decreased with the change of angle. The radial stress of surrounding rock reached the maximum value at θ=180°. As the angle was enlarged, the circumferential stress of surrounding rock was transformed from tensile stress into compressive stress and then into tensile stress again, and the circumferential tensile stress was the maximum (about 1.7 MPa) at θ=0°. Figure 8 displays the change laws of maximum tensile stress of lining under different internal pressures. It could be seen that as the internal pressure was increased, the circumferential stress of lining presented an approximate linear growth, and the maximum circumferential stress of lining was gradually transformed from compressive stress into tensile stress under about pi=5 MPa.

Finite element simulation
The finite element modeling and calculation were conducted via ABAQUS to simplify this into a 2D planar strain problem. The modeling was done within the xy plane, the boundary of axis x in the model ranged from -30 to 30, and the boundary of axis y ranged from -30 to 30, too. The mesh generation was realized within the whole model plane, the mesh elements were four-node bilinear plane stress quadrilateral elements, totaling 6,952 elements. The boundary and loading conditions were set as follows: The left and right sides of the model constrained its displacement in the direction x, and the model bottom constrained its displacement in the directions x and y; a vertical load of P=6. 3 MPa was set at the top of model, the internal pressure of q=10 MPa was applied inside the cavern, and the gravity of model was not considered, as shown in Figure 9. The idea and procedures of finite element calculation are as follows: (1) Establish the components representing the underground circular tunnel and lining, and generate the suitable meshes; (2) Set tangential smooth contact under the only existence of normal contact between the lining and surrounding rock; (3) Apply displacement boundary conditions and corresponding load; (4) Operate the calculation. The maximum principal stress distribution laws of surrounding rock and lining structure, which are calculated through the finite element method, are displayed in Figure 10. In ABAQUS, the tensile stress was positive and compressive stress was negative. It could be known that the maximum tensile stress borne by the lining structure was 28.37 Mpa, appearing at the side close to the cavern on the top of lining structure. Figure  11 shows the minimum principal stress distribution laws of surrounding rock and lining structure, which are also calculated through the finite element method. As shown in the figure, the maximum compressive stress of lining was about 9.98 Mpa, which accorded with the set initial loading condition. The minimum principal stress of surrounding rock presented a "butterfly"-like distribution, the radial compressive stress of surrounding rock was reduced obviously at upper and lower positions, and along the directions of positive and negative 45°, the radial compressive stress of surrounding rock changed uniformly.
The maximum tensile stress distribution of lining along the apex of arch is shown in Figure 12. It could be seen that both theoretical solution and numerical solution of tensile stress reached the maximum values (29.33 Mpa and 28.31 Mpa, respectively) on the internal surface of lining. As the distance from the center of cavern was enlarged, the maximum tensile stress of lining was gradually reduced, and the theoretical solution differed little from the numerical solution all the time, so they could be mutually verified.

Extended finite element simulation
According to the finite element modeling process in section 3.2, the extended finite element simulation was still done in ABAQUS. The extended finite element (XFEM) was formally proposed by Professor Belytschko of U.S. Northwestern University in 1999 [16][17] . It follows such idea: The discontinuous displacement in the crack propagation is described by introducing a nonlinear function nearby the cracks.
The XFEM method forms a unit decomposition function on the solution domain through the finite element shape function, defines the displacement mode nearby the node xi, and approximates the infinite dimensional problem into a finite dimensional problem, and the displacement on the solution domain is approximately expressed as follows: where u app is the displacement at any point in the domain, i is the unit domain of crack influence, N i (x) is defined as a group of shape function on the influence domain, β ij is the nodal displacement in the influence domain, and p j (x) is a group of base function in the influence domain. In the XFEM, the displacement jump (discontinuous) in the crack propagation is usually characterized using the Heaviside step function, and the singularity of stress field at crack tip is characterized using the Westergaard function.
The simulation parameters and process of ABAQUS XFEM referred to the analysis example in 3.1 and conventional finite element simulation in 3.2. Especially, in the simulation of lined crack tip, the fracture energy was G Ι f =G Ⅱ f =G Ⅲ f =80N/m according to the theory of concrete fracture mechanics. It was assumed that the cracks were propagated only in the concrete lining [16][17][18][19] . Figure 13 shows the displacement distribution laws of circular lining under the action of high-pressure gas and prefabricated cracks. It could be observed that the lining was mainly extruded in the vertical direction and pulled open in the horizontal direction, presenting a "flat duck egg" shape. The cracks at the top and bottom of lining were further propagated and stretched towards the deeper part along the already existing cracks, but no transverse cracks appeared, and the maximum opening displacement of the internal surface of lining was 0.93 mm.
The stress field distribution of circular lining under the effect of high-pressure gas and prefabricated cracks is shown in Figure 14. It could be known that due to the existence of cracks, the stress around the cracks reached as high as 94.77 MPa, with a strong stress concentration phenomenon, the stress at the part slightly away from the crack was approximate to the stress when the crack influence was not considered, and this was different from a viewpoint in the classical linear elastic fracture mechanics, namely, the stress at near end of crack tended to be infinitely great.

Conclusions
(1) The derived analytical solution to the lined circular pressure tunnel under the nonhydrostatic pressure field can be reduced to the existing solution under the internal pressure of p i =0, the maximum circumferential stress of lining is closely related to the internal pressure pi, and with the increase of p i , the circumferential tensile stress of lining is increased rapidly, and a critical pressure exists. When p i is greater than the critical value, the maximum circumferential stress of lining is tensile stress. When p i is smaller than the critical value, the maximum circumferential stress of lining is compressive stress.
(2) The analysis example and finite element numerical example both indicate that under the effect of high-internal-pressure gas, a large tensile stress will appear at the top of lining, where transverse cracks may be easily formed, which will result in the cracking failure of lining and shorten the service life of tunnel.
(3) The XFEM simulation results manifest that under the joint action of nonhydrostatic pressure field and high internal pressure, apparent cracking will occur to the concrete lining, but the cracks fail to run through the lining. The largest open crack is about 0.93 mm, being 0.2 mm larger than that specified in the seepage control requirement. It is suggested that the concrete cracking should be controlled by strengthening the reinforcement during the actual construction of underground internal pressure tunnel.

Date of reception:
About the author: Dezhu Tan (1971-), male, born in Shiyan, Hubei Province, Bachelor's Degree, senior engineer, mainly occupied in the scientific research on the highway tunnel engineering