Junhao Zhong, Xiaoli Yang*

School of Civil Engineering, Central South University, Changsha 410075, China

Keywords:Face stability Piecewise linear method Hoek-Brown strength criterion Multi-horn rotational mechanism Limit analysis

A B S T R A C T The face stability problem is a major concern for tunnels excavated in rock masses governed by the Hoek-Brown strength criterion.To provide an accurate prediction for the theoretical solution of the critical face pressure, this study adopts the piecewise linear method (PLM) to account for the nonlinearity of the strength envelope and proposes a new multi-horn rotational mechanism based on the Hoek-Brown strength criterion and the associative flow rule.The analytical solution of critical support pressure is derived from the energy-work balance equation in the framework of the plastic limit theorem; it is formulated as a multivariable nonlinear optimization problem relying on 2m dependent variables(m is the number of segments).Meanwhile,two classic linearized measures, the generalized tangential technique(GTT) and equivalent Mohr-Coulomb parameters method (EMM), are incorporated into the analysis for comparison.Surprisingly, the parametric study indicates a significant improvement in support pressure by up to 13%compared with the GTT,and as expected,the stability of the tunnel face is greatly influenced by the rock strength parameters.The stress distribution on the rupture surface is calculated to gain an intuitive understanding of the failure at the limit state.Although the limit analysis is incapable of calculating the true stress distribution in rock masses, a rough approximation of the stress vector on the rupture surface is permitted.In the end, sets of normalized face pressure are provided in the form of charts for a quick assessment of face stability in rock masses.

1.Introduction

After long-term development in the mountain tunnels and underground spaces of cities, both construction technology and equipment have reached new heights.As efficient and reliable equipment, the tunnel boring machine (TBM) currently serves an increasing number of tunnel engineering projects and plays an important role in most major projects.The employment of advanced equipment greatly reduces the risk,and an accurate prediction of the required support pressure retaining on the tunnel face becomes the main concern.Regarding the TBM-excavated tunnel,Kasper and Meschke[1]conducted a finite element simulation to evaluate the shield-driven tunnel involved in each excavation stage.The influence of underground water, grouting pressure at the shield tail and the interaction of the soil-slurry system on the face stability are considered.Kim and Tonon[2]carried out a finite element analysis of a slurry-shield-driven tunnel based on the assumption of an ideal membrane at the slurry-soil contact face and obtained the minimum face pressure from the relationship between the face pressure and the corresponding displacement.However, such a numerical model simplified several aspects of the TBM parameter, such as cutterhead inclination, opening ratio and gravity of the machine,and usually overestimated the required face pressure compared with the analytical solution.Ma et al.[3]focused on the stress problems resulting from the stress release of the surrounding rock or anisotropic stress and investigated the impact of different confining stresses on the TBM performance by conducting a full-scale linear cutting test in intact granite.This work provides important guiding significance for subsequent TBM-excavated tunnels under complicated confining stress conditions.Paltrinieri et al.[4]analyzed the performance of TBMs working in highly jointed rock masses and fault zones based on the TBM performance database and established the relationships between TBM performance and various advanced conditions including rock strength and disturbance degree.Tang et al.[5]accessed the stability of tunnel faces excavated by TBM with an advanced microseismic monitoring system by which the relationship between the number and energy of microseismic events and the excavation action is established and by which the detailed spatial distribution and occurrence time of microseismic events can be determined.Such a system seems promising in controlling face stability.A more detailed introduction to TBM-excavated tunnels under adverse geological conditions is provided in an overview reported by Gong et al.[6].

In theoretical analysis,numerous contributions concerning face stability are made by the limit equilibrium method [7–10] and limit analysis method [11–15].The referred works by the limit equilibrium method were based on the same wedge-prism failure mechanism and a critical support pressure is derived based on the equilibrium between the wedge and the prism.However, in practice, the actual failure shape differs from the simple standard geometry which is so simplified that it cannot describe the real failure region.Although the implementation of the limit equilibrium method requires a presupposition of the failure mechanism similar to the upper-bound theorem and assumes a possible stress distribution in the failure surface, the solution obtained by the limit equilibrium method is neither a rigorous upper bound nor a lower bound of the true solution,because it does not strictly satisfy the definitions of the upper bound theorem and lower bound theorem.Therefore,the investigation in this study is carried out in the framework of limit analysis.

