Coding the “New Implicit Method” in MatLab for preliminary tunnel design

Implementação no MatLab do “Novo Método Implícito” para pré-dimensionamento de túneis

Abstract

The tunnel design is a subject that demands a three-dimensional analysis dealing with two different players: lining and rockmass. The traditional methods of getting the necessary parameters for design do not consider a correct interaction between the two players, but the New Implicit Method (NIM) takes this characteristic in its core and develops some formulations for elastic, plastic and viscous rockmasses. Understanding that the stiffness of the lining and the distance between the lining and the tunnel face change the convergence of the tunnel, a code in MatLab for NIM is validated through Finite Element (FEM) with its results being presented. The validation of this method was compared with FEM analysis and the results obtained an average accuracy of 12% what represents a good approximation regarding geotechnical issues.

Keywords:  tunnel design; new implicit method; elastic; plastic and viscoplastic rockmass; convergence

Resumo

O dimensionamento de túneis requer uma análise tridimensional lidando com duas estruturas muito diferentes: o maciço e o revestimento. Os métodos simplificados tradicionais de dimensionamento de túneis não consideram a correta interação entre essas duas estruturas. O “Novo Método

Implícito” leva em consideração a correta interação e é desenvolvido para maciços elásticos, plásticos e viscoplásticos. O parâmetro chave desta interação é a convergência no instante da colocação do revestimento que é afetada, entre outros, pela rigidez do revestimento. Um código em MatLab do NIM foi desenvolvido e validado com resultados de cálculos numéricos em elementos finitos. Os resultados desta comparação fornecem uma precisão de 12%, demonstrando um ótimo resultado para obras geotécnicas.

Palavras-chave:  dimensionamento de túneis; Novo Método Implícito; maciços elásticos; plásticos e viscoplásticos; convergência

1. Introduction

The design of tunnel support is one of the most complex tasks in the tunnel design, and most of the projects remain empirically analyzed and designed. However, it is getting unsuitable, since it is being proved that the interaction of the soil and lining matters in a way to change the stresses taken by the concrete as well as the deformation caused by the rockmass. This deformation of the rockmass over the lining is the main parameter and it is going to be called convergence.

The Association Française des Tunnels et de L’Espace Souterrain (AFTES) [¹] classifies the tunnel support design methods into four types: (i) purely empirical methods indicating the most appropriate type of support for a situation defined from various geotechnical classification systems; (ii) methods for determining the loads acting on the support, regardless of support type and deformation; (iii) support design methods which consider loads exerted by the ground as input data with some integration of support stiffness, deformation and the reactions of the surrounding ground; (iv) and more recent methods taking full account of the ground/support interaction.

The first method consists on the reproduction of succeeded tunnel supports and does not provide an assessment of the safety factor. The second and third method derives from the strength of the materials and consider the tunnel a semi-circular beam reacting to external loads. The support reactions are given as Winkler springs to reproduce the strength of the soil.

The last classification is derived from the continuous media theory and this is the reason why it was chosen to be applied in this work. This classification integrates the traditional Convergence-Confinement Method (CV-CF) for analytical solution and some numerical solutions with Finite Element which can describe better the 3D effect around the tunnel face.

The New Implicit Method (NIM) was developed by BERNAUD and ROUSSET [³] as evolution of the Convergence-Confinement Method. The NIM fully integrates the influence of the lining on the final convergence at equilibrium, since it shows that stronger the lining is, the smaller is the convergence (dimensionless displacement). Regarding that, the paper’s main objectives are: (1) to present the New Implicit Method and (2) to introduce the NIM in MatLab for preliminary tunnel design.

Every relevant and recurrent parameter shown in the text are summarized on Appendix I at the end.

2. Interaction problem

The determination of the load transferred to the support requires an analysis of the interaction of the load-deformation characteristics of the elements comprising the system ground/support. The preliminary design of tunnels foresees the necessary strength of the lining since it is able of presenting the convergences of the rockmass.

Let us consider a circular deep tunnel (radius R_i), driven in a homogeneous and isotropic rockmass [Figure 1], initially subjected to a geostatic stress field (σ₀):

${\sigma }_0=-P_{\infty }\underset{\_}{\underset{\_}{1}}with{P}_{\infty }=\rho gz(z=meandepth)$ (1)

Figure 1
Simplification of the behavior of the problem. (BERNAUD and ROUSSET, 1994, [³])

