Stability Analysis of the Horseshoe Tunnel Face in Rock Masses

Accurately estimating the stability of horseshoe tunnel faces remains a challenge, especially when excavating in rock masses. This study aims to propose an analytical model to analyze the stability of the horseshoe tunnel face in rock masses. Based on discretization and “point-by-point” techniques, a rotational failure model for horseshoe tunnel faces is first proposed. Based on the proposed failure model, the upper-bound limit analysis method is then adopted to determine the limit support pressure of the tunnel face under the nonlinear Hoek–Brown failure criterion, and the calculated results are validated by comparisons with the numerical results. Finally, the effects of the rock properties on the limit support pressure and the 3D failure surface are discussed. The results show that (1) compared with the numerical simulation method, the proposed method is an efficient and accurate approach to evaluating the face stability of the horseshoe tunnel; (2) from the parametric analysis, it can be seen that the normalized limit support pressure of the tunnel face decreases with the increasing of geological strength index, GSI, Hoek–Brown coefficient, mi, and uniaxial compressive strength, σci, and with the decreasing of the disturbance coefficient of rock, Di; and (3) a larger 3D failure surface is associated with a high value of the normalized limit support pressure.


Introduction
With the development of shield technology, various tunnels, such as circular and non-circular tunnels (i.e., elliptical, rectangular, and horseshoe tunnels), are excavated in practical engineering. Compared with the circular tunnels, non-circular tunnels are more popular in practical engineering due to their higher section utilization and lower construction costs. To maintain the stability of the excavation face, the shield machine always provides a continuous supporting pressure at the tunnel face. A collapse failure will happen at the tunnel face when the supporting pressure is not sufficient to resist the movement of soil or rock towards the tunnel. Therefore, in the excavation of shield tunnels, one of the most important issues is to determine the required minimum supporting pressure to ensure the stability of the excavation face.
To solve this problem, the following three methods have been repeatedly adopted: (1) experimental tests [1,2]; (2) numerical simulations [3,4]; (3) analytical approaches [5][6][7][8][9][10]. Because, compared with experimental tests and numerical simulations, the analytical approaches are easier to implement and can provide a vehicle for engineers to directly design required minimum supporting pressures for the tunnel face [8], the analytical approaches, especially the upper-bound limit analysis approach, are the most popular methods to assess stability issues. Recently, many scholars have proposed various analytical models to investigate tunnel face stability [5][6][7][8]. For instance, considering the soil arching effect, Han et al. [5] proposed a failure model composed of five conical blocks to estimate

