Abstract
Natural soils are mostly in an unsaturated state and the corresponding mechanical properties differ significantly from a saturated one. In traditional stability analysis for a shield tunnel face, the soil mass is typically assumed to be dry or saturated for convenience of analysis. In this work, based on the upper bound theorem of classical plasticity theory and the log-spiral failure mechanism, face stability of the shield tunnel excavated in unsaturated soils under vertical steady flow is studied. The profile of shear strength is determined by virtues of the unified effective stress approach and the analytical solution of matric suction under unsaturated steady flow. On this basis, the analytical expression for the supporting pressure is deduced and the sequential quadratic programming is employed to search for the optimal upper bound solution. Through parametric analysis it is found the fitting parameters of SWCC and the groundwater level affect the supporting pressure notably. Besides, the relevant influence of vertical discharge is also significant when clayey soils are concerned.
1. Introduction
With the gradual maturity of shield tunnelling method, it has been widely applied in the construction of infrastructural facilities, especially in urban areas. One of the key problems that engineers concern during shield construction is the stability of tunnel face. As a result of excavation, stress release takes place and the original equilibrium of soil mass is disturbed. In order to prevent the possible failure of tunnel face, a supporting pressure is usually applied by the shield machine in various forms like bentonite slurry, earth pressure, and compressed air [1]. The magnitude of this supporting pressure should be reasonably determined, since a small pressure will lead to an active failure or a collapse, while an overlarge pressure may give rise to a passive failure or a blowout [2, 3].
In the literature, the stability of tunnel face has been extensively studied and various methods have been utilized, including experimental tests [4–6], numerical simulation [7, 8], and limit equilibrium method [9–11]. These approaches are all effective means to address the stability issue of tunnel face and they possess distinct advantages over others. For the limit equilibrium method, its simple principle and easy computation make it popular in the practical application. However, the result obtained from a limit equilibrium analysis is not rigorous since assumptions with respect to forces are required. An alternative approach is the limit analysis method which is based upon the upper and lower bound theorems of classical plasticity [12]. Apart from the merits of the limit equilibrium method, this method also encompasses a rigorous theory basis. In recent years, the limit analysis method has been widely employed in the corresponding researches on tunnel face stability. Tang et al. [13] extended the three-dimensional failure mechanism of tunnel face proposed by Leca and Dormieux [14] and investigated the face stability in layered soils. In view of deep rock tunnels, Zhang et al. [15] constructed a log-spiral failure mechanism and performed the corresponding upper bound analysis where the Hoek-Brown failure criterion was employed. Based on the same failure mechanism, Liang et al. [16] studied the stability of a shield tunnel face excavated in nonhomogeneous and anisotropic soils. Mollon et al. [17, 18] proposed a spatial discretization technique for generating the failure surface point by point, and the solutions of collapse pressure were significantly improved with the aid of the constructed failure mechanism. Due to the advantages of this technique, it was later utilized to study tunnel face stability in more complex settings [1, 19]. Besides, in the stability analysis of tunnel face, some scholars also combined the limit analysis theory with other approaches, including the finite element method [20, 21] and the reliability theory [22, 23].
It is noted that only dry or saturated condition is discussed in the existing research on tunnel face stability. However, natural soils are mostly unsaturated and the dry or saturated condition is a particular case [24, 25]. Due to the existence of matric suction, mechanical properties of unsaturated soils vary greatly from those of saturated soils. As a matter of fact, stability of slopes and bearing capacity of foundations have been extensively investigated in the context of unsaturated soils [26–30], and relevant research findings show the impact of unsaturated soil characteristics is notable. Therefore, in order to ensure the safe construction of shield tunnels, it is essential to understand the corresponding face stability in unsaturated soils.
In this work, upper bound analysis is conducted for the face stability assessment of shield tunnels excavated in unsaturated soils. The vertical unsaturated steady flow is considered and three different flow conditions are discussed, including evaporation, no-flow, and infiltration. Based on Darcy’s law and Gardner’s model, the matric suction profile along the vertical direction can be obtained [26, 27]. After that, by combining the unified effective stress approach and the closed-form equations for suction stress established by Lu et al. [31], the apparent cohesion profile along the depth is determined, which is later incorporated into the theoretical derivation for analytical solution of the required supporting pressure. For specified strength parameters, the corresponding upper bound solution is searched through optimization. Parametric analysis is finally presented and the obtained conclusions can provide a better guidance to the design and construction of shield tunnels.
2. Background Theory
2.1. Upper Bound Theorem of Limit Analysis
The limit analysis method has become an effective approach to address stability problems in soil mechanics since it was introduced by Chen [12]. The limit theorems on which limit analysis is based is established in the light of assumptions of perfect plasticity and associated flow rule. In contrast to the lower bound theorem, the limit analysis based upon upper bound theorem is mostly performed for stability evaluation of geotechnical structures. The upper bound theorem states that, for satisfied velocity boundary conditions, the actual collapse load is no larger than the load derived by equating the rate of work done by the external forces to the rate of internal energy dissipation [12]. This theorem can be mathematically formulated aswhere and are tensors of stress and strain rate, respectively, denotes the surcharge load on the surface , is the velocity along the velocity discontinuity surface, and is the body force within the volume .
2.2. Unified Effective Stress Approach
In Bishop’s effective stress model, the suction stress is zero when the soil mass is dry. However, both experimental studies and field observations show that the realistic suction stress with respect to dry condition is non-zero for clayey and some silty soils. This phenomenon is mainly attributed to the existence of physicochemical interparticle forces which is ignored in the Bishop’s model. In order to account for this effect in the existing framework, Lu and Likos [25] extended the Bishop’s effective stress formulation by introducing the concept of suction stress. For both saturated and unsaturated conditions, the unified effective stress is expressed as follows:where is the total stress, denotes the pore air pressure, and is the suction stress. By employing the two fitting parameters of SWCC, Lu et al. [31] established the following closed-form equations for the suction stress .where is the pore water pressure and and are fitting parameters of SWCC. More specifically, approximately represents the inverse of the air-entry value of unsaturated soil and relates to the corresponding pore size distribution. According to Lu and Likos [25], the typical value ranges of and with respect to different soils are illustrated in Table 1.
3. Suction Stress and Apparent Cohesion under Vertical Steady Flow
Unlike saturated soils whose permeability coefficient is constant, the hydraulic conductivity of unsaturated soils highly depends on the matric suction. Based on the Gardner’s model for this characteristic dependency and the Darcy’s law for steady state flow, Lu and Likos [25] deduced the analytical solution for matric suction, which is expressed as follows:where is the vertical discharge, is the saturated permeability coefficient, is the unit weight of water, and denotes the vertical distance above the groundwater level. Besides, the upward-positive sign convention is applied for parameter ; that is, evaporation is represented by a positive value of whilst infiltration is represented by a negative value. Combining (3b) and (4), the suction stress of unified effective stress approach under condition of vertical steady unsaturated flow yields
With specified parameters the suction stress at any distance above the groundwater line can be derived according to the above equation. In Figure 1, profiles of suction stress in four representative soils under steady flow conditions are shown. The detailed parameters of these four types of soil mass are given in the corresponding pictures. It can be observed that the distribution of suction stress is significantly influenced by the soil types or, more specifically, by the parameters and . In general, with the increase in and , the suction stress above the groundwater level gradually decreases and the impact of flow rate also tends to disappear. For clay, the suction stress keeps increasing, whereas there exists an obvious turning point for the other three types of soils.