The Equation (1) gives the initial stresses (before the excaThe Equation (1) gives the initial stresses (before the excavation). It is dependent on the P_∞ which is the pressure exerted by the soil mass above the tunnel. According to PARK [⁹] since it is a circular deep tunnel, it is accepted to consider the same pattern for vertical and horizontal axes. Parameters ρ and g are the unit mass and gravity acceleration respectively.

With the assumption that the lining is a ring of constant thickness e, made of a homogeneous and isotropic material; and, moreover, that the lining is set at a distance d₀ from the plane and vertical tunnel face, the problem becomes an axisymmetric problem dependent only on the radial distance r to the tunnel axis and the time t.

In the other hand, the deep section of the tunnel can be treated as a Plane Strain problem since the lining speed is constant (time dependent), so the parameters will depend only on the radial distance r as well.

Firstly, the most used method to predict the stresses upon the lining is the Convergence-Confinement method [¹] which is a procedure that allows the load imposed on a support installed behind the face of a tunnel to be estimated. Per CARRANZA and FAIRHURST [⁵], when a section of support is installed in the immediate vicinity of the tunnel face, it does not carry the full load to which it will be subjected eventually; a part of the load is redistributed around the excavation face itself. However, as the tunnel and face advance, this ‘face effect’ decreases and the support must carry a greater proportion of the load that the face had carried earlier.

The method explained in the previous paragraph is defined by two curves [Figure 2]: the convergence curve is the dashed line, also known as the ground reaction curve that gives the internal pressure exerted at the wall versus the convergence of the tunnel’s wall and it depends on the characteristics of the rockmass and the adopted plastic criteria; the second curve is the confinement one (or support characteristic curve), which illustrates the relationship between the pressure P_i versus the closure of the lining U_i ^s defined as the straight line in [Figure 2].

Figure 2
Tunnel equilibrium in a diagram Pi x Ui (a) and equilibrium of a tunnel (b) (BERNAUD and ROUSSET, 1996, [³])

The Confinement curve starts at the convergence of the installation of the lining U_o and its slope is defined as the lining strength K_s.

The parameters of Figure 2 A are:

U_o and P₀: the convergence (radial deformation) and the pressure at the distance d₀ from the tunnel face;

U_eq and P_eq: the convergence and the pressure at a section far from the tunnel face at equilibrium, where both pressure and convergence are stabilized;

R_p: the plastic radius is the limit between the elastic and plastic zone around the tunnel.

Figure 2 illustrates the relationship between the pressure P_i versus the closure of the lining U_i ^s (convergence). Since it is shown in Figure 2 (A) that the rockmass presents an elastoplastic behavior, there is a plastic zone around the tunnel (B). The convergences U_o and U_eq are shown in (B) and the pattern of initial increase and final stabilization of the rockmass displacement is shown as well.

Eventually, the Convergence-Confinement curve is drawn and the safety factor can be obtained by comparing the concrete strength with the stress at a distance far from the tunnel face that is imposed by the rockmass upon the lining structure.

3. The lining stiffness

PANET [⁸] brings equations for the lining stiffness (K_s). There are two main equations which describe the lining stiffness. They are derived from a thin shell and a thick pipeline. This work uses the normal solution only, considering that in deep tunnel, there is a uniform convergence only. For other cases, Panet [⁸] shows the lining stiffness for bending resistance. Considering a thin shell and R_i /_e>10, the lining stiffness is:

$K_s=\frac{E_L}{1-{\nu }_L{}^2}\frac{e}{R_i}$ (2)

In the other hand, when a thick pipeline is used (R_i /_e<10), the lining stiffness is:

$K_s=\frac{E_L\left[R_i^2-\left(R_i-e\right)²\right]}{\left(1-{\nu }_L\right)\left[\left(1-2{\nu }_L\right)R_i^2+\left(R_i-e\right)²\right]}$ (3)

Where:

R_i: radius of the excavated tunnel;

e: thickness of the lining;

E_L: Modulus of Young of the lining;

ν_L: Poisson coefficient of the lining.

4. The New Implicit Method (NIM)

The Convergence at equilibrium is the final displacement imposed on the lining ring and it has been calculated by the NIM, the Pressure at equilibrium is defined by the NIM too. The Method derived from the CV-CF method contains four main steps and for every rockmass behavior (elastic, elastoplastic and viscoplastic) the steps follow the same sequence (see Figure 3).

Figure 3
Sequence of solution of a tunnel convergence U_eq and Pressure P_eq

The New Implicit Method (NIM) was developed by BERNAUD and ROUSSET [³] and consists on deducing the curve of convergence with a geometrical transformation from unsupported to supported tunnel. The transformation is a function a⁰ (x) of the unsupported tunnel and it was approximated by:

$a^0\left(x\right)=\frac{{\text{U}}_i^f\left(x\right)-U_i^f\left(0\right)}{U_i^f\left(\infty \right)-U_i^f\left(0\right)}=1-{\left(\frac{0.84R_p}{x+0.84R_p}\right)}^2$ (4)

Where U_f = U_i ^f (0) is the convergence at the tunnel face of the unsupported tunnel and U_i ^f (∞)=U_∞ is the convergence of the unsupported tunnel far from the tunnel face. The NIM showed that a^s (x) depends on the K_s (lining stiffness) that had not been predicted by the CV-CF method; so, a function for the supported tunnel is proposed:

$a^s\left(x\right)=a^0\left(\alpha x\right),with\alpha =\alpha \left(K_s^{\text{'}}\right)$ (5)

Then, for a given constitutive law of rockmass, the application of the NIM can be summarized by the following steps:

• a) Choose the reduced lining function α^* (K_s' ).Then, for the supported tunnel, Equation (5) becomes:

$x^{\prime }=\frac{x}{R_i};{\alpha }^{\text{*}}=\alpha \frac{R_i}{R_p}$ (6)

$a^s\left(x\text{'}\right)=1-{\left(\frac{0.84}{{\alpha }^{\text{*}}{x}^{\prime }+0.84}\right)}^2$ (7)

• b) Solve the convergence of the tunnel face U_f according to the rockmass behavior.

• c) The solution of the problem is then obtained by solving the following equations:

$\begin{array}{l}{P}_{eq}=C_F\left(U_0,U_{eq}\right)\\ P_{eq}=C_V\left(U_{eq}\right)\end{array}$ (8)

• d) In the System (8), the convergence U₀ is a function of U_eq (that explains the name of ‘implicit method’). By (7) and (8) we can write:

$U_0=a^s\left(d_0\right)\left(U_{eq}-U_f\right)+U_f$ (9)

Figure 4 sums up how the NIM method works integrating the two charts which correlate the convergence, the pressure and the location of the studied section by setting the lining supported function and getting the value for U₀ and its respective pressure.

Figure 4
Calculation of a tunnel with the NIM

Another use of the simplified method is to determine the ratio λ₀ (or the convergence U₀) at the moment of lining installation and then use this value as a good approach to calculate the tunnel in a 2D plane strain condition with a more complex geometry and loading (anisotropic loading for instance).

Figure 5 shows how the pressure P_i is obtained by using the ratio λ₀.

Figure 5
Application of the ratio λ₀ in a plane strain model

4.1 The Numerical Code GEOMEC91

The tunnel face advance process is studied here by the ‘deactivation/reactivation’ method developed in the numerical code ‘GEOMEC91’ [⁴] It is a 2D axisymmetric analysis of tunneling that considers the sequences of excavation and lining placement as it is illustrated in Figure 6. The excavation process is simulated by the ‘deactivation’ of the stiffness of the elements representing the ground removed during the excavation step (face advance). This is accomplished by a significance reduction of their Young’s Modulus and Poisson’s ratio.

Figure 6
Typical mesh of GEOMEC91 (BERNAUD and ROUSSET, 1996 [⁴])

Tunnel support is placed using the reactivation process at a specified distance, d₀, from the advancing face. This is calculated by substituting the lining characteristics in the corresponding lining elements. At this moment, the lining elements are stress and strain free.

BERNAUD [⁴] validated this numerical code with the experimental values of the Boom clay in Belgium, and for a 230-m deep tunnel the in-situ convergence was about 2.20% and the numerical solution returned a convergence of 2.04% presenting an average difference of 0.16% (between GEOMEC91 and experimental measures) showing the accuracy of the CEOMEC91 used here as the validation tool of the New Implicit Method (see Figure 6).

QUEVEDO [¹⁰] and FIORE et al [⁶] compared the solution of the GEOMEC91 with the numerical solution given by software ANSYS and showed an accuracy of GEOMEC91 used here. The model used by Fiore presented an accuracy of 5% when compared to the GEOMEC91 solution. The values for the convergence and stress at equilibrium are two of the driving characteristics to build a reliable and safe tunnel, and as the text follows, the solution of the incognitos cited before are going to be given for different behavior of the rockmasses depending upon the geology, type and creep characteristic of the rock around the lining.

4.2 Elastic rockmass

Some dimensionless parameters are going to be introduced to reduce the number of parameters shown in text. For an elastic solution, the dimensionless parameters are:

$P_{\infty }^{\text{'}}=\frac{P_{\infty }}{E};K_s^{\text{'}}=\frac{K_s}{E};d_0^{\text{'}}=\frac{d_0}{R_i}$ (10)

Where:

P_∞: Initial geostatic stress as defined in Equation (1); E: Modulus of Young of the rockmass;

K_s: stiffness of the lining.

The steps to follow by the NIM to find the convergence and pressure at equilibrium are as follow:

• a) Get the entry parameters: P^∞, E, K_s, d₀, ν (Poisson coefficient of the rockmass) and R_i;

• b) For an unsupported tunnel, we have:

$U_{\infty }=\frac{1+\nu }{E}{P}_{\infty };U_f=0.27U_{\infty }$ (11)

• c) For elastic, and plastic solutions, the function a^s (x) is given from Equation (7) in function of a dimensionless lining function (α*), and it is illustrated in Figure 7 proving that this function is dependent on the lining stiffness:

${\alpha }^{\text{*}}=\alpha \frac{R_i}{R_p}with{R}_p=R_{i}forelasticrockmass$ (12)

Figure 7
(a) Lining function α* for N_s and (b) dimensionless average lining function α* for N_s (BERNAUD and ROUSSET, 1996 [⁴])

The approximate method by BERNAUD and ROUSSET[³] proposed an average function for the reduced lining function, so α* = α*_average. See Figure 7:

${\alpha }^{\text{*}}={\alpha }^{\text{*}}{}_{average}=1.82\sqrt{K_s^{\text{'}}}$ (13)

• d) So, the function a^s (x') from Equation (7) is going to be changed by replacing the average function of α* from Equation (13):

$a^s\left(x\text{'}\right)=1-{\left(\frac{0.84}{1.82\sqrt{K_s^{\text{'}}}{x}^{\prime }+0.84}\right)}^2$ (14)

$a^s\left(d_0^{\text{'}}\right)=1-{\left(\frac{0.84}{1.82\sqrt{K_s^{\text{'}}}{d}_0^{\text{'}}+0.84}\right)}^2$ (15)

• e) Calculate the convergence and stress at equilibrium:

$U_{eq}=\frac{P_{\infty }+K_s{U}_f\left(1-a^s\left(d_0^{\text{'}}\right)\right)}{K_s\left(1-a^s\left(d_0^{\text{'}}\right)\right)+\frac{E}{1+\nu }}$ (16)

$P_{eq}=P_{\infty }-\frac{E}{1+\nu }{U}_{eq}$ (17)

• f) The convergence at d₀ is given by Equation (9).

4.3 Plastic rockmass

The rockmasses subjected to plastic strain presenting cohesion and/or friction angle will be treated here with Tresca (ϕ = 0) and Coulomb (ϕ > 0) criteria. However, some other criteria like Hoek-Brown is valid and presented in the original text of NIM [³].

$P_{\infty }^{\text{'}}=\frac{P_{\infty }}{E};K_s^{\text{'}}=\frac{K_s}{E};d_0^{\text{'}}=\frac{d_0}{R_i};C^{\prime }=\frac{C}{E}$ (18)

4.3.1 Tresca criterion

This criterion is valid for rockmasses whose main geotechnical characteristic is given by the cohesion. The steps to follow by the NIM to find the convergence and pressure at equilibrium are as follow:

• a) Get the entry parameters: P_∞, E, K_s, d₀, C, ν = 0.5 and R_i.

• b) For Tresca, the Friction Angle (ϕ) is zero and a new parameter is introduced:

$N_s=\frac{P_{\infty }}{\text{C}}$ (19)

• c) The reduced dimensionless lining function is given in Equation (13).