Nonlinear Hoek-Brown Failure Criterion
The generalized nonlinear Hoek-Brown failure criterion was first introduced by Hoek [18] to analyze the stability of engineering structures in rock masses and now has been widely adopted to analyze stability issues [18,19]. The generalized nonlinear Hoek-Brown failure criterion can be described as follows: where σ 1 and σ 3 denote the maximum and minimum principal stresses, respectively, and σ ci denotes the uniaxial compressive strength of the rock material. The values s, a and m can be obtained by Equations (2)-(4): where m i is the nonlinear Hoek-Brown coefficient; GSI is the geological strength index; and D i is the disturbance coefficient of rock. According to Equations (1)- (4), there are four input parameters for the generalized nonlinear Hoek-Brown failure criterion, namely σ ci , GSI, m i , and D i . The values of these four parameters can be obtained from Prise [19]. In order to estimate the stability of geotechnical structures using the generalized nonlinear Hoek-Brown failure criterion, Yang et al. [20] proposed a tangent technique to transform the nonlinear failure criterion into a linear one, as shown in Figure 1. From Figure 1, it can be seen that the tangent technique can provide upper-bound solutions for stability issues because the tangent line is above the nonlinear envelope. In Figure 1, the tangent line at point M can be expressed as follows: where τ and σ n denote the shear stress and normal stress at the failure surface, respectively; c t and ϕ t are the equivalent cohesion and equivalent friction angle, respectively; and c t can be calculated as follows: where i m is the nonlinear Hoek-Brown coefficient; GSI is the geological strength index; and Di is the disturbance coefficient of rock. According to Equations (1)-(4), there are four input parameters for the generalized nonlinear Hoek-Brown failure criterion, namely ci σ , GSI , i m , and Di. The values of these four parameters can be obtained from Prise [19]. In order to estimate the stability of geotechnical structures using the generalized nonlinear Hoek-Brown failure criterion, Yang et al. [20] proposed a tangent technique to transform the nonlinear failure criterion into a linear one, as shown in Figure 1. From Figure 1, it can be seen that the tangent technique can provide upper-bound solutions for stability issues because the tangent line is above the nonlinear envelope. In Figure 1, the tangent line at point M can be expressed as follows: where τ and σn denote the shear stress and normal stress at the failure surface, respectively; t c and t ϕ are the equivalent cohesion and equivalent friction angle, respectively; and t c can be calculated as follows: It is worth noting that the equivalent cohesion t c is the mutative parameter with the change in equivalent friction angle t ϕ . In this study, the equivalent friction angle t ϕ is not a given value but an additional optimization parameter, which is different from the linear Mohr-Coulomb failure criterion.  Figure 2 shows the rotational failure model for the horseshoe tunnel face. From Figure 2, it can be seen that the proposed failure model rotates around a horizontal X-axis passing through point O with a uniform angular velocity ω. As shown in Figure 2, It is worth noting that the equivalent cohesion c t is the mutative parameter with the change in equivalent friction angle ϕ t . In this study, the equivalent friction angle ϕ t is not a given value but an additional optimization parameter, which is different from the linear Mohr-Coulomb failure criterion. Figure 2 shows the rotational failure model for the horseshoe tunnel face. From Figure 2, it can be seen that the proposed failure model rotates around a horizontal x-axis passing through point O with a uniform angular velocity ω. As shown in Figure 2, the failure model proposed in this study can be divided into two parts, S1 and S2. The proposed failure model is assumed to be bounded by two discrete log-spirals as follows:

Rotational Failure Model for the Horseshoe Tunnel Face
where r A and r B are the rotation radius of points A and B, respectively; and β A and β B are rotation angles of points A and B, respectively. Based on the discretization and the "pointby-point" techniques [6], the proposed failure model can be constructed by following five steps: (1) based on the discretization technique, the horseshoe tunnel face is first uniformly discretized into n points, namely points A 1 , A 2 , A 3 . . . A n/2 and B 1 , B 2 , B 3 . . . B n/2 (see Figure 2); (2) the failure part S1 is discretized into n 1 parts by a series of planes passing through the horseshoe tunnel face and the rotation center O (see Figure 2); (3) the failure part S2 is discretized into n 2 parts by a series of planes with angle dβ between adjacent planes passing through point O; (4) based on the "point-by-point" technique, the point P i,j+1 (i = 1, 2, 3, . . . N, j = 1, 2, 3, . . . N) at the failure surface is generated from the given points P i,j and P i+1,j with the starting points A 1 , A 2 , A 3 . . . A n/2 and B 1 , B 2 , B 3 . . . B n/2 at the horseshoe tunnel face (see Figure 3); and (5) the 3D failure surface for the horseshoe tunnel face can be obtained by connecting all adjacent points, P i,j , P i+1,j and P i,j+1 , as shown in Figure 4. More details about the discretization and the "point-by-point" techniques can be found in Mollon et al. [6]. β are rotation angles of points A and B , respectively. Based on the discretization and the "point-by-point" techniques [6], the proposed failure model can be constructed by following five steps: (1) based on the discretization technique, the horseshoe tunnel face is first uniformly discretized into n points, namely points 1 2 3 /2 , , ... n A A A A and     β are rotation angles of points A and B , respectively. Based on the discretization and the "point-by-point" techniques [6], the proposed failure model can be constructed by following five steps: (1) based on the discretization technique, the horseshoe tunnel face is first uniformly discretized into n points, namely points 1 2 3 /2 , , ... n A A A A and Figure 2); (2) the failure part S1 is discretized into 1 n parts by a series of planes passing through the horseshoe tunnel face and the rotation center O (see Figure  2); (3) the failure part S2 is discretized into 2 n parts by a series of planes with angle d β between adjacent planes passing through point O ; (4) based on the "point-by-point" technique, the point , 1 ) at the failure surface is generated from the given points Pi,j and Pi+1,j with the starting points

