Hydro-mechanical response of excavating tunnel in deep saturated ground

Excavating a tunnel in a deep and saturated ground affects the short- and long-term hydro-mechanical (H-M) response in the ground surrounding the opening. However, the interactions between transient pore pressure behavior and the corresponding deformation and stresses in the ground ahead of and behind the tunnel face are still not well understood. This paper investigates the transient H-M response of excavating a tunnel in a deep and saturated ground using a two-dimensional axisymmetric coupled model in the computer code Fast Lagrangian Analysis of Continua (FLAC). The tunnel was advanced in a stepwise excavation procedure consisting of undrained excavation and drained consolidation until the final tunnel face was reached. The final excavated face was then left to consolidate toward the steady-state condition. The main results of the paper are as follows: (1) when simulating a tunnel excavation in deep saturated ground using the convergence-confinement method, the unloading factors should be nonlinear and should consists of the mechanical unloading factor in the form of excavation force and the hydraulic unloading factor in the form of excavation pore pressure. These two unloading factors are necessary because the induced H-M response near the tunnel face is a rather transient response instead of an initial or final response. Moreover, it is observed that the pore pressure dissipation is not linear either with time or with distance to the tunnel face, (2) a relationship between the unloading factors and the distance to the tunnel face should then be established. This relationship is vital because it will provide the timing for tunnel support installation, and (3) the extrusion and the convergence of the advance core could be related through the proposed equations capturing the linear relationships between the face extrusion and its convergence as well as between the core extrusion and its pre-convergence. Through these relationships, the tunnel engineer may be able to estimate the magnitude of the deformation ahead of the face, which will subsequently allow control of the deformation behind the face.

a time-dependent seepage and consolidation induced by the transient nature of the coupled stress-pore pressure interaction. Even after the final excavated face is reached, ground consolidation in the form of slow dissipation of pore pressure, will continue during the so-called standstill period until the steady state condition is reached [11]. Thus, the advancing tunnel will result in two hydro-mechanical (H-M) responses. First, it will result in short-term response from when the ground is being progressively excavated (undrained loading). Second, it will result in long-term response from when the ground is consolidating (drained loading).
While short-term response is usually concentrated near the tunnel face and is dependent upon the techniques of excavation, the long-term response is usually noticeable long after the completion of the excavation. However, the degree to which the excavation affects these short-and long-term H-M responses of the surrounding ground remains unclear. How the transient pore pressure behavior will affect the corresponding deformation and stresses in the ground ahead of and behind the tunnel face is also poorly understood.

Importance of coupled H-M response in tunneling
In practice, an advancing tunnel is often simulated with an uncoupled approach, neglecting the transient coupled response. For instance, the so-called steady-state approach has been used to study the effects of fluid flow on tunnel deformations [9,12]. In this approach, the whole excavated length is reached in one step (instead of a progressive step-by-step excavation). The model is then brought to geomechanical equilibrium under the steady state condition (instead of under time-dependent pore pressure distribution). Consequently, the uncoupled approach cannot yield the transient deformation behavior that supposedly results from the advancing tunnel because the coupling effect is not properly captured. It is imperative to predict the transient deformation behavior correctly in order to determine when the support should be installed [13].
Furthermore, while it is important to know the magnitude of wall displacement behind the tunnel face, it is also vital to understand the coupled response that has already started in the ground ahead of the face. When a tunnel is being progressively excavated, a substantial volume of the ground ahead of the face is removed. This operation converts the triaxial state of stress to a plane stress condition that causes loss of confinement in the face. Consequently, radial deformations will be induced in the ground behind the tunnel face, called the convergence, and in the ground ahead of the face (the advance core), called the pre-convergence ( Fig. 1). At the same time, the loss of axial confinement will also trigger outward axial deformation of the face, called the extrusion. Unfortunately, in deep saturated ground, this type of deformation behavior of the advance core has not been fully explored. The relationships between wall convergence and core extrusion considering the transient coupling effect have not been clearly identified.