The aforementioned theoretical investigations by limit analysis are, however, limited to soil materials governed by the conventional linear Mohr-Coulomb(M-C)failure criterion.For the nonlinear failure criterion, such as the power-law criterion and Hoek-Brown strength criterion, relevant efforts are rarely involved.The Hoek-Brown strength criterion, after several rounds of improvement,has become a more generalized criterion in rock engineering and is widely accepted by scholars and engineers [16].However,this nonlinear criterion is difficult to directly introduce into the analysis in the form of principal stress, and an appropriate transformation is necessary.Thus, the problem of obtaining an analytical form for the Hoek-Brown strength envelope in τ-σnspace has concerned many researchers(τ is the shear stress;σnis the normal stress).Early,Ucar[17]attempted to develop an analytical expression for determining the shear strength from a specific case of the Hoek-Brown strength,namely a rock strength constant a=1/2.Such work was developed by Kumar [18], who extended the previous solution to the generalized Hoek-Brown strength criterion.Recently, Shen et al.[19] rearranged the formulas of Kumar and substituted the nonlinear function (σn/σc)1-afor a linear function(σcis the uniaxial compression strength).An analytical expression of shear strength parameters is captured from the strength envelope,but this procedure is solely applicable to specific ranges of the material properties.An alternative method referred to as the equivalent Mohr-Coulomb parameters method (EMM), performed early by Hoek [20], is widely accepted by scholars and engineers.The corresponding ctand δtare derived by the regression process that balances the areas above and below the equivalent M-C envelope,where ctand δtare the equivalent cohesion and the equivalent internal friction angle, respectively.Apart from the aforementioned methods, a numerical procedure originally developed by Zhang and Chen [21] is extensively applied to address the problems of geotechnical engineering with the nonlinear M-C criterion.This linearization measure is referred to as GTT and facilitates subsequent developments of the upper-bound analysis smoothly.Similar to the procedure by Zhang and Chen [21], the equivalent M-C parameters of the Hoek-Brown strength criterion are derived by Yang and Yin [22].Although a rigorous upper-bound solution is captured, the substitution of the nonlinear envelope for a single tangential line seems too rough.

In terms of the Hoek-Brown strength criterion, many works have contributed to the investigation of geotechnical engineering stability problems [23–29].Interestingly, only a few contributions are devoted to the analysis of face stability in the Hoek-Brown media [30–33].Therefore, the face stability of tunnels excavated in weakened surrounding rock remains a major concern.

Therefore, in this work, the face stability problem of the TBMexcavated tunnel in heavily fractured and low-quality rock masses governed by the Hoek-Brown strength criterion is revisited.Based on the parametric form of the Hoek-Brown strength criterion, the Hoek-Brown envelope is substituted by n linear segments, and an improved multi-horn rotational mechanism determined by 2m dependent variables is generated in light of the normality of the plastic flow rule.The critical face pressure is derived from the energy-work balance equation and quickly optimized by an original marine predators algorithm (MPA) developed by Faramarzi et al.[34].Two more linearized measures (EMM and GTT) for the Hoek-Brown strength criterion are introduced to capture the equivalent M-C parameters and the corresponding results calculated by the EMM and GTT are compared with the proposed method.

2.Methodology

2.1.Hoek-Brown strength criterion

To characterize the strength of geomaterials, several proposals of strength criteria have been developed and extensively applied to geotechnical engineering,among which the two most prevalent strength criteria are widely employed, namely, the M-C strength criterion and the Hoek-Brown strength criterion.The former is usually applied to solve the problems in linear or nonlinear soil materials while the latter is typically utilized to address the problems in rock materials that show a palpable nonlinearity in strength.The superiority of the Hoek-Brown criterion in building a direct connection between material parameters and the practical state of rock masses has contributed to widespread acceptance by engineers compared with the other strength criteria.Below,a brief review of the Hoek-Brown strength criterion is given.

In 1980, the inceptive form of the Hoek-Brown strength criterion, as presented in Fig.1, was proposed, and after long-term application and development, a generalized form of the Hoek-Brown strength criterion was represented in a relationship between major stress and minor principal stress:

Fig.1.Schematic diagram of the Mohr circle and Hoek-Brown strength envelope with reference to Hoek [20].Note: σn is the the normal stress; σt the uniaxial tensile strength;T the stress vector on failure surface;δ the rupture surface of rock masses; v the velocity vector on failure surface;vn the normal velocity vector;vt the tangential velocity; and E* the isotropic tensile strength.

where σcis the uniaxial compressive strength of the intact rock;σ1the maximum principal stress; σ3the minimum principal stress;and mb, a, and s the material parameters depending on the rock state/quality.

where miis a rock type-dependent parameter; GSI the geological strength index, derived from a comprehensive evaluation of the geological environment,rock structural features and surface characteristics; and Dca coefficient that describes the disturbance situation of rock.Dc=0 indicates an undisturbed in situ state while Dc=1 indicates a fully disturbed state,and in this study,an undisturbed state is taken.For the situation where data are unmeasurable, Hoek [20] summarized a few typical values of miwith regard to different rock types for possible reference,and some more specific and interesting investigations on the parameter mican be found in recent works [35].

2.2.Applicability of limit analysis to rock mass

Practical experience indicates that rock and concrete exhibit brittleness under tensile stress and perform a certain ductility under compression.Relevant scholars have qualitatively proven that the plastic flow of rock and concrete would happen when enough confining pressure is provided, and clearly concluded that the transformation from brittleness to ductility is closely related to the confining pressure and temperature.Indeed, this ductility deformation of rock materials makes the application of the upper-bound limit analysis feasible.Chen [36] stated that when rock is assumed to be in zero tension strength or in limited compressive strength, it can be treated as a real material.The infinite plasticity assumption is precise and conservative in zero tension for the employment of the limit analysis theorem; however, this applicability is doubtful under compressive stress because the obvious brittleness under compression and the apparent dropping of stress–strain curve after ultimate strength observed in experiments do not conform with the assumption of the limit theorem.However,if the strain of rock and concrete in compression is small and merely appears, the ductility of rock and concrete shown in compression before the obvious dropping of stress may be enough to meet the requirements of limit analysis.