Work Rate Calculations
Based on the proposed failure model, the limit support pressure of the tunnel face can be derived by equating the work rate exerted by external force to the internal energy dissipation rate using the upper-bound limit analysis method [21][22][23][24][25], which described

Work Rate Calculations
Based on the proposed failure model, the limit support pressure of the tunnel face can be derived by equating the work rate exerted by external force to the internal energy dissipation rate using the upper-bound limit analysis method [21][22][23][24][25], which described by Equation (8). Note that the external forces in this study mainly include the supporting pressure at the tunnel face and the rock gravity.
where . W T denotes the work rate exerted by the supporting pressure; . W D denotes the internal energy dissipation rate; . W γ denotes the work rate exerted by the rock gravity. The details of the work rate calculation will be presented in following subsections.

Internal Energy Dissipation Rate
As shown in Figure 5, the internal energy dissipation rate of the proposed failure model can be expressed by: where c t is the equivalent cohesion, which can be derived by Equation (6); R i,j is the rotation radius of the center of the discretization surface F i,j ; and S i,j is the area of the discretization surface F i,j , as shown in Figure 5.

Work Rate Exerted by the Rock Weight
The work rate of the rock weight can be calculated by: where , i j V is the volume of each discretization block of the proposed failure model and  Figure 5).

Work Rate Exerted by the Supporting Pressure
The work rate of the supporting pressure acting on the horseshoe tunnel face can be calculated by:

Work Rate Exerted by the Rock Weight
The work rate of the rock weight can be calculated by: where V i,j is the volume of each discretization block of the proposed failure model and β i,j is the rotation angle of the discretization block V i,j corresponding to the discretization surface F i,j (see Figure 5).

Work Rate Exerted by the Supporting Pressure
The work rate of the supporting pressure acting on the horseshoe tunnel face can be calculated by: where R i,0 is the rotation radius of the center of the discretization surface F i and β i,0 is the rotation angle of the center of the discretization surface F i .

Determination of the Limit Support Pressure
Combining Equation (8) up to Equation (11), the limit support pressure can be calculated by: where N γ and N c t are the nondimensional coefficients, representing the influence of the rock gravity and equivalent cohesion on the limit support pressure. With the combination of Equations (8)-(11), the expressions of N γ and N c t can be written as: Combining Equation (12), a nonlinear minimum optimization algorithm named fminsearch is adopted to calculate the limit support pressure with respect to β E , r E /H, and ϕ t in this study. To ensure a responsible failure model of the tunnel face, the optimization for the limit support pressure is under the constraint expressed in Equation (15).