(a) Clay

(b) Silt

(c) Loess

(d) Sand
The existence of suction stress provides additional cohesion to the soil mass and it is termed as apparent cohesion. In the frame work of unified effective stress approach, the apparent cohesion is defined as follows [26, 27]:where is the effective friction angle of soils. As a result, the total cohesion of unsaturated soil is the sum of effective cohesion and apparent cohesion , namely, . Due to the suction stress above the groundwater level in clay and silt being larger than in loess and sand, as shown in Figure 1, the corresponding total cohesion also differs significantly from the effective cohesion. Therefore, for tunnel construction in clay and silt, it is particularly essential to take mechanical properties of unsaturated soils and unsaturated steady flow into consideration.
4. Log-Spiral Failure Mechanism of Shield Tunnel Face
Due to the disturbance induced by tunnel excavation, the original equilibrium of the soil mass in front of tunnel face is destroyed. When there is not enough supporting pressure, the soil mass tends to move into the tunnel and collapse accidents may occur. For dry or saturated soils, the log-spiral failure mechanism can be employed to depict this failure of tunnel face [15, 16], where the associated flow rule is satisfied on the slip surface. This mechanism is also applicable to face stability analysis in unsaturated soils, considering the fact that the degree of saturation has no impact on their internal friction angle [25]. As shown in Figure 2, the rigid block , which is bounded by the two log-spiral curves and and the tunnel face, rotates about point at the incipience of collapse. In terms of the polar coordinate system, which is centered at point and starts from the vertical direction, distances and for any point on and to point can be expressed aswhere and are lengths of and , respectively, and are corresponding angles with the vertical direction. By designating as diameter of the tunnel and as angle of and the vertical direction, the following equations can be derived according to geometrical relationship.