Prior to the application of kinematic analysis,a plastic shearing layer with zero thickness should be introduced,so that the velocity discontinuities become permissible in the limit state.The kinematics-based limit analysis theorem stipulates that in any kinematically admissible velocity field,the total work rate of external loads does not exceed the internal energy dissipation rate.

where the first term on the left is the internal energy dissipation rate occurring in the failure body; σijand ˙εijthe tensors of stress and strain rate,respectively;the second term on the left the internal energy dissipation rate occurring along the velocity discontinuity surface; [v]jthe velocity jump vector in the velocity discontinuity surface; nithe unit normal vector;vithe distributed loads; Xithe distributed loads; Tithe force acting on the failure surface; V the volume of failure block; L the velocity discontinuity surface; and S the boundary of failure block.To respect the normality rule of the plasticity theorem, the vector of the strain rate must be perpendicular to the limit state surface;in analogy,the velocity vector in the vn-vtplane should be perpendicular to the envelope, as depicted in Fig.1.The first term on the right is the external work rate occurs along the boundary surface S;and the second term in right the work rate produced by distributed loads Xisuch as gravity or inertial force.Predicated upon the relationship of Eq.(3), an upper-bound failure load that is no less than the true failure load is obtained.

2.3.Linearization of the Hoek-Brown criterion by the PLM

Kinematic analysis can be readily conducted in the framework of the plastic limit analysis theorem, however, the representation of the Hoek-Brown strength criterion in the form of principal stress is not convenient for further development.Thus, using the analytical procedure for the strength envelope of the generalized Hoek-Brown criterion proposed by Kumar [18], an alternative form of the formula in Eq.(1) is derived from stress space τ,σn() by introducing an angular parameter δ.The specific derivations can be found in reference [18]; here, only the eventual expressions of shear and normal and stresses are presented:

As indicated in Fig.2,it is also the angle δ that the external normal of the strength envelope and the horizontal line will be inclined at.Readily, the criterion in Eq.(1) on the principle stress plane is transformed onto the τ,σn() plane, which implies that for a random angle δ,stress components τ and σncan be calculated from Eqs.(4) and (5).

Fig.2.Piecewise linear method for approximation of Hoek-Brown strength envelope.Note: δi is the ith rupture angle; δm the mth rupture angle.

Based upon the transformation of the Hoek-Brown strength criterion, a piecewise linear method (PLM) is employed in this paper to substitute the strength envelope for several tangential lines with different rupture angles,as illustrated in Fig.2.This approximation makes the analysis more feasible and allows accounting for the nonlinearity of the strength envelope depending on the variable angle.The rupture angle gradually increases when the stress state on the rupture surface varies from high level to low level.For each tangential line,the corresponding rupture angle is considered constant, but the rupture angle between adjacent tangential lines is different.In addition, each tangential line corresponds to a unique point on the strength envelope,and the stress state determined by that point is considered identical.Although the upper-bound limit analysis is incapable of computing the real stress state of the rock mass, it is promising to provide a rough approximation of the stress distribution on the rupture surface.

3.Typical approximations on the Hoek-Brown criterion

3.1.The generalized tangential technique

The linear M-C strength criterion is able to deal with the most geotechnical engineering problems,however,the Hoek-Brown criterion due to its distinct nonlinearity fails to be directly introduced into the theoretical analysis.Therefore, a numerical procedure named GTT was first developed by Zhang and Chen [21] to investigate the problem of slope stability with nonlinear strength criteria.Thislinearizationmeasurefacilitatessubsequent developments of the upper-bound analysis smoothly and has been widely applied to many geotechnical structures.Using the same procedure employed by Zhang and Chen [21], the equivalent M-C parameters of the Hoek-Brown strength criterion are derived by Yang and Yin [22].

As presented in Fig.3,passing through a point P at the strength envelope, a unique tangential line is determined, and the corresponding tangential equation reads:

Fig.3.Linearization of the Hoek-Brown criterion by GTT with reference to Yang and Yin [22].

where ctand δtare the intercept and slope of the tangential line,respectively.Combining the parametric form of the Hoek-Brown criterion in stress space τ,σn() and the tangential equation in Eq.(6), the equivalent M-C strength parameters ctand δtyield the following relationship:

It has been proven that the limit failure loads obtained from the tangential technique is a rigorous upper bound of the actual limit loads because the material strength provided by the tangential line is always greater than or equal to that of the actual failure envelope.In fact, this single tangential line method is a specific case of the PLM elaborated in the former section, namely m=1.

3.2.Equivalent Mohr-Coulomb parameters method