Objective of the study
The main objective of this paper is to perform a numerical study on the transient H-M response of excavating a tunnel in a deep and saturated ground. The H-M responses of the ground will be analyzed as well as their transient behaviors during the face excavation (referred as the short-term response) and the standstill periods (referred as the long-term response). The practical application of the results for tunneling in deep saturated ground are discussed, including the relationships between wall convergence and core extrusion.
A flowchart of the steps taken to fulfill the objective of this study can be seen in Fig. 2. As shown, the conceptual numerical model (before tunnel advance) is built first. The conceptual model can be seen later in Fig. 4b. The hydraulic and geomechanical equilibrium conditions for the model are then established with the prescribed boundary conditions and in situ stresses. The excavation is then started by removing the elements in the excavation section. After successfully performing model calibration against the existing analytical solutions, the numerical model is then used to assess the H-M response of an advancing tunnel in deep saturated ground.

Numerical modeling of an advancing tunnel
An advancing tunnel is represented as an axisymmetric coupled model that is created in a 2-D explicit finite difference code called Fast Lagrangian Analysis of Continua or FLAC [5]. The program can perform fluid flow and geomechanical calculations based on Biot's theory of poroelasticity [2]. The computation starts from geomechanical equilibrium. From this state, the solution to a coupled H-M problem for the next time step is achieved simultaneously by first performing the flow calculation to evaluate the increased pore pressure and then the geomechanical problem to evaluate the volumetric strain contribution. To evaluate the effective stress and observe ground failure, the total stress correction is then applied [6]. The total stress correction is performed to accommodate pore pressure rising from fluid flow calculation and from the mechanical volumetric strain. Figure 3 shows the numerical steps carried out in FLAC to incrementally solve Biot's coupled H-M equations. At any time t, FLAC solves the dynamic equilibrium equations in terms of the total stresses σ ij , the displacement u i , applied boundary loads F i , and the body loads ρg i z, as where ρ is the unit weight, g i is the gravitational acceleration, and z is the vertical coordinate of a mass. A superposed dot indicates derivative with respect to time. From the velocities u i , the strain rates ε ij are determined from the velocity gradient through the following compatibility relation:  From the strain rates, the effective stress increments σ ′ ij are calculated using the constitutive relation D ijkl as The volumetric strain increment ε vol =ε kk is also calculated from the strain rates and is then passed onto the fluid diffusion equation as where p is the pore pressure, k ij is the permeability tensor, µ, ρ f and K f are the fluid viscosity, unit weight and bulk modulus, respectively, S is the storativity, h is the fluid elevation, φ is the porosity, and q is the fluid source or sink. The pore pressure increment ṗ and the effective stress increments σ ′ ij are then added to get the total stress increments σ ij according to Terzaghi's effective stress law defined as where δ ij is the Kronecker delta tensor. In the next time step t + dt, the new total stress increments σ ij are then used in the dynamic equilibrium in Eq. (1) and the calculation procedure is repeated for the desired length of time.
In this paper, the tunnel is modeled as an axisymmetric and unsupported circular tunnel excavated at a depth of H = 225 m and with a radius of R = 2.5 m. An axisymmetric model is used in this paper so that the effect ground surface can be disregarded. Therefore, the analysis performed in this paper is intended for unsupported, full-face circular The problem layout and the H-M boundary conditions can be seen in Fig. 4a. The model prescribes pore pressure p = 0 or a permeable boundary at the tunnel perimeter to represent drainage, pore pressure p = p(t) or an impermeable boundary on the left and on the bottom boundaries to allow the pore pressure to change with time, and pore pressure p = p o on the upper and on the right boundaries to represent far-field pore pressure.
The advance of the tunnel is modeled by gradually deactivating 40 finite difference zones of thickness x = 1.25 m (R/2) for over 10 days of excavation period, indicating the excavation rate of v a = 5 m/day (Fig. 4b). When each zone x is being deactivated, an instantaneous tunnel excavation (excavation period is assumed to be t = 0) is carried out to represent the undrained loading. To represent the drained loading, the excavation is then followed by letting the consolidation takes place with a period of t = x/v a . Both loadings are carried out under the permeable wall and tunnel face (p = 0). This stepwise excavation is continued until the excavation arrives at the final tunnel face, that is at x/R = 0 (see Fig. 4b). At this point, the response of the ground surrounding the tunnel is called the short-term response. The ground consolidation is then allowed to occur until the standstill period of 930 days is achieved. At this point, the response of the ground surrounding the tunnel is called the long-term response, meaning the ground has reached its steady state condition.
The characteristic time of the diffusion process t c is 5 days, which is calculated as where k H is the ground conductivity, γ w is the unit weight of water, and S is the ground storativity calculated as S = n K w + 1 K + 4G 3 , where n is the ground porosity, K w is the water bulk modulus, and K and G are the ground's bulk and shear moduli, respectively. The H-M parameters of the ground are given in Table 1.