Validation
Since no analytical models have been proposed to investigate the face stability of the horseshoe tunnel in rock masses, the numerical results obtained from the numerical software FLAC 3D were employed to validate the proposed method. Figure 6 shows the established numerical model based on the FLAC 3D . As a symmetrical tunnel is considered, the assessment of the tunnel face stability using the numerical approach is only based on half of the horseshoe tunnel (see Figure 6). It is worth noting that the numerical model is large enough (the length, width, and height of the numerical model, respectively, equaling 130 m, 100 m, and 60 m) to eliminate the effect of boundary conditions. As the rock deformation occurs mainly near the tunnel face, the mesh in the vicinity of the tunnel face is densified in order to improve the calculation efficiency (see Figure 6). As also shown in Figure 6, the boundary conditions of the numerical model are as follows [3,4]: the bottom of the numerical model is full-fixed, the vertical faces of the numerical model are fixed in the normal directions, and the top of the numerical model is free to displace. An elastic-plastic constitutive model following the nonlinear The Hoek-Brown failure criterion is adopted to simulate the rock masses, and the shell structure is used to simulate the tunnel lining element with the thickness being 22 mm, Poisson's ratio being 0.2, and Young's modulus being 15 Gpa. The rock properties adopted in the comparisons are as follows: γ = 25 kN/m 3 , m i = 5 − 35, GSI = 20 − 30, σ ci = 1 MPa, and D i = 0.5 [19]. The Young's Modulus E m and passion ratio V m of the rock mass can be obtained by [26,27]: (16) merical model are fixed in the normal directions, and the top of the numerical model is free to displace. An elastic-plastic constitutive model following the nonlinear The Hoek-Brown failure criterion is adopted to simulate the rock masses, and the shell structure is used to simulate the tunnel lining element with the thickness being 22 mm, Poisson's ratio being 0.2, and Young's modulus being 15 Gpa. The rock properties adopted in the comparisons are as follows:  Figure 6. Numerical model. Figure 6. Numerical model.
Based on the established numerical model, by gradually reducing the supporting pressure until the tunnel face is in a critical state, the corresponding supporting pressure can be recognized as the limit support pressure σ T [28]. A comparison between the limit support pressure obtained by the proposed method and the results from the numerical model is presented in Figure 7. As shown in Figure 7, the value of limit support pressure obtained by the proposed method greatly decreases with the increasing of m i , which is consistent with the numerical results. Additionally, the limit support pressures obtained in this study are slightly lower than the numerical results, with a maximum difference less than 8%. This observation illustrates that the proposed method can be used to accurately estimate the stability of the horseshoe tunnel face under the nonlinear Hoek-Brown failur criterion. Based on the established numerical model, by gradually reducing the supporting pressure until the tunnel face is in a critical state, the corresponding supporting pressure can be recognized as the limit support pressure σT [28]. A comparison between the limit support pressure obtained by the proposed method and the results from the numerical model is presented in Figure 7. As shown in Figure 7, the value of limit support pressure obtained by the proposed method greatly decreases with the increasing of i m , which is consistent with the numerical results. Additionally, the limit support pressures obtained in this study are slightly lower than the numerical results, with a maximum difference less than 8%. This observation illustrates that the proposed method can be used to accurately estimate the stability of the horseshoe tunnel face under the nonlinear Hoek-Brown failure criterion. Furthermore, the numerical simulation method is time-consuming and usually requires some redundant input variables, some of which are difficult to specify [28]. Compared with the numerical simulation method, the proposed method is easier to implement and can provide a vehicle for engineers to directly design required limit support pressures for the horseshoe tunnel face. In summary, one can conclude that the proposed method is an efficient and accurate approach to determining the face stability of the horseshoe tunnel in rock masses.

Effect of the Geological Strength Index, GSI , on the Tunnel Face Stability
In this section, the effect of the geological strength index, GSI , on the limit support pressure of the horseshoe tunnel face is discussed (see Figure 8). The curve of normalized limit support pressure, σΤ/γΗ, with GSI is presented in Figure 8. The adopted param- 10 60 GSI = − . It can be seen in Figure 8 that the normalized limit support pressure, σT/γH, decreases with the geological strength index, GSI , increases, and the decreasing rate Furthermore, the numerical simulation method is time-consuming and usually requires some redundant input variables, some of which are difficult to specify [28]. Compared with the numerical simulation method, the proposed method is easier to implement and can provide a vehicle for engineers to directly design required limit support pressures for the horseshoe tunnel face. In summary, one can conclude that the proposed method is an efficient and accurate approach to determining the face stability of the horseshoe tunnel in rock masses.