Apart from the GTT, another typical method to address the Hoek-Brown strength criterion is the EMM.In the process of ceaselessly improving the Hoek-Brown strength criterion, Hoek [20]introduced the relationship among GSI,mb,s,and a,and considered the disturbance effects of stress release and blast on rock mass by using the disturbance coefficient Dc.Moreover,for the convenience of investigating the underground tunnel excavation and slope stability, the expressions of equivalent cohesion ctand internal friction angle δtin the M-C criterion are derived from the regression process that balances the areas above and below the equivalent M-C envelope.

As shown in Fig.4,the dashed line represents the equivalent MC envelope which is captured by curve fitting.The whole stress space is separated into three regions.Noted that the cohesion and shear stress tend to be overestimated and that the friction angle will be underestimated for the stress tensor in Region 1,while an underestimated cohesion and an overestimated internal friction angle and shear stress are obtained in Region 3.The equivalent shear strength parameters in the form of M-C take the following expressions:

Fig.4.Approximation of the Hoek-Brown strength criterion by the EMM with reference to Hoek [20].Note:φt is the equivalent internal friction angle.

where σ3n=σ3max/σc, and σ3maxis expressed as:

where γ is the unit weight of the rock mass; C the buried depth of the tunnel;and σcmthe compressive strength of the rock mass.σcmis normally calculated by Eq.(11):

The aforementioned procedure provides an approximated substitution for the Hoek-Brown strength criterion and the corresponding equivalent M-C parameters are also incorporated into the face stability analysis for comparison.

4.Kinematic analysis by piecewise linear method

4.1.Modified rotational mechanism of the tunnel face

In this section,a modified multi-horn rotational mechanism for a tunnel face excavated in rock masses is generated in detail.As presented in Fig.5,a rotational mechanism of a deeply buried tunnel of overburden C and diameter D is considered.This failure mechanism is composed of a completely collapsed block that revolves around the fixed axis passing through O performing a rotational motion.The upper boundary A1AiAmand lower boundary BBiBmrepresent rupture surfaces which are also explained as a narrow transition layer where velocity discontinuity is allowed.The boundaries of the failure mechanism consist of several log-spiral fragments where the stress reaches the yield level and the stress vector T presented in Fig.1 achieves the Hoek-Brown strength envelope.The boundary BBiBmis discretized into m segments while the number of segments for the boundary A1AiAmis m-1 resulting from the geometric feature of the failure mechanism.The generation of rupture surfaces BBiBmand A1AiAmstarts from the tunnel face.In terms of the first segment of the lower boundary, one can rotate the anticlockwise radial line OB of length r0in angle α1,determine the expression of arc BB1based upon the associative flow rule,and then generate the upper boundary from the crown of the tunnel face A1with the initial radial being.Adopting the same procedure, one can define the next segment until the radial line intersects at the same apex Am(Bm).Each segment is determined by two angular variables δiand αi, as shown in Fig.5.The rupture angle δiwithin segments between successive radial lines is constant.Notably, the log-spiral curves at points Aiand Biare not necessarily smooth but continuous.

Fig.5.Modified failure mechanism of a deep tunnel face with m segments.

In the ith sector, the upper and lower boundaries between radial lines OBi-1and OBiare defined by the following log-spiral expressions:

where θ is the polar angle of the rupture surface; θithe polar angle of Ai,θi=θi-1＋αi(i=1, 2, ∙∙∙, m);r′and r the lengths of the rotation center to the velocity discontinuity surface; rithe lower polar radius of the ith segment; and r′i the upper polar radius of the ith segment.(length of OAi) and ri(length of OBi) are calculated as:

where the subscript k is an index.To ensure that the failure mechanism is closed,the upper and lower boundaries should intersect at the same point Am(Bm); hence, setting rm=, Eq.(16) is derived.

In the triangle OA1B, the following relationship is readily derived based on the sine theorem.

where θ0is the polar angle of B.Combining Eqs.(16) and (17), the angle αmof the last segment yields

To facilitate the latter derivation, rsiis defined as the distance from the origin point O to A1B and the upper boundary A1AiAmand is analytically represented in the polar coordinate system as follows:

4.2.Establishment of the energy-work balance equation

Predicated on Eq.(3), the equilibrium of the internal work rate and external work rate is established.The internal energy is solely dissipated in the narrow transition layer while the external work rates are contributed by the gravity of the rock and uniformly dis-tributed support pressure.Then, an upper bound solution of the limit support pressure that is no less than the true support pressure is obtained by the following balance equation:

where Dsis the rate of energy dissipation;Wγthe work rate done by gravity; and Wσsthe work rate resulted from the required support pressure σs.The energy dissipated in the volume of rock is integrated by the following expression:

In practice, due to the rigid body assumption on the failure block, the plastic strain rate in the failure body is zero; thus, the energy dissipated in the volume of rock yields zero.Assuming that the thickness of the narrow layer is zero,the internal energy dissipation rate along the thin transition layer is calculated by the following integral formula:

Alternatively,the integrand function in Eq.(22)is rewritten as a product of stress and velocity on the rupture surface per unit area.

where v is equal to the length of velocity vector V; dithe energy dissipation rate per unit area;σnithe normal stress of the ith segment; and τithe shear stress of the ith segment.σniand τi, in the limit state surface, are two mutually perpendicular components of stress vector T.The summation of the internal energy dissipation rate occurring in all segments yields:

where ω is the angular velocity of the collapse block during the incipient of failure.

As illustrated in Fig.5, a micro element dA in the ith sector shown light gray is selected to perform a double integration for calculating the work rate produced by the weight of rock masses.

where ρ is equal to the length of micro element dA to the rotation center; and θi-1and θithe polar angles of successive boundaries in the ith sector with light gray.

The work rate produced by surface pressure uniformly applied to the tunnel face is presented as

where r0is the length of the first radial line.

Substituting the terms in Eq.(20)for work rates in Eqs.(24)–(26)and rearranging them,a dimensionless group σs/γD is derived:

Eq.(27) is taken as a multivariate nonlinear function σs/γD=f (αi,δi) where the number of dependent variables that need to be optimized is 2m.This function is difficult to solve by the conventional enumeration method; thus, an alternative algorithm named MPA is adopted as a solving tool.This algorithm is a new meta heuristic optimization algorithm proposed by Faramarzi et al.[34].Its inspiration stems from the survival of the sea fittest theory,that is,marine predators choose the best foraging strategy between the Lévy walk and the Brownian walk.The algorithm has the characteristics of a stronger optimization ability than the conventional genetic algorithm and sequential quadratic programming.The optimization should be conducted subjected to the following constraint conditions.

5.Kinematic analyses by classic approximations on the Hoek-Brown criterion

5.1.Classic failure mechanism of the tunnel face

Concerning the common measures for approximating the nonlinear strength criteria, the GTT proposed by Zhang and Chen[21] is widely used by scholars.In this section, the GTT is adopted to substitute the Hoek-Brown strength envelope for a tangential line.The corresponding results are calculated as a reference to the accuracy of the PLM.A similar rotational mechanism like Fig.5 is generated respecting the associative flow rule.Because the approximation of the Hoek-Brown strength criterion is a single tangential line, only one rupture angle, namely the equivalent internal friction angle, will be determined.Therefore, the velocity discontinuity of the mechanism will not be divided into several segments like that shown in Fig.5.A constant rupture angle exists for all rupture surfaces.

Fig.6 presents the active failure mode of a tunnel with a buried depth of C and a diameter of D.The horn-like mechanism ABE rotates with the axis passing through point O and is bounded by log-spirals AE and BE.The log-spirals start from the crown and convert of the tunnel face and intersect at apex E with an angle of 2δtdue to the normal flow rule.The expressions of the upper boundary and lower boundary of the horn-like mechanism are written as:

Fig.6.Failure mechanism for linear substitution of the Hoek-Brown envelope.

where rAis the polar radius of crown;rBthe polar radius of invert;θAthe polar angles of crown;θBthe polar angles of invert;and θEthe polar angle of apex E.rA(length of OA) and rB(length of OB) satisfy the following relationships:

To make the failure mechanism closed, the upper and lower boundaries should intersect at the same point E; hence, letting r =r′, θEreads:

For the convenience of developing the latter derivation,defining rsas the radial coordinate of the upper boundary EAB, thus, rsis analytically represented as follows in two distinct regions.

5.2.Establishment of the power balance equation

Analogous to the previous section, the rate of internal energy dissipation occurring in the velocity discontinuity takes the following form:

where ω is the angular velocity of the collapse block at the incipient failure state.

In the micro element dA,the integral expression of the external work rate induced by the weight of rock masses yields

The work rate produced by the support pressure uniformly acting on the tunnel face is presented as:

By substituting the terms in Eq.(20)for work rates in Eqs.(34)–(36) and rearranging them, the dimensionless group σs/γD is derived:

Eq.(37) is formulated as the objective function σs/γD=f (θA,θB,δt) which is a special case of Eq.(27), namely the number of segments is set to m=1.Additionally, the above objective functions are applicable to the EMM by substituting the corresponding shear strength parameters into Eqs.(27) and (37).Also, the optimization of Eq.(37) should be implemented under the following constraint conditions.

The specific optimization process of the required support pressure by the three methods is presented in Fig.7 where the objective function obtained by the PLM depending on 2m variables is optimized by MPA.The objective functions obtained by the EMM and GTT depending only on three variables are calculated by an enumeration algorithm.

Fig.7.Flowchart of the calculation process.

6.Results and discussions

6.1.Comparisons

As shown in Eq.(27), by providing a set of angle variables(αi,δi),a unique failure mechanism is correspondingly determined,based on which a unique support pressure can be calculated.By optimizing the σs/γD=f αi,δi(), the maximum support pressure is obtained.Existing works concerning the two-dimensional face stability of rock tunnels are insufficient, especially in theoretical research and numerical simulation.Thus, the correctness of the results in this work is validated from two aspects.On the one hand,based on the relationship between the GTT and the PLM, namely GTT is identical to the PLM when m=1.Therefore, the corresponding support pressure predicted by the GTT should be equivalent to that of the PLM with m=1.Such mutual verification between the GTT and the PLM has been carried out.On the other hand, to further improve the credibility of the GTT, results calculated by the GTT are compared with the solution by Zhang et al.[37], where a deep-buried tunnel in rock masses with D=10 m and C=20 m is considered, and the same GTT applied in this work is used to address the nonlinearity of the Hoek-Brown criterion.The comparative results in Table 1 show that the difference between the twoworks is small and results from the error of data extracted from Zhang et al.[37].In a word,the validity of the proposed procedure is ensured to some extent.