Validation against analytical solution
Before analyzing the transient H-M response, the axisymmetric model in this paper is first validated against the analytical solutions in Bobet [3]. For the purpose of validation and direct correspondence with the analytical solutions, the so-called uncoupled simulation is used. In this approach, the excavated region (L = 50 m) in Fig. 4b is excavated in one step and is done without calculating the contribution of volumetric strain toward pore pressure generation. Thus, there is no excess pore pressure buildup during the short-term consolidation. The model is run in the "flow-only simulation" for a total period of t = 930 days and then is brought to mechanical equilibrium. In the short term, because the excavation does not introduce excess pore pressure (Fig. 5a), changes in stresses and deformation (Fig. 5b, c) are only affected by the mechanical response of the ground after the excavation. As expected, the radial deformation increases in magnitude toward the opening. In the long term, as the seepage force toward the tunnel boundary decreases, the impact of effective stresses on the radial displacement occurs not only around the opening but also at a great distance from the opening. Consequently, radial displacement increases with the increase in radial distance from the tunnel (Fig. 5c). This result validates the axisymmetric model used in this paper.

Contours of the induced H-M response
The contours of pore pressure, radial displacement, and effective stresses are plotted in Fig. 6. These contours present wide-ranging views of the induced short-term H-M response shortly after the excavation (at the start of the standstill period) and long-term H-M response after the drained consolidation (at the end of the standstill period). From the contours of the pore pressure field (Fig. 6a), there appears to be a transition from the undrained condition shortly after the final face excavation to the drained condition after the steady state is reached. In the short-term, the undrained condition is illustrated by a steep gradient of pore pressure concentration surrounding the tunnel face at x/R = 0, and the gradient extends far behind the tunnel face at x/R = 20 with gradually decreasing magnitude. In the vicinity of the tunnel face, the instantaneous pore pressure develops from 0 to 2 MPa within a short distance of Comparison of a pore pressure p, b radial displacement u r , and c effective stresses σ ′ in the short-term and long-term consolidations between the present axisymmetric model and the uncoupled analytical solutions from Bobet [3] only − 1R from the face (0.8 MPa/m). At x/R = 20, the pore pressure gradient gently decreases to 0.2 MPa/m within 4R from the tunnel axis.
In addition to the steep pore pressure gradient, the undrained loading develops a region of excess pore pressure further inside the advance core (from − 1 to − 4R), which increases the in-situ pore pressure by 11%, from 2.25 to 2.50 MPa. This highlyconfined region may be a favorable condition, as the confinement will temporarily strengthen the ground against shearing during the excavation. In the long term, the confinement disappears with time as the excess pore pressure also dissipates, transforming the ground toward the steady-state condition. This transient hydraulic response consequently affects the corresponding mechanical response of the ground in the form of radial displacement (Fig. 6b). Even though the radial displacement at the tunnel wall only slightly increases with pore pressure dissipation, the radial displacement in the ground above the tunnel is greatly affected by the seepage flow toward the opening. In the ground near the top boundary, the radial displacement increases from about 1 cm in the short term to a maximum of 6 cm in the long term.
The increase in deformation with time is somewhat expected. Due to the seepage flow in the vicinity of the face, the effective radial stress increases from 1.5 MPa in the short term to more than 2.0 MPa in the long term (Fig. 6c). Similarly, due to the dissipation of pore pressure, the contour of the effective tangential stress in Fig. 6d expands from the proximity of the opening in the short term (within 4R from the tunnel axis and − 2R from the tunnel face) to be more distributed toward the advance core in the long term (up to − 8R from the tunnel face). Consequently, long-term core stability will be affected. If the load transfer exceeds the strength of the core, the tunnel face will collapse, followed by the collapse of the opening.