Effect of the Geological Strength Index, GSI, on the Tunnel Face Stability
In this section, the effect of the geological strength index, GSI, on the limit support pressure of the horseshoe tunnel face is discussed (see Figure 8). The curve of normalized limit support pressure, σ T / γH, with GSI is presented in Figure 8. The adopted parameters are L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , σ ci = 1 MPa, D i = 0.5, m i = 5-35, and GSI = 10-60. It can be seen in Figure 8 that the normalized limit support pressure, σ T /γH, decreases with the geological strength index, GSI, increases, and the decreasing rate decreases with the increasing of GSI. For instance, for m i = 10, the normalized limit support pressure, σ T /γH, decreases from 0.35 to 0.03 when GSI increases from 10 to 50. This indicates that a lower limit support pressure is associated with a high value of GSI, and a high GSI is beneficial for the tunnel face stability.  Figure 9 that the 3D failure surface is large for high values of GSI. Combining this with Figure 8, one can conclude that a larger 3D failure surface is associated with a high value of the normalized limit support pressure.     Figure 9 that the 3D failure surface is large for high values of GSI. Combining this with Figure 8, one can conclude that a larger 3D failure surface is associated with a high value of the normalized limit support pressure. Figure 10 shows the curve of the normalized limit support pressure σ T /γH with the Hoek-Brown coefficient m i . The analytical parameters are L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , σ ci = 1 MPa, D i = 0.5, GSI = 10-50, and m i = 5-35. Figure 11 presents the effect of the Hoek-Brown coefficient, m i , on the 3D failure surface for L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , σ ci = 1 MPa, D i = 0.5, GSI = 10, and m i = 10-35. As shown in Figure 10, the normalized limit support pressure, σ T /γH, decreases with the Hoek-Brown coefficient, m i , increases, and the decreasing rate decreases significantly with the increasing of m i . For example, for GSI = 10, the normalized limit support pressure σ T /γH decreases from 0.35 to 0.11 when m i increases from 10 to 35. This observation can be illustrated by the change of the 3D failure surface from a larger size to a smaller one with the increasing of the Hoek-Brown coefficient, m i , from 10 to 35 (see Figure 11).    Figure 11). Figure 10. Effect of i m on the normalized limit support pressure σT/γH. Figure 11. Effect of i m on the 3D failure surface. Figure 10. Effect of m i on the normalized limit support pressure σ T /γH.

Effect of the Uniaxial Compressive Strength, σ ci , on the Tunnel Face Stability
The effect of the uniaxial compressive strength, σ ci , on the normalized limit support pressure, σ T /γH, is discussed in this section, and the curve of the normalized limit support pressure, σ T /γH, with different uniaxial compressive strength, σ ci , values is presented in Figure 12. The adopted parameters are as follows: L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , GSI = 10, m i = 5, D i = 0-1.0, and σ ci = 1-10 MPa. Figure 13 presents the effect of the uniaxial compressive strengths, σ ci , on the 3D failure surface for L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , GSI = 10, m i = 5, D i = 0.4, and σ ci = 2-10 MPa. From Figure 10, it can be seen that the normalized limit support pressure, σ T /γH, decreases with the σ ci increases. For instance, for D i = 0.4, the normalized limit support pressure, σ T /γH, decreases from 0.33 to 0.11 when σ ci increases from 2 MPa to 10 MPa. This is consistent with Senent et al. [29] for the stability analysis of the circular tunnel face. This observation is illustrated by the changes of the 3D failure surface from a larger size to a smaller one with the increasing of the uniaxial compressive strength, σ ci , from 2 MPa to 10 MPa (see Figure 13).

Effect of the Uniaxial Compressive Strength, ci σ , on the Tunnel Face Stability
The effect of the uniaxial compressive strength, ci σ , on the normalized limit support pressure, σT/γH, is discussed in this section, and the curve of the normalized limit  σ , values is presented in Figure 12. The adopted parameters are as follows: 11.6m L = , the normalized limit support pressure, σT/γH, decreases from 0.33 to 0.11 when ci σ increases from 2 MPa to 10 MPa. This is consistent with Senent et al. [29] for the stability analysis of the circular tunnel face. This observation is illustrated by the changes of the 3D failure surface from a larger size to a smaller one with the increasing of the uniaxial compressive strength, σci, from 2 MPa to 10 MPa (see Figure 13).   This section discusses the effect of the disturbance coefficient of rock, D i , on the normalized limit support pressure, σ T /γH, and the curve of the normalized limit support pressure, σ T /γH, with different rock disturbance coefficients, D i , is presented in Figure 14. The analytical parameters are as follows: L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , m i = 5, GSI = 10, σ ci = 1-10 MPa, and D i = 0-1.0. Figure 15 presents the effect of the disturbance coefficient of rock, D i , on the 3D failure surface for L = 11.6 m, H = 10.3 m, γ = 25 kN/m 3 , m i = 5, GSI = 10, σ ci = 10 MPa, and D i = 0.2-0.8. It can be seen that the normalized limit support pressure, σ T /γH, increases with the increasing of the rock disturbance coefficient, D i , and the increasing rate of the normalized limit support pressure, σ T /γH, increases significantly with the rock disturbance coefficient, D i , increases. For instance, for σ ci = 10 MPa, the normalized limit support pressure, σ T /γH, increases from 0.07 to 0.40 when D i increases from 0.2 to 0.8. This observation can be illustrated by the changes of the 3D failure surface from a smaller size to a larger one with the increasing of the rock disturbance coefficient, D i , from 0.2 to 0.8 (see Figure 15).