• d) Calculate the convergence U_∞ and at the tunnel face U_f, as follows:

$U_{\infty }=\frac{1.5\text{C}}{E}{e}^{\left(N_s-1\right)}forincompressiblematerial\left(v=0.5\right)$ (20)

$\frac{U_f}{U_{\infty }}=0.413-0.0627N_{s}for1<N_s\le 5$ (21)

• e) The a^s (d₀') is the same as the shown one in Equation (14);

The solution of this criterion valid for P_eq < P_∞ - C; if P_eq > P_∞ - C the rockmass is elastic. So, solving the System (8), the convergence and the stress at equilibrium are as follow:

$a\ln U_{eq}+b{U}_{eq}+c=0$ (22)

$\{\begin{array}{c}a=-\frac{C}{K_s}\\ b=a^s\left(d_0^{\text{'}}\right)-1\\ c=\left(1-a^s\left(d_0^{\text{'}}\right)\right)U_f+\frac{P_{\infty }}{K_s}+\frac{C}{K_s}\left(\ln \left(\frac{1.5C}{E}\right)-1\right)\end{array}$ (23)

$P_{eq}=\left(U_{eq}-U_0\right)K_s$ (24)

• g) The plastic radius is:

$R_p=R_i{e}^{-\left(\frac{P_{eq}-P_{\infty }}{2C}+1\right)}$ (25)