Table 1 Critical face pressures compared with Zhang et al.[37].

Prior to the stability analysis, the number of segments is of necessity to be defined because an increase in segments indeed gains a lot in accuracy but is also accompanied by a burden in computing time.Therefore, a compromise between accuracy and efficiency should be achieved.Based on the following parameters:D=10 m,Dc=0,mi=5,GSI=10,σc=1 MPa and γ=25 kN/m3,Table 2 shows the variation in the required support pressure with the number of segments.The required support pressure gradually improves with the increasing number of sectors, which is also accompanied by a linear increase in computing time.The optimization in one time ended with the conditions that the difference between the results of two consecutive calculations is smaller than 1×10-6with the increment in constraint variables being 0.001°.Considering insignificant improvement in the support pressure(less than 0.1%),the number of sectors is selected as 10 and applies to subsequent content.

Earlier investigations on face stability based on the Hoek-Brown strength criterion assume a linear substitution for the strength envelope [27,31].An outer tangential line of a certain point in the nonlinear envelope is adopted to approximate the entire nonlinear envelope.The corresponding equivalent shear strength parameters are derived from the relationship between the tangential line and envelope.The optimal tangential line is selected when the critical face pressure is the best.In this work, such an approximation indicates a constant rupture angle in the mechanism,which is a specific case of PLM, i.e., m=1.Before analysis, it is necessary to declare that the solutions calculated by the GTT and PLM are rigorous upper bounds while the results obtained by the EMM are probably not a rigorous upper bound since the linear envelope of the EMM is derived by the regression process rather than the outer tangential line by the GTT.

In terms of different diameters, Table 3 provides the required support pressure calculated from m=10 and m=1 for comparison.It is worthy of noting that the difference between the two methods reaches approximately 13%.This improvement is significant andimparts confidence in the method presented in this paper.Fig.8 compares the results calculated from the PLM, GTT and EMM.As expected,the support pressure given by the PLM is always greater than that of the GTT.The results by the EMM maintain consistency with the other two methods when miand GSI remain in a small range, whereas the difference becomes significant for the parameters making the support pressure less than 40 kPa.This difference will be discussed in the parametric study section.

6.2.Influences of rock strength properties

This subsection aims at studying the influences of rock strength properties, including the geological strength index GSI, uniaxial compressive strength σc/γD and strength parameter mi, on the required face pressure.The disturbance coefficient is taken as Dc=0, indicating undisturbed in situ rock masses.

The dimensionless support pressure σs/γD varying with σc/γD is presented in Fig.9a for a particular dataset, namely mi=7 and GSI=10.It is obviously noted from the results of Fig.9a that the face stability sharply decreased with a significant increase in the required support pressure for decreasing values of σc/γD (approximately σc/γD<5).The required support pressures predicted by the PLM are greater than those of the GTT and EMM,and the variation trend shows a well consistency with each other.Notably,the required support pressures calculated by the EMM merely show agreement with the PLM and GTT for σc/γD ≤10, approximately,but for σc/γD>10, the difference gradually becomes significant.The specific reason is explained in Fig.9b, where the equivalent M-C strength parameters approximated by the EMM and GTT are presented.The cohesion and friction angle estimated by the GTT shows opposite variation tendencies because the strength envelope appears to be an evident nonlinear feature in the low stress level.Compared with the GTT,the EMM tends to estimate a larger cohesion and a smaller internal friction angle.Although the differences in the equivalent M-C shear strength parameters between the two linearized methods are exacerbated,the disparity in cohesion estimated by the two methods is more significant because the variation trend of cohesion in the GTT is negatively correlated with the uniaxial compressive strength of rock,while in EMM,it is positively correlated with the uniaxial compressive strength of rock.According to the limit analysis theorem, an overestimated cohesion consequently leads to an overestimated internal energy dissipation rate, which produces a stable tunnel face and a smaller prediction of support pressure.

Fig.11.Impact of strength parameter mi on the required support pressure and equivalent M-C parameters with σc/γD=4 and GSI=10.

Figs.10 and 11 represent the variation in dimensionless support pressure with the GSI and the strength parameter mi,respectively.Similar to the dependence on the uniaxial compressive strength,σc/γD, the face stability decreases with the reduction in strength parameters,which is also expressed by an increase in the required support pressure.The kinematic method of limit analysis states that the solution derived in any kinematically admissible velocity field is an upper bound of the actual failure load.If a kinematically admissible velocity field derives the least upper-bound solution,we conclude that the corresponding mechanism is better.In this study, the support pressure produces negative work and is not a failure load.Therefore,a maximum upper-bound is deemed better.The PLM always predicts support pressure better than the other two methods, which demonstrates that the proposed multi-horn mechanism is better.Based upon a rough observation of Figs.10 and 11, the required support pressure predicted by the EMM becomes distinct from those of the PLM and GTT at GSI=20 and mi=15.Interestingly, these illustrations indicate a rather moderate dependence of the support pressure on the strength parameters GSI and micompared with the dependence on the parameter σc/γD.