Effect of the Disturbance Coefficient of Rock, i D , on the Tunnel Face Stability
This section discusses the effect of the disturbance coefficient of rock, i D , on the normalized limit support pressure, σT/γH, and the curve of the normalized limit support pressure, σT/γH, with different rock disturbance coefficients, i D , is presented in Figure   14. The analytical parameters are as follows: 11.6m L = . It can be seen that the normalized limit support pressure, σT/γH, increases with the increasing of the rock disturbance coefficient, i D , and the increasing rate of the normalized limit support pressure, σT/γH, increases significantly with the rock disturbance coefficient, i D , increases.
For instance, for 10MPa ci σ = , the normalized limit support pressure, σT/γH, increases from 0.07 to 0.40 when i D increases from 0.2 to 0.8. This observation can be illustrated by the changes of the 3D failure surface from a smaller size to a larger one with the increasing of the rock disturbance coefficient, i D , from 0.2 to 0.8 (see Figure 15). Figure 14. Effect of i D on normalized limit support pressure, σT/γH. Figure 14. Effect of D i on normalized limit support pressure, σ T /γH.

Effect of the Disturbance Coefficient of Rock, i D , on the Tunnel Face Stability
This section discusses the effect of the disturbance coefficient of rock, i D , on the normalized limit support pressure, σT/γH, and the curve of the normalized limit support pressure, σT/γH, with different rock disturbance coefficients, i D , is presented in Figure   14. The analytical parameters are as follows: 11.6m L = . It can be seen that the normalized limit support pressure, σT/γH, increases with the increasing of the rock disturbance coefficient, i D , and the increasing rate of the normalized limit support pressure, σT/γH, increases significantly with the rock disturbance coefficient, i D , increases.
For instance, for 10MPa ci σ = , the normalized limit support pressure, σT/γH, increases from 0.07 to 0.40 when i D increases from 0.2 to 0.8. This observation can be illustrated by the changes of the 3D failure surface from a smaller size to a larger one with the increasing of the rock disturbance coefficient, i D , from 0.2 to 0.8 (see Figure 15).

Conclusions
In the framework of the upper-bound limit analysis theory, this study analyzes the stability of horseshoe tunnel faces in rock masses. Compared with previous analytical solutions, the following improvements have been achieved: (1) In this study, the nonlinear Hoek-Brown failure criterion is first incorporated into the stability analysis of horseshoe tunnel faces in rock masses. The comparisons between the results from the proposed method with the numerical results illustrate that the proposed method is an efficient and accurate approach to assessing the face stability of horseshoe tunnels. (2) Based on the proposed method, the effect of rock properties on the normalized limit support pressure and the 3D failure surface are presented. It is shown that, for selected cases, the normalized limit support pressure of the tunnel face greatly decreases with the increasing of GSI, m i , and σ ci , and decreasing of D i ; the 3D failure surface become larger with the decreasing of GSI, m i , and σ ci , and increasing of D i ; a larger 3D failure surface is associated with a high value of the normalized limit support pressure; and high GSI, m i , and σ ci are beneficial for the tunnel face stability.
A possible extension of this work could be further investigation into the active and passive stability of deep/shallow horseshoe tunnel faces [30] considering the seismic loads [31] and the seepage forces [32].