For Mises, the solution is given by the System (8) with C replaced by $\raisebox{1ex}{$2C$}\!\left/ \!\raisebox{-1ex}{$\sqrt{3}$}\right.$.

4.3.2 Coulomb criterion

This criterion is valid for rockmasses whose main geotechnical characteristic is given by cohesion and friction angle when $P_{eq}<\frac{2}{K_p+1}\left\{P_{\infty }+\frac{H}{2}\left(1-K_p\right)\right\}$ . The steps to follow by the NIM to find the convergence and pressure at equilibrium are as follow:

• a) Get the entry parameters: P_∞, E, K_s, d₀, C, ϕ, ν = 0.5 and R_i;

• b) For Coulomb, two new parameters depending on the friction angle (ϕ) are introduced:

$R_c=\frac{2Ccos\varphi }{1-sin\varphi }$ (26)

$N_s=\frac{2P_{\infty }}{R_c}$ (27)

When the Friction Angle is zero, the solution is for Tresca criterion shown in Equation (22).

• c) For Coulomb, the reduced lining function (13) is changed by the insertion of the influence of the friction angle on it:

${\alpha }^{\text{*}}={\alpha }^{\text{*}}{}_{average}+0.035\varphi =1.82\sqrt{K_s^{\text{'}}}+0.035\varphi$ (28)

• d) The Equation (28) changes the Equation (7) for a^s (d₀') as follows:

$a^s\left(d_0^{\text{'}}\right)=1-{\left(\frac{0.84}{\left(1.82\sqrt{K_s^{\text{'}}}+0.035\varphi \right)d_0^{\text{'}}+0.84}\right)}^2$ (29)

• e) Calculate the parameters K_p and H for Coulomb:

$K_p=\frac{1+\sin \varphi }{1-\sin \varphi }$ (30)

$H=C\frac{\cos \varphi }{\sin \varphi }$ (31)

• f) Calculate the convergence U_∞ and at the tunnel face U_f from Equation (21), as follows:

$U_{\infty }=\frac{1.5}{E}\frac{K_p-1}{K_p+1}\left(P_{\infty }+H\right){\left\{\frac{2}{K_p+1}\frac{P_{\infty }+H}{H}\right\}}^{2/\left(K_p-1\right)}$ (32)

• g) The solution of the System (8) gives the convergence and stress at equilibrium:

$a{U}_{eq}^{n_1}+b{U}_{eq}^{n_2}+c=0$ (33)

$\left\{\begin{array}{c}a={\left\{\frac{\left(K_p+1\right)E}{1.5\left(K_p-1\right)\left(P_{\mathrm{\infty }}+H\right)}\right\}}^{\frac{K_p-1}{2}}\left(\frac{K_p+1}{2}\right)\frac{K_s}{E}\left(1-a^s\left(d_0\right)\right)\\ b={\left\{\frac{\left(K_p+1\right)E}{1.5\left(K_p-1\right)\left(P_{\mathrm{\infty }}+H\right)}\right\}}^{\frac{K_p-1}{2}}\left(\frac{K_p+1}{2}\right)\left\{U_f\frac{K_s}{E}\left(a^s\left(d_0\right)-1\right)+\frac{H}{E}\right\}\\ c=-\frac{\left(P_{\mathrm{\infty }}+H\right)}{E}\\ n_1=\frac{K_p+1}{2}{;n}_2=\frac{K_p-1}{2}\end{array}\right.$ (34)

$P_{eq}=\left(U_{eq}-U_0\right)K_s$ (35)

4.4 Viscoplastic rockmass

The NIM for viscoplastic rockmasses was developed as the same way as the elastic and plastic rockmasses, so the solution is given by getting a transformation equation from unsupported to supported tunnel. The reduced velocity of excavation V* is the new parameter for viscoplasticity and it is written as:

$V^{\text{*}}=\frac{\eta V}{R_{i}C}$ (36)

Where:

η: viscosity constant.

The transformation function obtained from the unsupported tunnel is changed and it is given as follows:

$a_{vp}^0\left(x\right)=\frac{U_i^f\left(x\right)-U_i^f\left(0\right)}{U_i^f\left(\infty \right)-U_i^f\left(0\right)}=1-{\left(\frac{0.84}{\frac{x}{A{R}_p}+0.84}\right)}^{2B}$ (37)

The parameters A and B in Equation (37) are:

$A\left(V^{\text{*}}\right)=0.018+0.982e^{\left(-{\left(\frac{\log \left(1+V^{\text{*}}\right)}{2.8}\right)}^{6.4}\right)}$ (38)

$2B\left(N_s,V^{\text{*}}\right)=b_{02}(p_1+\left(1-p_1\right)e^{\left(-{\left(\frac{\log \left(1+V^{\text{*}}\right)}{2.3}\right)}^{3.3}\right)}$ (39)

$b_{02}\left(V^{\text{*}}\right)=0.195+1.805e^{\left(-{\left(\frac{\log \left(1+V^{\text{*}}\right)}{2.71}\right)}^{4.84}\right)}{p}_1\left(N_s\right)=e^{\left(1.756-0.878N_s\right)}$ (40)

However, for the supported tunnel, as it was developed for the elastoplastic rockmasses, the relationship between the unsupported and supported equation is given as:

$a_{vp}^s\left(x\right)=a_{vp}^0\left(\alpha x\right)$ (41)

$a_{vp}^s\left(x\right)=1-{\left(\frac{0.84}{\frac{\alpha x}{A{R}_p}+0.84}\right)}^{2B}$ (42)

• a) The calculation of the lining function α is different from the ones presented before for elastoplastic solution and it is as follows:

${\alpha }^{\text{*}}=\alpha \frac{R_i}{R_p}=K_s{{}^{\text{'}}}^{C_1}{e}^{C_2}$ (43)

$\begin{array}{l}{C}_1\left(N_s,V^{\text{*}}\right)=p_{c1}-\left(p_{c1}-0.5\right)e^{\left(-{\left(\frac{\log \left(1+V^{\text{*}}\right)}{3.54}\right)}^{m_{c1}}\right)}\\ p_{c1}\left(N_s\right)=e^{\left(-0.118+0.817N_s\right)}\\ m_{c1}\left(N_s\right)=e^{\left(2.094+0.0467N_s\right)}\end{array}$ (44)

$\begin{array}{l}{C}_2\left(N_s,V^{\text{*}}\right)=p_{c2}-\left(p_{c2}-0.6\right)e^{\left(-{\left(\frac{\log \left(1+V^{\text{*}}\right)}{3.70}\right)}^{m_{c2}}\right)}\\ p_{c2}\left(N_s\right)=e^{\left(-0.713+1.172N_s\right)}\\ m_{c2}\left(N_s\right)=-0.564+2.51N_s\end{array}$ (45)

• b) Calculate the convergence at the supported tunnel face (U_f), as follow:

$\frac{U_f}{U_{\infty }}=p_2+\left(0.413-0.0627N_s-p_2\right)e^{\left(-{\left(\frac{\log \left(1+V^{\text{*}}\right)}{2.70}\right)}^{2.9}\right)}$ (46)

$p_2\left(N_s\right)=\frac{0.27N_s}{e^{N_s-1}};1<N_s\le 5$ (47)

• c) The $a_{vp}^s\left(d_0\text{'}\right)$ gets:

$a_{vp}^s\left(d_0\text{'}\right)=1-{\left(\frac{0.84}{\frac{d_0\text{'}}{A}{\alpha }^{\text{*}}+0.84}\right)}^{2B}=1-{\left(\frac{0.84}{\frac{d_0\text{'}}{A}{K}_s{{}^{\text{'}}}^{C_1}{e}^{C_2}+0.84}\right)}^{2B}$ (48)

Since it is a viscoplastic solution, a criterion for the plasticity surface is required, the solution shown in the next topics is given for Tresca criterion presented before.

4.4.1 Tresca solution

• a) The convergence U^∞ is:

$U_{\infty }=\frac{1.5C}{E}{e}^{\left(N_s-1\right)};\frac{U_f}{U_{\infty }}from\left(45\right)$ (49)

• b) The solution of the System (8) is:

$a\ln U_{eq}+b{U}_{eq}+c=0$ (50)

$\{\begin{array}{c}a=-\frac{C}{K_s}\\ b=a_{vp}^s\left(d_0\right)-1\\ c=\left(1-a_{vp}^s\left(d_0\right)\right)U_f+\frac{P_{\infty }}{K_s}+\frac{C}{K_s}\left(\ln \left(\frac{1.5C}{E}\right)-1\right)\end{array}$ (51)

$P_{eq}=\left(U_{eq}-U_0\right)K_s$ (52)

• c) The plastic radius is given by Equation (25).

Figure 8 illustrates the influence of the reduced velocity on the convergence profile of the tunnel, showing that the faster it gets, the smaller is the convergence. According to BERNAUD and ROUSSET [³], when the velocity trends to infinity, the solution trends to the elastic convergence, when there is no enough time to generate residual strain.

Figure 8
The influence of the reduced velocity on the convergence. (BERNAUD and ROUSSET, 1994 [³])

5. Results and validation

The code is developed in MatLab and gives the .gui executable. The code has been developed with the equations of the New Implicit Method shown before, and afterwards some response charts and solutions are given to compare the solution of NIM with GEOMEC91, the finite element tool described in Item 4.1.