Longitudinal profiles of the induced H-M response
Profiles of pore pressure, convergence, and effective stresses taken at 1 m above the tunnel wall are plotted along the longitudinal axis of the tunnel (Fig. 7). The comparison between the coupled and uncoupled simulation of these H-M responses is also shown.
From Fig. 7a, the coupled simulation clearly induces the non-uniform pore pressure profile at the tunnel wall. As can be seen previously from the contour plot in Fig. 6a, the excavation induces a sudden or short-term increase in pore pressure up to 2.65 MPa in the region just ahead of the tunnel face (from − 1R to − 4R). In this region, due to the nature of ground diffusivity, the pore pressure ahead of the face dissipates at a slower speed than the excavation progress. This condition results in a higher stress concentration and a stiffer core, which is then responsible for the excess pore pressure buildup in the corresponding region inside the core. The pore pressure build-up gradually subsides with time towards its steady-state condition as shown by its long-term profile. This profile corresponds to the long-term pore pressure contour in Fig. 6a.
The behavior of the non-monotonic pore pressure buildup beyond its initial value is called the Mandel-Cryer effect [4,8]. This phenomenon has also been observed at the tunnel crown when tunneling in shallow ground [10]. This behavior occurs due to the stress concentration in the tunnel crown that is still stiff and undrained. This pore pressure characteristic also highlights the difference between the coupled simulation in this paper and the uncoupled simulation that has been performed by previous researchers [9,12]. In the uncoupled simulation, the Mandel-Cryer effect does not exist (Fig. 7a), and the pore pressure diffuses faster than in the coupled simulation. In the short term, at the tunnel face x/R = 0, the pore pressure drops to 0.7 MPa in the uncoupled simulation but is still maintained at 1.5 MPa in the coupled simulation.
The difference in the pore pressure behavior between the coupled and uncoupled simulations has direct consequences for the corresponding mechanical responses of the ground in terms of displacements and stresses. The coupled simulation results in smaller convergences than the uncoupled one does, particularly for displacement in long-term consolidation (Fig. 7b). In this state, the convergences from the coupled simulation are only 0.1 cm, 1.4 cm and 4.5 cm at locations ahead of the face (x/R = − 20), at the face (x/R = 0) and behind the face (x/R = 20), respectively. These amounts are 80%, 18% and 12% lower than those from the uncoupled simulation: 0.5 cm, 1.7 cm and 5.1 cm, respectively. In the short term (Fig. 7b), the effect of core stiffening due to the Mandel-Cryer effect can be seen at x/R = − 2, which is the location of the highest pore pressure buildup (see Fig. 7a). At this location, the amount of displacement from the coupled simulation is only 0.2 cm, which is 50% lower than that from the uncoupled one (i.e., 0.4 cm), reconfirming the confinement effect resulting from excess pore pressure buildup (see explanation for Fig. 6a).
Compared to the displacement behavior, the direct consequences of the coupling mechanism on the radial and tangential stress behaviors are rather complex and dependent on the amount of consolidation time. In this tunneling problem, there are two competing effects that contribute to the induced H-M response. The first is the fluid flow effect: the seepage flow toward the excavation induces seepage forces, increasing the Fig. 7 Longitudinal profiles of a pore pressure p, b radial displacement u r , and c effective stresses σ ′ : comparison between the coupled and uncoupled simulations effective stresses of the surrounding ground. The second is the geomechanical effect: at the same time as the seepage flow, the mechanical response of the ground will induce an increase in pore pressure, reducing the effective stresses. In the advancing tunnel problem, the degree of each competing effect is time-dependent in nature, consequently inducing transient H-M responses of the surrounding ground as the tunnel progresses. In short-term consolidation (Fig. 7c), the coupling interaction is dominantly governed by the geomechanical effect. The degrees of the effective stresses induced by the coupled and uncoupled simulations are rather similar, despite the lower pore pressure value induced by the latter simulation. The opposite occurs in long-term consolidation, that is, the fluid flow effect outweighs the geomechanical effect. The increase in pore pressure from the geomechanical response is now not as significant as it was in the short-term consolidation. In the coupled simulation, the increase of effective stresses due to fluid flow toward the excavation is now lower than in the uncoupled one, which now corresponds to the lower pore pressure value induced by the uncoupled simulation. The difference is particularly obvious for the effective stresses in the ground ahead of the tunnel. The comparison of the induced effective stresses at different consolidation times from both simulations is presented in Table 2.

