Abstract
Excavation of tunnels produces a redistribuition of stresses and induces deformations in the rock mass around the tunnel’s cross section. In the case of elasto-plastic behavior of rock mass, plastic zones may appear. It is important to quantify the influence of this zone on the overall response of the tunnel. In this paper, we deduce a fully analytical solution in terms of displacements and stresses around a circular deep tunnel. The aim here is not to replace a 3D numerical calculation. This kind of analytical calculation are only useful to have a good understanding of the tunnel behavior in the preliminary phases of the project. For example, to perform parametric studies useful to choosing good parameters to introduce in a 3D numerical calculation. A homogeneous and isotropic rock mass is considered. For elasto-plastic behavior, the Tresca’s constitutive model with associate flow rule and Mohr-Coulomb’s constitutive model with non-associate flow rule are considered. For both, the idealized stress-strain curve presents a linear istropic hardening law. A geostatic-hydrostatic state of initial stresses and infinitesimal strains is assumed. The analytical solutions are compared with the FEM solutions demonstrating excellent agreement.
Keywords:
Deep tunnels; analytical solution; constitutive model; elastoplasticity
1 INTRODUCTION
During the analysis of a deep tunnel it is essential to understand the behavior of the rock mass after the excavation. The excavation induces a disturbance zone around the tunnel that depends on the properties of the rock mass, the excavation method, in situ stress field, the geometry of the cross section of the tunnel and the type and stiffness of the support and the timing of its installation. This plastic zone is defined as a region in the rock mass where its stiffness parameters have changed due to the process of excavation and support of the tunnel. One of the tools commonly used to understand the interdependence of these factors is analytical solutions.
The analytical solutions consist of the use of equations deduced from the assumptions of continuous mechanics. However, for the closed-form solution, it is necessary to adopt some simplifying hypotheses, such as circular section, infinite homogeneous medium, axissimetry boundary conditions and geostatic-hydrostatic stress field. In any case, these solutions have the potential of providing a quick overview of the problem and are very useful to identify the relationship between the most important variables. The aim here is not to replace a three-dimensional numerical calculation. This kind of analytical calculation are only useful to have a good understanding of the tunnel behavior in the preliminary phases of the project. In addition, as they are quick calculations, they allow parametric studies that can be useful when choosing the parameters of a three-dimensional numerical calculation that takes more time.
The first analytical solution employed in deep tunnels derives from theoretical development to determine the field of stresses and deformations around openings in elastic media. Lamé (1866Lamé, G. (1866). Leçons sur la théorie mathématique de l’elasticité des corps solides par G. Lamé. Gauthier-Villars.) proposed the first solution (in plane state of deformations) for cylindrical openings in linear elastic medium submitted to an initial state of hydrostatic geostatic stresses. Afterwards, Kirsch (1898Kirsch, C. (1898). “Die theorie der elastizitat und die bedurfnisse der festigkeitslehre.” Zeitschrift des Vereines Deutscher Ingenieure, 42, 797-807.) proposed a solution in elastic media generalizing the initial stresses field. In the following century, along with the development of the phenomenological theory of plasticity, several analytical solutions were proposed. In general, considering different yield criteria (Mohr-Coulomb, Hoek-Brown, Tresca for example), stress-strain behavior (elasto-plastic perfect, hardening/softening or brittle) and strategies for considering volumetric reduction during plastic deformation (such as using a non-associated flow rule). Studies such as Fenner (1938Fenner, R. (1938). Untersuchungen zur erkenntnis des gebirgsdrucks. Glückauf, 74, 681-695 and 705-715.), Morrison and Coates (1955Morrison, R. G. K., & Coates, D. F. (1955). Soil mechanics applied to rock failure in mines. The Canadian Mining and Metallurgical Bulletin, 48(523), 701-711.), Salencon (1969Salencon, J. (1969). Contraction quasi-statique d’une cavite a symetrie spherique ou cylindrique dans un milieu elastoplastique. In Annales des ponts et chaussées (Vol. 4, pp. 231-236).), Peck (1969Peck, R. B. (1969). Deep excavations and tunneling in soft ground. Proc. 7th ICSMFE, 1969, 225-290.), Daemen and Fairhurst (1970Daemen, J. J. K., & Fairhurst, C. (1970). Influence of failed rock properties on tunnel stability. In The 12th US Symposium on Rock Mechanics (USRMS). American Rock Mechanics Association.), Panet (1976Panet, M. (1976). Analyse de la stabilité d'un tunnel creusé dans un massif rocheux en tenant compte du comportement après la rupture. Rock mechanics, 8(4), 209-223.), Berest et al. (1978Berest, P., Minh, D. N., & Panet, M. (1978). Contribution a l’étude de la stabilité d’une cavité souterraine dans un milieu avec radoucissement. Revue Française de Géotechnique, (4), 63-72), Detournay and Fairhurst (1987Detournay, E., & Fairhurst, C. (1987). Two-dimensional elastoplastic analysis of a long, cylindrical cavity under non-hydrostatic loading. In International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts (Vol. 24, No. 4, pp. 197-211). Pergamon.), Corbetta (1990Corbetta, F. (1990). Calculs Analytiques et Numériques de Tunnels Profonds (Doctoral dissertation, Ecole Nationale Supérieure des Mines de Paris)) are exemples of this contribution. Brown et al. (1983) present a good review of these solutions until 1980s and presents its own analytical solution.
In recent analytical solutions, part of the researchers is studying the influence of a non-hidrostatic initial stress field considering the stress along the axis of the tunnel the major principal. Studies such as Reed (1988Reed, M. B. (1988). The influence of out‐of‐plane stress on a plane strain problem in rock mechanics. International Journal for Numerical and Analytical Methods in Geomechanics, 12(2), 173-181.), Pan and Brown (1996Pan, X. D., & Brown, E. T. (1996). Influence of axial stress and dilatancy on rock tunnel stability. Journal of Geotechnical Engineering, 122(2), 139-146.), Lu et al. (2010Lu, A. Z., Xu, G. S., Sun, F., & Sun, W. Q. (2010). Elasto-plastic analysis of a circular tunnel including the effect of the axial in situ stress. International Journal of Rock Mechanics and Mining Sciences, 47(1), 50-59.), Zhou and Li (2011Zhou, X. P., & Li, J. L. (2011). Hoek-Brown criterion applied to circular tunnel using elastoplasticity and in situ axial stress. Theoretical and applied fracture mechanics, 56(2), 95-103.), Wang et al. (2012Wang, S., Wu, Z., Guo, M., & Ge, X. (2012). Theoretical solutions of a circular tunnel with the influence of axial in situ stress in elastic-brittle-plastic rock. Tunnelling and Underground Space Technology, 30, 155-168.), including finite strain, such Vrakas and Anagnostou (2014Vrakas, A., & Anagnostou, G. (2014). A finite strain closed‐form solution for the elastoplastic ground response curve in tunnelling. International Journal for Numerical and Analytical Methods in Geomechanics, 38(11), 1131-1148.) are examples of this contribution.
This article presents two analytical solution: the first considers a Tresca’s constitutive model with linear hardening behavior with associated flow rule. The second case, more complex, considers a Mohr-Coulomb’s constitutive model with linear hardening behavior with non-associated flow rule, in wich post-peak dilatancy occurs at a constant rate with major principal strain. In both cases, an initial geostatic-hydrostatic stress field and infinitesimal deformation is adopted.
2 DEFINITION AND HYPOTHESIS OF THE PROBLEM
The Figure 1 shows the geometry and boundary conditions of the problem. A circular deep tunnel of radius in a homogeneous and isotropic rock mass with an elastic-plastic hardening behavior. A plane state of deformation is assumed, i.e., the section is far from the three-dimensional effects of the excavation face, toghether with axial symmetry boundary conditions for uniform loading . is the pressure on the wall of the tunnel that simulates the decompression caused by the advance of the excavation and placement of the lining. For example, it is used for generate the characteristic curve of rock mass in the Convergence-Confinement method (Panet et al., 2001Panet, M., Givet, P. D. C. O., Guilloux, A., Duc, J. L. D. G. N. M., Piraud, J., & Wong, H. T. S. D. H. (2001). The convergence-confinement method. AFTES-recommendations des Groupes de Travait.).
When decreases from and the yield criterion is achived, a cylindrical plastic zone begins to develop around the tunnel. In an elastic-plastic rock mass with hardening, this zone can have two internal domains: a hardening elastic-plastic zone and a perfect plastic zone. In the absence of lining, when , depending on the properties of rock mass, the advancement of this plastic zone can lead to high deformations, indicating the possible failure of the rock mass.
3 TRESCA CRITERION
3.1 Equations of the problem
Considering the axial symmetry condition of the problem in Figure 1, the resulting stress state at a distance is defined by the radial stress and circumferential stress, which are the minor and major principal stresses, and respectively, when the section of interest is far from the excavation face and with initial state of geostatic-hydrostatic stresses. The resulting displacement is defined by the radial displacement .
The strain tensor in polar coordinates is
and the stress tensor is written as
Applying the condition , i. e., the problem is statically admissible, have the following differential equation:
Considering the boundary conditions of Figure 1 which are
In linear elasticity behaviour, where (Hooke’s law), the solution of this boundary value problem is
In infinitesimal plasticity the strain deformation is the sum of the tensor of elastic deformations and the tensor of the plastic deformation , i.e.,
This deformation will occur when the stresses reach the Tresca yield criterion written as
It should be noted that the out-of-plane stress is taken as the intermediate principal stress . This assumption is valid with geostatic-hydrostatic in situ stress conditions and sections far from the excavation face. In this solution the plastic deformation was assumed to be associated, i.e.,, and given by the following expression
Here we are interested in the particular dependence of the cohesion , shown in Figure 2, which introduces the linear hardening law.
Variation of the cohesion: A -linear hardening behavior and B - perfect plasticity behavior
And therefore,
where
and the internal variable is the plastic strain magnitude .
Introducing (12) in (10) and replacing in Hook's law gives
Making the addition of (15) and (16):
So, using the (17) in (3) we obtain two equations valid in all the domain: (18) and (19) that we will use in the following deductions.
Replacing (18) in (19) the last one becomes:
where
Considering a monotonic decreasing function in starting from , we have the different behavior zones in the rock mass, show in Figure 3. A first phase, when remains close enough to the , the entire rock mass remains in elasticity, having only zone 1. However, when is below a limit value, necessary for the stress state to reach the yield criterion at , the plastic zone 2 with linear hardening behavior begins. When decreases further, a third zone, of perfect plasticity, may appear.
3.2 Solution in the elastic zone 1
The solution in this zone is classic (Mandel, 1966Mandel, J. (1966). Cours de mécanique des milieux continus. Cours de mécanique des milieu continus. Gauthier-Villars.). We define the radius those where the elastic zone begins. The conditions and in allows obtaining the value of the constant using (18) and (19):
Using the condition in the elastic domain and the boundary condition in we obtain the complete solution in this zone:
3.3 Solution in the plastic zone 3
We define the radius as the limit between the plastics zones 3 and 2. In this perfect plastic zone, the yield criterion is given by
As far as the yield criterion is zero in this zone, the equilibrium equation (3) becomes
The integration of the (19) together with the boundary condition of (4) leads to the following equation of the radial stress:
Replacing (30) into (18) we obtain the expression of the radial displacement given by
Afterwards, using the equilibrium equation and the plane state condition of deformations, the rest of the solution is given by
3.4 Solution in the plastic zone 2
In this zone of plasticity with hardening behaviour the yield criterion is written as
As far as this criterion is null in this zone, the integration of the equilibrium equation with the condition in leads to
where is the value of the radial stress in . Due to the continuity of the radial stress between zone 3 and 2, equation (30) provides
In the same way that for the zone 3, the expression of the displacement is obtained using (18), which gives:
Using the equilibrium equation and the plane state condition of deformations, the rest of the solution is given by:
3.5 Calculation of the plastic radius and
In order to finish the resolution of the problem we must calculate the plastic radius and .In we have
Substituting (40) in (20) we obtain
The value of is finally obtained writing the continuity of the radial stress in , i.e.,
Together with (36) and (41) the value of is obtain in an explicit way:
with the following constants:
3.6 Comparison between the numerical and analytical solution
The numerical solution was obtained through the elastoplasticity algorithm implemented in the GEOMEC 91 finite element code. The local elastoplasticity integration algorithm of the stresses and internal variables and the global equilibrium iteration algorithm can be obtained in Bernaud (1991Bernaud, D. (1991). Tunnels profonds dans les milieux viscoplastiques: approches expérimentale et numérique (Doctoral dissertation, Ecole Nationale des Ponts et Chaussées).). The finite element mesh with 100 isoparametric elements with 9 nodes per element, together with their boundary conditions is shown in Figure 4.
The characteristics of the following example are those of the Boom clay (Bernaud, 1991Bernaud, D. (1991). Tunnels profonds dans les milieux viscoplastiques: approches expérimentale et numérique (Doctoral dissertation, Ecole Nationale des Ponts et Chaussées).) obtained from samples at 240m deep: , , , , , and . Figure 5 illustrates the convergence curve, , for a given application. The point marked in Figure 5 corresponds in FEM numerical calculation and we can see that numerical and analytical solutions are identical.
The Figure 6 shows the plastic radius for the internal pressure range.
For a given pressure the Figure 7 illustrates the stress versus the radius for analytical and numerical solutions.
As before, we can show the excellent agreement between the analytical and numerical solution. In Figure 7, through , the three zones that delimit the elastoplastic regime can be clearly seen.
4 MOHR-COULOMB CRITERION
4.1 Equations of the problem
In this case, the same strain and stresses tensor given by (1) and (2), equilibrium equation (3), the boundary conditions (4) and (5) and the superposition of the deformations elastic and plastic (10) are used. But the yield function is written as
where
and is the friction angle. In this model the hardening behaviour was adopted only in the cohesive parameter . Assuming a non-associated flow rule, i.e., , where
and
obtained for plastic deformation
where is the material dilatation angle. The internal variable associated with the isotropic hardening behavior is given by:
Here we are interested in the particular dependence of the in (55), shown in Figure 8, which introduces the linear hardening law.
Variation of the function: A - linear hardening behaviour and B - perfect plasticity behavior
And therefore,
where
The equilibrium equation (3) and the Hooke’s law give the following relations (that are valid in all zones of behavior)
Considering a monotonic decreasing function in starting from , we have the different behavior zones in the rock mass. A first phase, when remains close enough to the , the entire rock mass remains in elasticity, having only zone 1. However, when is below a limit value, necessary for the stress state to reach the yield criterion at , the plastic zone 2 with linear hardening behavior begins. When decreases further, a third zone, of perfect plasticity, may appear.
In the following, we describe the deductions corresponding to the two cases of plasticity. For the first, we have only the plastic with linear hardening zone. The second is more complex because exist two plastics zones in the rock mass: the zone with perfect plasticity and the zone of plasticity with hardening behaviour.
4.2 First case: only one plastic zone
In this case exists only the plastic with hardening zone corresponding to the zone 2 illustrates in the Figure 9.
4.2.1 Solution in the elastic zone 1
The solution in this zone is obtained writing the conditions:
These conditions applied to the (58) and (59) gives the following solution:
4.2.2 Solution in the plastic zone 2
Considering the continuity of and in :
the displacements are obtained by integration of (59) between and . Expressed as a function of we have
with
Due to the nullity of the yield criterion in this zone, the equilibrium equation writes
Due to the fundamental (58) we can eliminate . So we obtain
Equation where the second member can be written as a function of , given in the (63). So, an integral-differential equation in is expressed as
The method of resolution (74) uses technics purely mathematics. Therefore, we expose in the following directly the solution obtained:
with
The constants and are obtained taken into account the boundary and continuity conditions of stresses, given by:
Therefore:
Finally, the displacement is a function only of the variable
with the follow constants:
4.2.3 Calculation of the plastic radius
In order to complete the solution of the problem, we must calculate the plastic radius . This is performed writing the nullity of and of the yield criterion in in the (58):
The deduction can be performed using the (67) and (80) and obtaining an implicit equation of :
with
This solution is only valid as long as the stress values maintain the condtion .
4.3 Second case: two plastics zones
In this case, the second zone (perfect plasticity) appears in (see Figure 10).
4.3.1 Solution in the plastic zone 3
In this zone the yield criterion writes:
Due to the nullity of in this zone, the equilibrium equation is:
The resolution of the differential equation above with the limit condition gives the classical solution:
The calculation of the displacements is performed by the integration of the fundamental (59) between the radius and :
Where representes the displacement in the point . Due to the continuity of the displacements can be calculate by the equation of in the second plastic zone.
4.3.2 Solution in the plastic zone 2
In this zone, the stresses and displacements are written as the same way that in the first calculation, (67) and (80). The difference concerns the limit condition of the first calculation that is replaced by the condition . So, the radial stress in this zone is written with the (80) but with others values of and that are:
All the others constants are the same as of the first calculation. We remember:
The displacements in this zone are given by the (90) with the constants of (91) to (95). The constants and must be calculated with the (106) and (107). We remember this equation:
4.3.3 Calculation of the plastic radius and
As before, we need to calculate the plastic radius and to complete the solution of the problem. With this goal we write the conditions about in and :
Considering the (58) the first condition is
The derivation in function of for the (67) written in and the value of allow to write the above equation as a function of and :
where and the constants:
The second condition ( in ) gives:
The derivation of the (67) is written in and we define in the (118). Finally, we can express as a function of .
with the following constants:
So, the (113) and (119) are a system of two implicit equations for and that must be solve numerically. First we replace the expression of in function of given by the (119) in the (113). With a numerical calculation (calculation of the zero of a function) we find the value of for a given . This value introduced in the (119) gives the value of and consequently the value of .
4.4 Comparison between the numerical and analytical solution
The numerical solution uses the same program as the the previous solution (GEOMEC 91) whose mesh and boundary conditions are shown in Figure 4. The characteristics of the following example are those of the Boom clay (Bernaud, 1991Bernaud, D. (1991). Tunnels profonds dans les milieux viscoplastiques: approches expérimentale et numérique (Doctoral dissertation, Ecole Nationale des Ponts et Chaussées).) obtained from samples at 240m deep: , , , , , , , and . Figure 11 illustrates the convergence curve for a given application.
The Figure 12 shows the plastic radius for spectrum of internal pressures.
For a given pressure (with ) Figure 13 illustrates the stress versus the radius r for analytical and numerical solutions.
Stress versus the radius for analytical and numerical solutions for the Mohr-Coulomb criterion
5 ANALYTICAL AND NUMERICAL CONVERGENCES - COMPARISON
Table 1 summarizes the result of analytical and numerical convergences of six different cases. Two considering the Tresca model (T) and five considering the Mohr-Coulomb model (MC). The characteristics of rock mass correspond to that of Boom’s clay (Bernaud, 1991Bernaud, D. (1991). Tunnels profonds dans les milieux viscoplastiques: approches expérimentale et numérique (Doctoral dissertation, Ecole Nationale des Ponts et Chaussées).).
6 CONCLUSIONS
In this paper, simple closed-form solutions of the displacements and stresses for an elastic-plastic rock mass of the circular deep tunnel with Tresca’s and Mohr-Coulomb’s linear hardening model was presented. Three zones appear around the tunnel: an elastic zone, a plastic hardening zone and a perfect plastic zone. The analytical solution for these three zones is shown. The analytical results were compared with the FEM results and a perfectly agreement was obtained.
References
- Berest, P., Minh, D. N., & Panet, M. (1978). Contribution a l’étude de la stabilité d’une cavité souterraine dans un milieu avec radoucissement. Revue Française de Géotechnique, (4), 63-72
- Bernaud, D. (1991). Tunnels profonds dans les milieux viscoplastiques: approches expérimentale et numérique (Doctoral dissertation, Ecole Nationale des Ponts et Chaussées).
- Brown, E. T., Bray, J. W., Ladanyi, B., & Hoek, E. (1983). Ground response curves for rock tunnels. Journal of Geotechnical Engineering, 109(1), 15-39.
- Corbetta, F. (1990). Calculs Analytiques et Numériques de Tunnels Profonds (Doctoral dissertation, Ecole Nationale Supérieure des Mines de Paris)
- Daemen, J. J. K., & Fairhurst, C. (1970). Influence of failed rock properties on tunnel stability. In The 12th US Symposium on Rock Mechanics (USRMS). American Rock Mechanics Association.
- Detournay, E., & Fairhurst, C. (1987). Two-dimensional elastoplastic analysis of a long, cylindrical cavity under non-hydrostatic loading. In International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts (Vol. 24, No. 4, pp. 197-211). Pergamon.
- Fenner, R. (1938). Untersuchungen zur erkenntnis des gebirgsdrucks. Glückauf, 74, 681-695 and 705-715.
- Kirsch, C. (1898). “Die theorie der elastizitat und die bedurfnisse der festigkeitslehre.” Zeitschrift des Vereines Deutscher Ingenieure, 42, 797-807.
- Lamé, G. (1866). Leçons sur la théorie mathématique de l’elasticité des corps solides par G. Lamé. Gauthier-Villars.
- Lu, A. Z., Xu, G. S., Sun, F., & Sun, W. Q. (2010). Elasto-plastic analysis of a circular tunnel including the effect of the axial in situ stress. International Journal of Rock Mechanics and Mining Sciences, 47(1), 50-59.
- Mandel, J. (1966). Cours de mécanique des milieux continus. Cours de mécanique des milieu continus. Gauthier-Villars.
- Morrison, R. G. K., & Coates, D. F. (1955). Soil mechanics applied to rock failure in mines. The Canadian Mining and Metallurgical Bulletin, 48(523), 701-711.
- Pan, X. D., & Brown, E. T. (1996). Influence of axial stress and dilatancy on rock tunnel stability. Journal of Geotechnical Engineering, 122(2), 139-146.
- Panet, M. (1976). Analyse de la stabilité d'un tunnel creusé dans un massif rocheux en tenant compte du comportement après la rupture. Rock mechanics, 8(4), 209-223.
- Panet, M., Givet, P. D. C. O., Guilloux, A., Duc, J. L. D. G. N. M., Piraud, J., & Wong, H. T. S. D. H. (2001). The convergence-confinement method. AFTES-recommendations des Groupes de Travait.
- Peck, R. B. (1969). Deep excavations and tunneling in soft ground. Proc. 7th ICSMFE, 1969, 225-290.
- Reed, M. B. (1988). The influence of out‐of‐plane stress on a plane strain problem in rock mechanics. International Journal for Numerical and Analytical Methods in Geomechanics, 12(2), 173-181.
- Salencon, J. (1969). Contraction quasi-statique d’une cavite a symetrie spherique ou cylindrique dans un milieu elastoplastique. In Annales des ponts et chaussées (Vol. 4, pp. 231-236).
- Vrakas, A., & Anagnostou, G. (2014). A finite strain closed‐form solution for the elastoplastic ground response curve in tunnelling. International Journal for Numerical and Analytical Methods in Geomechanics, 38(11), 1131-1148.
- Wang, S., Wu, Z., Guo, M., & Ge, X. (2012). Theoretical solutions of a circular tunnel with the influence of axial in situ stress in elastic-brittle-plastic rock. Tunnelling and Underground Space Technology, 30, 155-168.
- Zhou, X. P., & Li, J. L. (2011). Hoek-Brown criterion applied to circular tunnel using elastoplasticity and in situ axial stress. Theoretical and applied fracture mechanics, 56(2), 95-103.
Edited by
Editor:
Publication Dates
-
Publication in this collection
18 Sept 2020 -
Date of issue
2020
History
-
Received
12 Mar 2020 -
Reviewed
30 June 2020 -
Accepted
31 July 2020 -
Published
05 Aug 2020