5.1 Elastic solution

The first results are shown for elasticity used for instantaneous strain pattern, we adopted:

$E=500MPaand\nu =0.498$ (53)

The solution obtained for the convergence through the Finite Element Analysis (U_eq ^FE (%)) was taken from the GEOMEC91 solution introduced by [⁴]. The errors between the FEM and NIM are shown and presented good approximation as shown in Table 1.

Table 1
Validity of the method for elastic rockmass

Constitutive law	P∞'	d0'	Ks'	Ueq FE (%)	Ueq NIM (%)	$\frac{U_{eq}{}^{NIM}-U_{eq}{}^{FE}}{U_{eq}{}^{FE}}$
Elasticity	0.008	0	0.72	0.87	0.76	12.64%
	0.008	0.67	7.2	0.90	0.93	-3.33%
	0.008	1.67	72	1.02	1.11	-8.82%

U_eq ^FE,NIM: solution for the convergence at a far distance from d₀ through Finite Element Method and NIM method analysis respectively.

Figure 9 shows the results obtained through the .gui executable for the second d₀' shown in the Table above, showing the charts developed automatically by the executable for stiffness varying from 0 to 20 000 MPa.

Figure 9
Solution for an elastic rockmass using .gui MatLab

The executable generates three charts which are the solutions for many supporting stiffness, it also generates a .xlsx file with the numerical solution for each supporting stiffness, the dashed line represents the solution for P_∞' = 0.008, d₀' = 0.67 and K_s' = 72.