5. Upper Bound Analysis for Required Supporting Pressure
In this section, analytical solution for the required supporting pressure is deduced based on the upper bound theorem of limit analysis and the log-spiral failure mechanism. The external forces considered herein include the gravity and the supporting pressure. The rate of work done by gravity is calculated by means of the superposition method, where the corresponding work rate is expressed as the sum of four subareas of , , , and . For the subarea , the rate of work of gravity yieldswhere is the unit weight of soil mass, denotes the velocity of a point in area , is the angular velocity of collapsing block, and is a function related to parameters and and its detailed expression is given in Appendix. For the triangular subarea , the rate of work of gravity isAnalogously, the rates of work of subareas and , respectively, are
The overall work rate of gravity can then be derived based on the solutions of these four subareas and the following relation:As for supporting pressure , it is typically regarded as a uniform pressure applied on the tunnel face, as shown in Figure 2. The rate of work of in the revolving velocity field is deduced as follows:where the minus sign indicates that the negative work is done by the gravity.
Due to the fact that the collapsing block is assumed to be rigid, the internal energy only dissipates along the slip surfaces and . Different from the saturated condition, the total cohesion of unsaturated soils under vertical steady flow is height-dependent. To determine the rate of internal energy dissipation, the vertical distance of a point on log-spiral curve above the groundwater level is obtained aswhere is the vertical distance of groundwater level and tunnel invert. Combining (18) and the apparent cohesion expressed in (6), the total cohesion of soil mass at this point, , can be derived. Consequently, the rate of internal energy dissipation along the log-spiral curve is deduced asSimilarly, the vertical distance of a point on the log-spiral curve and the groundwater level isBy designating as the total cohesion at this point, the corresponding rate of internal energy dissipation yieldsAs a result, the analytical solution for the required supporting pressure is deduced by equating the rate of work of external forces to the rate of internal energy dissipation, which is expressed as
For specified and and given values of other parameters, a numerical solution of supporting pressure can be calculated according to (22). The optimal upper bound solution of corresponds to the maximum of these values, and thus the following mathematical optimization model is established. In this paper, the sequential quadratic programming algorithm is employed to determine the required supporting pressure.
6. Numerical Results and Discussions
6.1. Validation of the Log-Spiral Failure Mechanism
Although the stability of shield tunnel face is in effect a three-dimensional issue, this paper treats it as plane strain for the convenience of incorporating the strength features of unsaturated soils. Prior to investigating the influence of diverse parameters, it is essential to examine the validity of this mechanism and to discuss the resulting difference to actual condition. In view of tunnel face stability in dry soils, Chambon and Corté [32] conducted a series of centrifugal-model tests. The strength parameters of experimental materials are as follows: unit weight ; internal angle of friction ; and cohesion . The variation of values of shear strength parameters is attributed to some uncertainty in the measurement. Corresponding to four combinations of strength parameters, the collapse pressure is solved in the light of the log-spiral failure mechanism, and its comparison with experimental result is shown in Table 2. It is found that, for the four combinations of strength parameters, the analytical solutions are larger than the experimental results, albeit to different extents. This is because the plane strain analysis is more conservative to a three-dimensional one as the beneficial effect of intermediate principal stress is ignored. Besides, as regards plane strain analysis of tunnel face stability, some other failure mechanisms are also available in the literature, e.g., the mechanism based on a discrete technique proposed by Mollon et al. [17] and the joined mechanism of translation and rotation proposed by Xu et al. [33]. Corresponding to a tunnel face in dry sand () with and , the collapse pressures with respect to various friction angles are derived in the light of these three mechanisms, as shown in Table 3. It can be found that the results are on the whole close. In fact, the discrete mechanism degenerates to the log-spiral mechanism when the internal friction angle is uniform and the discretization parameter is reasonably small. Therefore, the minor deviation between discrete and log-spiral mechanisms is attributed to . On the other hand, the results of log-spiral mechanism are slightly larger than those of joined mechanism. This is due to the log-spiral mechanism underestimates arching, as was pointed out by Chambon and Corté [32]. Despite this minor defect, the log-spiral mechanism is easy to use and can yield reasonable results, and hence it is suitable for the analysis of effect of unsaturated soil strength.
6.2. Effect of Fitting Parameters of SWCC
The fitting parameters and of SWCC are important parameters affecting the mechanical properties of unsaturated soils. The physical significances of these two parameters are as follows: approximates the inverse of the air-entry value of soil and reflects the corresponding pore size distribution. As shown in Figure 1, the profile of suction stress varies notably when the values of and change, thereby indicating a significant variation in the total cohesion . In order to investigate their influence on the required supporting pressure, the parameters are chosen as follows: , , , , and . The variation of with and is presented in Figure 3 in the context of different groundwater levels of , , and . In general, the discrepancy of the curves enlarges when increases. This can be interpreted as follows: on one hand, in soil mass with large values of and (e.g., sand shown in Figure 1(d)), the suction stress quickly diminishes with the increase in vertical distance above the water table; on the other hand, in soil mass with small values of and (e.g., clay shown in Figure 1(a)), the suction stress keeps increasing with the height. Consequently, these two opposite effects are magnified when the vertical distance above the ground water level increases. For specified , the required supporting pressure generally increases with the increase in and . An interesting case is that, for relatively large values of and , the required supporting pressure versus is nearly a horizontal line, as shown in Figure 3(d). This is also because the influence of matric suction mostly vanishes for the given parameters and the resulted solution approximates the one in dry condition. It can be concluded that the existence of matric suction provides a beneficial effect on face stability, especially for tunnels constructed in clayey and silty soils (i.e., soil mass with small values of and ). In addition, particular attention should be devoted to the possible variation of groundwater level, as the mechanical performance of unsaturated soils vary significantly when groundwater level changes.