Transient H-M response
To capture the transient response experienced by the tunnel wall before and after the excavation passes, a monitoring point is installed at x/R = 10 and at a location 1 m above the tunnel crown during the excavation period. Figure 8a shows that the pore pressure build-up at the monitoring point starts to develop at some distance before the excavation face arrives at the point.
At t = 3 days (x/R = 14), the pore pressure starts to deviate from p o = 2.25 MPa to its peak value of p = 2.65 MPa at t = 4.25 days (x/R = 11.5). The pore pressure then drops to p = 1.75 MPa when the tunnel face arrives at the monitoring point at t = 5 days (x/R = 10). As the excavation passes the monitoring point, the pore pressure falls sharply before it steadily decreases at t = 6 days (x/R = 8). This Mandel-Cryer effect can also be observed in Fig. 7a.
Similarly, at t = 3 days (x/R = 14), the convergence also starts to deviate from u r = 0 cm to the value of u r = 3.0 cm at t = 5 days and x/R = 10 (Fig. 8b). As the tunnel face continues to advance, the convergence keeps increasing before flattening at u r = 4.2 cm at t = 7 days and x/R = 6. At the same time, the effective stresses start to deviate from σ ′ = 2.25 MPa at t = 3 days (x/R = 14) and remain steady at σ ′ θ = 4.65 MPa and σ ′ r = 1.25 MPa at t = 7 days and x/R = 6 (Fig. 8c).

Consolidation
Simulation At x/R = − 20 (ahead) At x/R = 20 (behind)  The significance of the results is threefold. First, during tunnel advance, the induced H-M response near the face is a rather transient response instead of an initial or final response. Second, the pore pressure dissipation is not linear either with time or with distance to the tunnel face. Rather, the induced pore pressure is seen to increase beyond its initial value before the excavation has arrived at that location. This pore pressure buildup occurs due the transient coupling effect. Third, due to pore pressure dissipation in the advance core, the ground continues to converge during the standstill period, even after the excavation has reached the final tunnel face.
In terms of their potential field application, the first and second results could be used to formulate the unloading factors in the convergence-confinement method for tunneling in deep saturated ground. In this situation, the unloading factors should now be nonlinear and should consists of: (1) the mechanical unloading factor in the form of excavation force and (2) the hydraulic unloading factor in the form of excavation pore pressure.
Combined with the third result, a relationship between the unloading factors and the distance to the tunnel face can be developed. This relationship is vital because it will provide the timing for tunnel support installation.