5.2 Plastic solution

The results shown in Table 2 concerns the plasticity with a Tresca and a Coulomb criterion. The plastic solution is used in over consolidated soils which present residual strain. For this analysis, we adopted:

$E=500MPaand\nu =0.5$ (54)

P_eq ^FE,NIM: solution for the pressure at a far section from d₀^ through Finite

Table 2
Validity of the method for plastic rockmass

Constitutive law	P∞'	d0'	Ks'	Peq FE	Peq NIM	$\frac{P_{eq}{}^{NIM}-P_{eq}{}^{FE}}{P_{eq}{}^{FE}}$
Tresca
C = 1 MPa C = 1.33 MPa	0.008 0.008	0.33 1.00	7.2 72	1.48 0.90	1.50 0.85	-1.35% 5.56%
Coulomb
C = 1 MPa C = 1 MPa C = 1.33 MPa	0.008 0.008 0.01	0.67 1.67 1.67	0.24 2.4 24	0.76 0.35 0.309	0.75 0.31 0.31	1.32% 11.43% -0.32%

Element and NIM method analysis respectively.

The results for Coulomb and Tresca criteria present some errors smaller than 12% which is acceptable if compared to the results of the Finite Element analysis, indicating that the NIM is an adequate method for predicting the convergence and stress at equilibrium.

5.3 Viscoplastic solution

The results for viscoplastic cases are using the Tresca criterion as the plasticity surface. Some soils, mainly clays, present a creep behavior what increases the strain after some time and this behavior is important to assess. The results shown here have been compared to the ones given by GEOMEC91 as well. The results for many viscoplastic examples are given in Table 3.

Table 3
Validity of the method for viscoplastic rockmass

Ns	d0'	V*	Ks'	Peq FE	Ueq FE (%)	Peq NIM	Ueq NIM (%)	$\frac{P_{eq}{}^{NIM}-P_{eq}{}^{FE}}{P_{eq}{}^{FE}}$	$\frac{U_{eq}{}^{NIM}-U_{eq}{}^{FE}}{U_{eq}{}^{FE}}$
2	0.33	1000	0.72	0.86	1.06	0.77	1.11	-10.43%	4.72%
	1		36	0.72	0.38	0.79	0.37	10.33%	-2.63%
3	0.33	2000	18	1.78	0.39	1.77	0.39	-0.68%	0.00%
	1	150	0.24	0.92	2.8	0.92	2.79	0.51%	-0.36%
			2.4	1.59	2.0	1.74	1.86	9.54%	-7.00%
		500	0.72	1.49	2.1	1.44	2.15	-3.09%	2.38%
			7.2	1.97	1.66	2.19	1.49	11.02%	-10.24%
5	1	500	0.72	1.53	1.93	1.60	1.77	4.67%	-8.29%

As it is shown above, the method works well for every type of condition.

Figure 10 shows the results of the executable for a viscoplatic rockmass, the three charts given by the .gui are in function of the reduced velocity instead of the lining stiffness as it was for the elastic and plastic solutions shown before, the dashed line shows the solution for N_s = 5, d₀' = 1, V* = 500 and K_s' = 0.72 matching with the solution shown is Table 3.

Figure 10
Solution for a viscoplastic rockmass through .gui MatLab

6. Conclusion

A supported tunnel behaves as a three-dimensional structure for which strain and stress fields of the surrounding ground are strongly influenced. This paper develops an easy tool executable with MatLab to solve the problem of convergence and stress at equilibrium for elastic, elastoplastic and elastoviscoplastic tunnels.

It is important to know beforehand some geotechnical characteristics of the rockmass and the loading pattern to select the stress-strain criteria correctly. It is known that if the rockmass does not present permanent strain after being stressed it is going to behave as an elastic model, however the existence of creep and permanent deformations implies a plastic or a viscoplastic trend which give different solutions.

The problem was proved to be slightly different from the one solved by the Convergence-Confinement method, since the New Implicit Method highlights that the stiffness of the support is a great player on the convergence U₀ (at a distance d₀ from the tunnel face).

A point to be added is to take the permanent deformation of the supporting, for instance the concrete one, since it may change the curvature of the confinement curve.

The accuracy of the new method is satisfactory for a geotechnical preliminary study, regarding the numerical limitations. The average error stays at a maximum of 12%, with most of the solutions under 8% when compared to the numerical solution. Generally, most of the differences are being for the safety, since it is increasing the value of the stress at equilibrium which is imposed over the lining.

The results are good for values of K_s' greater than one. In addition to that, the results for smaller d₀' present better convergence than the ones of greater distance of lining placement. These smaller values represent most of the excavated tunnels, such as the value of d₀' ? 0,16 used in Tunnel Paraíso in Sao Paulo (Mafra, 2011, [⁷]). So, for conditions of K_s' > 1 and d₀' < 1, most of the solutions converge better to the numerical solution.