The PLM substitutes the Hoek-Brown strength envelope for several segments.For a unique failure mechanism,one can derive a set of equivalent M-C parameters,namely,ciand δion the rupture surfaces.In each segment, the stress distribution along the rupture surface is constant.To better understand the approximation of the PLM on the nonlinearity of the strength envelope, two sets of mechanical parameters of rock masses are selected for illustration.As presented in Fig.12, the left two subgraphs are the strength envelope and the right two subgraphs are the corresponding M-C strength parameters.The blue results are calculated based on the stress state of every point in the Hoek-Brown envelope and are continuous.The pink results are estimated by the PLM and exhibit a stepped distribution near the equivalent M-C parameters.The solid line and dashed line represent the variation in the friction angle δtand cohesion ctwith normal stress, respectively.Through optimizing Eq.(27),δtis the value of the dependent variable δithat corresponds to the optimal support pressure.Also, using Eqs.(4)and (5), the stress states τiand σniused to calculate the cohesion by ci=τi-σnitanδiare obtained.In addition, it clearly shows that the nonlinearity of the rock strength results in a decrease in the equivalent friction angle and an increase in the equivalent cohesion when the normal stress increases.

Fig.12.Hoek-Brown strength envelope and corresponding equivalent M-C strength parameters for two heavily fractured rock masses.

As stated in Section 2,the Hoek-Brown strength criterion in Eq.(1) can be transformed onto the τ,σn() plane from the principal stress plane by introducing a rupture angle δ as a parameter.It is found that the stress components on the strength envelope are uniquely dependent on rupture angle δ, which makes it possible to roughly approximate the range of stresses acting on the rupture surface.Although the stress distribution so determined may not be the real stress because the stress fields associated with the kinematic field may not be necessarily in equilibrium,some meaningful reasons for the variation trend of solutions may be revealed.

Given the parameters D=10 m, mi=5, GSI=10, σc=1 MPa and γ=25 kN/m3, the stress components on the rupture surface are drawn on the normalized Hoek-Brown strength envelope, as shown in Fig.13.The bullets in blue, green and white correspond to the stress states calculated for tunnels with diameters of 6, 8 and 10 m, respectively.These results clearly show that the stress state on the rupture surface maintains a low compression level for tunnels with small diameters, whereas the stress distribution of the large-diameter tunnel remains at a relatively high level,and the stress range is more scattered.Interestingly, the tensile stress always exists at the apex of the failure mechanism regardless of the diameter of the tunnel.Additionally,the collapse region of the rock tunnel expands with the increasing diameter, and the stress distribution on the rupture surface involves a greater scale.The corresponding critical failure mechanisms are plotted in Fig.14 to provide an intuitive illustration of the stress distribution on the kinematic mechanism.The dashed pink lines and solid blue lines represent the rupture surfaces.The bullets are placed in the middle of each segment in the same colors shown in Fig.13.All the segments are numbered from 1 to 8, which is consistent with the numbering of the stress points in Fig.13.

Fig.13.Calculated stress states on the rupture surface mapped on Hoek-Brown strength envelope for tunnels with different diameters.

Fig.14.Calculated stress states on the rupture surface mapped on the corresponding failure mechanisms.

After parametric analysis,a set of design charts calculated based on the newly proposed failure mechanism is provided for a quick assessment of the required support pressure of a tunnel excavated in heavily fractured rock masses.Fig.15 shows the dimensionless support pressure of the tunnel face σs/γD as a function of the dimensionless uniaxial compressive strength σc/γD for different values of the rock strength parameters miand GSI.As previously mentioned, the face stability of a TBM-excavated tunnel in undisturbed surrounding rock masses (Dc=0) is considered.The ranges of parameters are selected to represent poor-quality rock masses,where the tunnel face is prone to instability in practice.Specifically,σc/γD varies from 1 to 100, mivaries from 5 to 20 in equal intervals, and GSI is set lower than 25, which represents very poor-quality rock masses.As expected, the decrease in support pressure when the rock strength properties increase is illustrated in Fig.15.In most cases, the values of the support pressure are low,and in some cases,the predicted support pressure is negative,which indicates that the tunnel face has the ability to self-stabilize without external pressure to support it.

Fig.15.Design charts of the required support pressure in Hoek-Brown rock masses.

6.3.Discussions

The purpose of this study aims to further improve the accuracy of previous linear substitution for the Hoek-Brown criterion and to analyze the difference in the different approximated methods,including the tangential technique, equivalent M-C parameters method and piecewise linear method.