Convergence and extrusion
As described in "Introduction", the loss of radial and axial confinements at the tunnel face induces convergence and face extrusion as the excavation progresses. In fact, previous results have shown that the pore pressure continues to dissipate during the standstill period, indicating that the core ahead of the tunnel will continue to consolidate. Knowing the magnitude of the convergence and extrusion of the core during this standstill period may help determine the relationship between these types of deformation.
During the standstill period, radial and axial displacements at various locations ahead of the face continue to occur to considerable magnitudes (Fig. 9a, b), respectively), with those located at the face (δ r0 and u ex0 ) being the largest, as expected. During the standstill period (from t = 10 days to t = 930 days), the face convergence increases by 20% from 1.2 cm to 1.4 cm. The percentage increase of face extrusion with time is even larger than that of convergence. The face extrudes up to 37% by the end of the standstill period, from 3.2 to 4.4 cm (Fig. 9c). The fact that the percentage increase of face extrusion (37%) is larger than that of the convergence at the face (20%) demonstrates that the axial displacement of the tunnel wall is as important as its radial displacement. Moreover, compared to the trend at the face, the trend of pre-convergence, δ r1 and δ r2 , and the core extrusion, u ex1 and u ex2 , at some distance ahead of the face are similar but to a much lesser degree. By plotting a graph of extrusion versus convergence from Fig. 9a, b, the relationship between the convergence and the extrusion of the core can be established as shown in Fig. 10.
From Fig. 10, it can be seen that at all locations inside the core (at 0R, − 2R and − 4R), the trends of the relationships are similar in that they are all linear. The proposed linear equations are presented in the insets in Fig. 10. These equations may be used as first estimates of the magnitude of either the extrusion or the convergence when it is not possible to measure it directly. Specifically, knowing the behavior of face extrusion from convergence (Fig. 10a) may help the tunnel engineer to determine the method of face

Conclusions
This study investigated the transient H-M response when excavating a tunnel in a deep and saturated ground using a 2-D axisymmetric FLAC model. The induced short-term and long-term H-M responses of the ground were analyzed as well as their transient Fig. 10 Equations relating the linear relationship between extrusion and convergence of the advance core behaviors including the three deformation types, i.e., convergence, extrusion and preconvergence of the advance core.
The results showed that during the excavation period, the pore pressure build-up in the region ahead of the face has increased beyond its in-situ value even before the tunnel face has arrived. This so-called Mandel-Cryer effect would stiffen the advance core and provide temporary confinement to the face, increasing its shearing resistance. As the excavation passed, the excess pore pressure dissipated nonlinearly with time and the confinement disappeared.
The same behavior was also observed during the long-term drained consolidation period. After the excavation of the final face, the coupled H-M response of the ground did not stop shortly but continued with time, transitioning the ground toward its steadystate condition. Consequently, the pore pressure along the axis of the tunnel was higher in the coupled simulation than that in the uncoupled one, causing less convergence in the former than in the latter. This result can be potentially applied to improve the formulation of the unloading factors in the convergence-confinement method, as well as to develop the equation for constructing transient longitudinal displacement profiles considering the transient nature of the consolidating ground. Thus, a relationship between the transient deformation and the distance to the tunnel face could be developed. In tunnel engineering, this relationship is vital because it helps determine the timing for the installation of support.
In addition to the displacement, the ground also responded mechanically by transferring the excavation load toward the advance core during the long-term consolidation. The load transfer was indicated by the increase of the effective radial stress and the expansion of the effective tangential stress in the advance core ahead of the face. If the load transfer exceeds the strength of the core, the tunnel face will collapse, followed by the collapse of the opening. The transient deformation of the advance core also suggests that the axial displacement of the ground is equally important as its radial displacement.
To summary, the potential field application of the study is as follow: 1. When simulating a tunnel excavation in deep saturated ground using the convergence-confinement method, the unloading factors should be nonlinear and should consists of: (1) the mechanical unloading factor in the form of excavation force and (2) the hydraulic unloading factor in the form of excavation pore pressure. These two unloading factors are necessary because based on the result of this paper, the induced H-M response near the tunnel face is a rather transient response instead of an initial or final response. Moreover, it is observed that the pore pressure dissipation is not linear either with time or with distance to the tunnel face. 2. A relationship between the unloading factors and the distance to the tunnel face should then be established. This relationship is vital because it will provide the timing for tunnel support installation. 3. The extrusion and the convergence of the advance core could be related through the proposed equations capturing the linear relationships between the face extrusion and its convergence as well as between the core extrusion and its pre-convergence. Through these relationships, the tunnel engineer may be able to estimate the magnitude of the deformation ahead of the face, which will subsequently allow control of the deformation behind the face. It is for this purpose that these equations have been developed.
Even though the study is limited to 2-D problem and under limited number of simulations, it provides a new path to capturing the transient nature of tunnel excavation in deep saturated ground. Moreover, future research may use the framework that has been established in this paper to simulate transient H-M responses due to tunneling in more complex ground behaviors.