The numerical code GEOMEC91 presents an approximation of actual behavior. However, considering the validation of GEOMEC91 with clays shown before and that it is a realistic tool to preview the behavior of deep tunnels, it was used to evaluate the solutions given by the New Implicit Method encoded here in MatLab.

7. Acknowledgements

The authors recognize and would like to thank the support of the Federal University of Rio Grande do Sul as well as the scholarship provided by Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES), which incentives good researchers to keep the good knowledge and the improvement of technology in Brazil.

8. References

• ^[1] ASSOCIATION FRANÇAISE DES TUNNELS ET DE L’ESPACE SOUTERRAIN. Recommendations on the Convergence-Confinement Method. AFTES - 2001.
• ^[2] ASSOCIATION FRANÇAISE DES TUNNELS ET DE L’ESPACE SOUTERRAIN. Guidelines on diagnosis methods for lined tunnels. AFTES - 2005.
• ^[3] BERNAUD D.; BENAMAR, I. and ROUSSET, G. The “New Implicit Method” for tunnel in elastoplastic and viscoplastic fields. Revue Française de Géotechnique, n 68, 3-19, 1994.
• ^[4] BERNAUD, D. and ROUSSET, G. The “New Implicit Method” for tunnel analysis. International Journal for Numerical and Analytical Methods in Geomechanics, 20, pp. 673-690, 1996.
• ^[5] CARRANZA, C. and FAIRHURST, C. Application of the Convergence-Confinement Method of Tunnel Design to Rockmasses that satisfy the Hoek-Brown Failure Criterion. Tunneling and Underground Space Technology, vol 15, nº 2, pp. 187-213, 2000.
• ^[6] FIORE, P.V., MAGHOUS, D.B.; CAMPOS FILHO, A. A Tridimensional Finite Element approach to model a Tunnel with shotcrete and precast Concrete. Ibracon Structures and Materials Journal, vol. 9, n. 3, pp 403-413, 2016.
• ^[7] MAFRA, G. M. “Finite Element Method Analysis of Tunnel Excavated Rockmasses”. Master Thesis. University of São Paulo, Brazil, 2011.
• ^[8] PANET, M. “The calculus of tunnels through Convergence-Confinement Method”. Presses de l’école nationale des Ponts et Chauseés, Paris, 1995.
• ^[9] PARK, K. H. “Elastic Solution for Tunneling-Induced Ground Movements in Clays”. International Journal of Geomechanics, 4(4), 310-318, 2004.
• ^[10] QUEVEDO, F. P. M. “Long term behavior of Concrete Lined Deep Tunnels: Finite Element Model”. Master Thesis. Federal University of Rio Grande do Sul, Brazil, 2017.

^{9. APPENDIX I}

The most relevant parameters used in this paper are summarized in the following list. The parameters subscripted with (L used in E_L, ν_L and so on) are defined regarding the lining characteristics. The superscript (' used in Ks',d₀' and so on) indicates that such parameter is dimensionless.

• P_∞   Initial geostatic stress as defined in Equation (1);
• K_s   Stiffness of the lining as defined in Equations (2) and (3);
• d₀   Distance between the tunnel face and the lining during the excavation;
• U_o   Convergence of the tunnel at distance d₀; [dimensionless]
• U_f   Convergence of the tunnel at its face; [dimensionless]
• U_∞   Convergence of the unsupported tunnel very far from the face; [dimensionless]
• E   Modulus of Young;
• R_P   Plastic radius as Equation (25);
• C   Cohesion of the rockmass;
• ϕ   Friction angle;
• υ   Poisson coefficient;
• a,b,c   Parameters for the plastic solution system as Equations (23) and (34);
• H   Parameter on Coulomb Criterion;
• a^s _vp   Transformed function of unsupported to supported tunnel in viscoplastic field;
• R   Radius of the excavated tunnel;
• N_S   Parameter used in plasticity that indicates the level of plastification as Equations (19) and (27);
• U_eq   Convergence of the supported tunnel very far from the face; [dimensionless]
• P_eq   Pressure of the supported tunnel very far from the face;
• K_p   Coefficient of horizontal stresses at rest as Equation (30);
• a^s   Transformed function of unsupported to supported tunnel;
• V   Velocity of excavation;
• η   Viscosity constant;
• V*   Reduced velocity;
• P_eq   Pressure at equilibrium;
• R_c   Parameter on Coulomb Criterion;
• α*   Average solution which includes the influence of the lining strength on the convergence;
• A,B   Parameters for the viscoplastic solution as Equations (38) and (39);
• C₁, C₂   Parameters for the viscoplastic average solution as Equations (44) and (45).

The other remaining parameters are listed on work presented before. They are referred as soon as they show up, and the list above are with the most recurrent and important parameters to take note.

Publication Dates

• Publication in this collection
  May 2018

History

• Received
  03 Apr 2017
• Accepted
  10 Oct 2017