The previous applications of the Hoek-Brown strength criterion mainly adopted a single straight line to substitute for the nonlinear envelope.Such a linearization method consists of the outer tangential technique and linear regression fitting.Therefore,the nonlinear rock materials are accordingly substituted with the linear materials characterized by the M-C criterion.The equivalent parameters of the M-C criterion are chosen to attain a theoretically best upper bound of the face pressure.The piecewise linear method adopted in this paper ensures that stability analysis can be conducted with consideration of the true nonlinear pressure dependency of the failure envelope.The previous works analyzed by the GTT,namely,a specific case of the PLM, are a compromise to the computation and complicated derivation to some extent, and loss of accuracy is therefore inevitable.The EMM is a simplified substitution for the nonlinear strength envelope in a distinct stress region.The results calculated from the EMM show agreement with other methods for low-quality rock materials while it shows poor applicability to rock masses with high quality.

The application of the piecewise linear method on the Hoek-Brown strength criterion requires a modification of the rotational mechanism of the tunnel face.Because in the conventional face stability analysis where the tangential technique is utilized to linearize the Hoek-Brown strength criterion, the boundaries of the rotational failure mechanism are continuous and smooth curves that are a function of a single internal friction angle, and such an equivalent friction angle is the slope of the tangential line of Hoek-Brown envelope.However, when the piecewise linear method is employed,the Hoek-Brown envelope is substituted with n tangential lines; simultaneously, n equivalent internal friction angles are obtained.Based on the associative flow rule,the boundaries of the failure mechanism should be divided into n segments,and each segment corresponds to a single equivalent internal friction angle.

In the derivation process, some basic assumptions and prerequisites are made as same as the conventional rotational mechanism: (1) rock masses are assumed homogeneous and isotropic,joints and fissures are uniformly distributed, and the strength of rock materials is governed by the Hoek-Brown criterion; (2) the failure block is considered a rigid body rotating around the polar origin O; and (3) the plastic flow in the limit failure state respects the associative flow rule.According to the kinematic theorem of limit analysis, the energy-work balance equation is established based on the principle of virtual work, the small deformation hypothesis and the ideal elastic–plastic hypothesis.The small deformation hypothesis is the basis of the principle of virtual work.The ideal elastic-plastic hypothesis is a prerequisite for the small deformation hypothesis.In practice, after the geomaterial reaches the yield stress level, it will enter the infinite plastic flow state,and the strain will accumulate with increasing time.The size of the failure area will change greatly after failure, which is not convenient for analysis and calculation.Therefore, the limit analysis method assumes that the research objects are in a critical state,and the stress–strain curve just enters the flow plastic stage from the elastic stage.At this time,the plastic flow just occurs,the strain is of elastic magnitude,and for the limit analysis,the elastic deformation is normally disregarded.Thus, the geometric relationship of geotechnical materials slightly changes after failure.

A necessary improvement of this work consists of extending the face stability problems from two-dimensional space to threedimensional space.Also,some complicated situations encountered during underground excavation, such as pore water pressure and dynamic loads deserve to be investigated.

7.Conclusions

The stability of TBM-excavated tunnel faces is investigated in the Hoek-Brown media.Three kinds of substitutions are introduced to account for the nonlinearity of the strength envelope by providing linear and piecewise linear approximations for the nonlinear strength envelope.In the case of linear substitution,the conventional curvilinear cone mechanism is adopted in the analysis.The failure mechanism corresponding to the PLM is modified with the multi-horn mechanism to ensure its compatibility with the piecewise linear envelope.The closed-form solution of the critical face pressure is obtained in the framework of the plastic limit theorem.The critical face pressure is formulated as a multivariable nonlinear optimization problem (2m dependent variables), which can be quickly searched by the MPA.The appropriate number of segments (namely, m=10) is suggested in terms of accuracy and efficiency.The influence of rock strength parameters on the critical face pressure is investigated in detail, and dimensionless charts of the support pressure are provided for possible guidance in practical design.

Based on the results of the parametric study, some meaningful conclusions are drawn:

The accuracy of the critical face pressure predicted by the PLM presents a significant improvement up to 13% compared with the conventional GTT, which imparts confidence to the PLM.

The results calculated by the PLM, EMM and GTT show consistency with each other when miand the geological strength index remain within a small range, whereas the difference becomes significant for the parameters making the support pressure less than 40 kPa.In other words,the solution obtained by the EMM is accurate to a limited extent.

The influence of rock strength parameters on the required support pressure shows obvious nonlinearity, and the face stability seems more sensitive to the uniaxial compressive strength than GSI and mi.

In this study,the multi-horn solution is interpreted as a solution for the true nonlinear strength envelope but with piecewise constant distributions of stress on the rupture surface.Although such an assumption does not ensure the stress equilibrium, the stress equilibrium may not be necessarily required in the kinematic method of limit analysis, thus it is still a rigorous upper bound to the support pressure.The presented PLM,in a word,is an extension of the GTT, in contrast, and the GTT is a special case of the PLM,namely,m=1.Furthermore,note that an isotropic rock mass is presumed to satisfy the Hoek-Brown strength criterion, which is an essential prerequisite that makes it applicable to heavily fractured rock masses,while for the problems associated with high deformations of material subjected to high stresses,the solutions derived in the present paper are not applicable.

Acknowledgements

This work was supported by Fundamental Research Funds for the central universities of Central South University (No.2022ZZTS0153).The financial support is greatly appreciated.