(a) z0=0m

(b) z0=5m

(c) z0=10m

(d) z0=20m
6.3. Effect of Tunnel Diameter and Effective Angle of Internal Friction
The tunnel diameter and effective angle of internal friction are also key parameters affecting the face stability. In order to study their influence on the required supporting pressure, the parameters are chosen as follows: , , , , and . The variation of with and is shown in Figure 4 for two representative soils of clay ( and ) and silt ( and 0). It can be observed that the variation patterns are similar in these two types of soils, except for the magnitude of in clay that is much smaller than that in silt. This is attributed to the fitting parameters and of clay being smaller than those of silt, and thus the beneficial effect of matric suction is more obvious in clay. For specified soil type, increases with the increase in tunnel diameter , whilst it decreases with the increase in the effective frictional angle .

(a) Clay

(b) Silt
6.4. Effect of Vertical Discharge
As illustrated in Figure 1, the vertical discharge may exert a great impact on the distribution and magnitude of suction stress. In order to study its influence on the required supporting pressure, the parameters are chosen as follows: , , , , and . The variation of with different vertical discharges of , and , which, respectively, stand for infiltration, no-flow, and evaporation, is shown in Figure 5 for clay and silt. The values of and of these two types of soils are the same as those in the previous discussion, and the relevant saturated permeability coefficients are and , respectively. Besides, for comparison, the required supporting pressure with respect to the dry condition is also given, which is referred to as the case of no suction. It can be found that the solution of under dry condition is largest among all conditions, as there is no additional cohesion being considered. For clayey soil, the influence of vertical discharge is obvious: the infiltration increases from the no-flow condition, whereas the evaporation leads to the decrease in . This is a consequence of the corresponding influence of vertical discharge on the suction stress, as shown in Figure 1(a). However, the variation of with respect to is not evident for silt. It can be concluded that the face stability analysis based upon dry or saturated properties of soil mass is conservative. Besides, for tunnel excavation in clay, the case of infiltration (i.e., those induced by rainfall) should be particularly cared as the required supporting pressure increases notably from the no-flow condition.

(a) Clay

(b) Silt
7. Conclusions
The stability of shield tunnel face is studied in the light of the upper bound theorem of limit analysis. Different from existing research based upon dry/saturated soil strength, the shield tunnel is considered to be excavated in unsaturated soils. The unified suction stress approach is employed for shear strength representation and the existing analytical solution for profile of matric suction is used to consider the vertical steady flow. These formulae are then incorporated into the upper bound formulation and the analytical solution of the required supporting pressure of tunnel face is deduced. The corresponding upper bound solution is searched through optimization and the parametric analysis is conducted to discuss the influence of diverse parameters on face stability. The fitting parameters and groundwater table are found to have great impact on the supporting pressure. As regards the vertical discharges, they affect the supporting pressure significantly for clay, whilst their influence is negligible for other types of soil.
Appendix
The detailed expressions of are given as follows:
Data Availability
No data were used to support this study.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The preparation of the paper has received financial supports from National Natural Science Foundation of China (51804113, 51674115, 51434006, and 51874130). The financial supports are greatly appreciated.