J. Cent. South Univ. Technol. (2010) 17: 1357-1363

DOI: 10.1007/s11771-010-0642-x

3D FEM analysis for layered rock considering anisotropy of shear strength

ZHANG Yu-jun(张玉军), ZHANG Wei-qing(张维庆)

Institute of Rock and Soil Mechanics, Chinese Academy of Sciences, Wuhan 430071, China

? Central South University Press and Springer-Verlag Berlin Heidelberg 2010

Abstract: An empirical expression of cohesion (C) and friction angle  for layered rock was suggested. This expression was compared with a test result made by the former researchers. The constitutive relationship of a transversely isotropic medium and Mohr-Coulomb criterion in which C and  vary with directions were employed, and a relative 3D elasto-plastic FEM code was developed, in which the important thing was to adopt a search-trial method to find the orientation angle (ρ) of shear failure plane (or weakest shear plane) with respect to the major principal stress as well as the corresponding C and  Taking an underground opening as the calculation object, the numerical analyses were carried out by using the FEM code for two cases of transversely isotropic rock and isotropic rock, respectively, and the computation results were compared. The results show that when the rock is a transversely isotropic one, the distributions of displacements, plastic zones and stress contours in the surrounding rock will be non-axisymmetric along the tunnel’s vertical axis, which is very different from that of isotropic rock. The stability of the tunnel in transversely isotropic rock is relatively low.

Key words: layered rock mass; shear strength; anisotropy; three dimensional finite element method; computation analysis

1 Introduction

It is well-known that a layered rock possesses obvious anisotropic behavior of deformation and  strength, compared with a homogeneous rock, its conditions of stability and failure are more complicated [1]. At present a layered rock can be treated theoretically as a transversely isotropic medium, its stress-strain relationship includes five constants, and the strength parameters also exhibit some spatial distribution modes. In the past many years the scholars at home and abroad have done constantly researches on anisotropic shear strength of layered rock masses. For instance, LEE and PIETRUSZCZAK [2] assumed that the crack initiation of transversely isotropic rock masses is governed by the Hoek-Brown failure criterion, and described the anisotropy of the strength through the orientation dependence of strength parameters m and s. SAROGLOU and TSIAMBAOS [3] modified the Hoek-Brown failure criterion by incorporating a new parameter (k_b) to account for the effect of strength anisotropy. XU and ZHANG [4] obtained the orthotropic equivalent shear strength parameters for fractured rock mass by applying the simulation method of joint network, damage mechanics and theory of joint tensor. YANG et al [5] established an implicit anisotropic yield criterion for jointed rock masses through the Mohr-Coulomb failure condition for plane in every direction. YU [6] set up a strength criterion that can describe the inherent anisotropy as well as two failure modes of cracking and yielding for geotechnical material. COLAK and UNLU [7] undertook a study to determine the effect of transverse anisotropy of sedimentary rocks on the Hoek-Brown strength parameter (m_i). LYDZBA et al [8] considered a laminate consisting of interchanging layers of sandstone and claystone to be a particular formation and dealt with the specification of transversely isotropic failure criterion for stratified rocks. PIETRUSZCZAK  et al [9] described the inherent anisotropy of sedimentary rocks by employing a second-order microstructure  tensor, in which the distribution function specified the directional dependence of strength parameters. USHAKSARAEI and PIETRUSZCZAK [10] presented a macroscopic transversely isotropic failure criterion for structural masonry based on the critical plane approach. PIETRUSZCZAK and MROZ [11] formulated the failure criteria for the cohesive-frictional materials using two different approaches.

The above-mentioned scholars inspected and verified the anisotropic failure criteria suggested by themselves using laboratory tests and/or computation analyses, but they did not apply the criteria to FEM and other numerical simulations with practical background of  geotechnical engineering.

There should be τ=C(θ)+σtan φ(θ) if Mohr-Coulomb criterion is employed to describe the failure characteristics of a layered rock. That is, both cohesion (C) and internal friction angleare not constants, but are functions of the angle between the failure surface and the foliated plane (θ) (or the direction of maximum principal stress σ₁) [12-13]. For this reason, an empirical expression of C and  for a layered rock was suggested. On the basis of this constitutive relationship of a transversely isotropic medium and Mohr-Coulomb criterion in which C and  vary with directions were introduced. Therefore a relative 3D elasto-plastic FEM code was developed. Taking an underground opening as the calculation object, the numerical analyses were carried out for both cases of a transversely isotropic rock and an isotropic rock, respectively.

2 Anisotropic shear strength of layered rock

2.1 Expression of C and  for layered rock

Generally speaking, the shear strength of a layered rock (S) has a minimum value of S_min and a maximum value of S_max in directions parallel to and perpendicular to the foliated plane, respectively, and among them S varies with angle θ with respect to the foliated plane according to some rules (Fig.1) as:

S=f(S_min, S_max, θ)                              (1)

Fig.1 Inclination of plane with C andconcerned θ with respect to foliated plane

JAEGER [14] also gave an expression of S:

S=S₁+S₂cos[2(α-β)]                           (2)

where S₁ and S₂ are constants; α and β are angles between the failure surface as well as the foliated plane and the direction of maximum principal stress σ₁, respectively (see Fig.2).

As a kind of approximation, an empirical expression for C and  in any direction with θ was suggested [15]:

                   (3)

where 0≤θ≤π/2; C_min,  C_max and are the cohesion and the internal friction angles parallel to and perpendicular to the foliated plane, separately.

Fig.2 Inclination α and β of failure plane and foliate plane with respect to σ₁

2.2 Comparison with earlier test result

Here, the test results made by ATTEWELL and SANDFORD [16] were used. For Penrhyn slate through changing β of the foliated plane with respect to the axis of maximum principal stress, ATTEWELL and SANDFORD [16] obtained C and  values for different β magnitudes under the condition of triaxial compression. A curve-fit formula of C and  acquired by ATTEWELL and SANDFORD [16] is expressed as

                     (4)

where C₁, C₂,  and  are constants; γ and δ are inclinations of the failure surface and the foliated plane with respect to σ₁.

ATTEWELL and SANDFORD [16] did not give the concrete values of C₁, C₂, γ,and δ. In this work, a curve-fitting was made on the basis of results in Ref.[16] (directional influence curves of cleavage anisotropy upon C (Pa) and (?) at failure of the rock specimens) and produce:

                  (5)

From Eq.(5) we have

Thus

                    (6)

Let β=90?, that is, the foliated plane is perpendicular to the axis of σ₁ and the change of γ, δ and θ in Eqs.(4), (6) from 0? to 90?, the relative curves for C and  are demonstrated in Fig.3. From Fig.3 it can be seen that Eqs.(4) and (6) have a good agreement.

3 Formulae of 3D FEM

3.1 Constitutive relationship

For the transversely isotropic rock shown in Fig.4,

Fig.3 Relationship between shear strength and θ: (a) Cohesion C; (b) Friction angle

Fig.4 Local and global coordinate systems

the foliated plane is taken as x′y′ plane, and the normal direction of the foliated plane is taken as z′ axis, then the elastic matrix in the local coordinate system x′y′z′ is written as [17].

         (7)

where M=E₁/(1+μ₁)(1-μ₁-2nμ₂²); n=E₂/E₁; G₁=E₁/ 2(1+μ₁); E₁ and μ₁ are the elastic modulus and Poisson ratio in the foliated plane, respectively; and E₂, G₂ and μ₂ are the elastic modulus, the shear modulus   and Poisson ratio in the direction perpendicular to the foliated plane, respectively.

Stresses and strains are transformed from the local coordinate system x′y′z′ into the global coordinate system xyz, we get

                                    (8)

where σ, ε and D are stress tensor, strain tensor and elastic matrix in the global coordinate system, respectively, and

;

L=

;

L is the coordinate transformation matrix; l_i, m_i and n_i (i=1, 2, 3) are cosines of the angles between the axes of the local and global coordinate systems.

3.2 Failure criterion

For the layered rock and stress state presented in Fig.5, if Mohr-Coulomb criterion is applied, a shear surface will occur in σ₁-σ₃ plane, we have

                          (9)

and

                       (10)

            (11)

where τ and σ are shear and normal stresses, respectively; σ₁ and σ₃ are major and minor principal stresses, respectively; and ρ is angle between the failure surface and σ₃ direction.

Fig.5 Stress state and foliated planes

Substituting Eq.(10)-(11) into Eq.(9) leads to a yield function as:

     (12)

After a transformation, Eq.(12) becomes [18]

        (13)

and

where J₁ and J₃ are the first and third stress invariants, respectively; J′₂ is the second invariant of the deviatoric stresses, σ_x, σ_y and σ_z are normal stresses in x, y and z directions, respectively; σ′_x, σ′_y and σ′_z are deviatoric stresses in x, y and z directions, respectively; τ_xy, τ_yz and τ_xz are shear stresses in xy, yz and xz planes, respectively; and η is Lode angle.

When angle ρ between the failure surface and σ₃ axis and the corresponding C(θ) and  are known, the plastic matrix can be computed as

                      (14)

where H is a hardening function.

It is necessary to point out that if the rock is an isotropic one, its C and  are constants, so there is ρ=π/4+ when the failure happens [19]; but here the rock is a transversely isotropic one, C and  are functions of the orientations with respect to the foliated plane, so angle ρ and the corresponding C and  cannot be determined explicitly when the rock breaks. For this reason a search-trial method was put forward as follows.

In the local coordinate system shown in Fig.6, it is assumed that a shear failure surface is produced, and its normal cosines are l_x, m_x and n_x. As the normal is in σ₁-σ₃ plane, the angles between the normal and axes are T₁, T₂ and 0, respectively, so this results in:

                      (15)

where A=cos T₁ and B=cos T₂.

Fig.6 Layered rock and x′y′z′ coordinate system

                      (16)

where   b₂=

and n_x is the cosine of angle between the

normal of the shear failure surface and z′ axis (rotationally elastic symmetry axis), so it leads to

                                (17)

Then, according to the relationship among C,  and θ in Eq.(3), C and  may be calculated. But usually, {l_x, m_x, n_x} are unknown beforehand. For this purpose a search-trial method is employed, that is, let angle ρ change from 0? to 90?, and for every ρ_f, A_i=cos T_1i and B_i=cos T_2i are known, so n_xi and θ_i can be determined, and C_i and  are obtained.

Using  to find angle ρ_f corresponding to the maximum value of f_i, the surface with ρ_f is the real or potential shear failure surface to be evaluated. This work can be done by the computer easily.

4 Calculation examples

It is assumed that the longitudinal axis of a tunnel, whose orientation is β₀=0?, is z axis of a global coordinate system, and the strike and dipping angles of the foliated plane of a layered rock around the tunnel are β=45? and α=45?, respectively, so we have [20]:

The cross-section of the tunnel is a rectangular one with a size of 10 m×10 m. The thickness of the overburden is 20 m. For the cases of taking the rock as a transversely isotropic one and an isotropic one respectively, the numerical analyses by using the 3D elasto-plastic FEM code developed by the authors were carried out. The mesh for FEM and the mechanical parameters of rock are shown in Fig.7 and Table 1, respectively.

Fig.7 FEM mesh of rock

The displacement, plastic zones and contours of stresses in the plane perpendicular to z axis σ_x, σ_y obtained from the calculations are presented in Table 2 and Figs.8-9. From these data it can be seen that when the rock is a transversely isotropic one, the distributions of the displacements, plastic zones, σ_x and σ_y contours in the surrounding rock are non-axisymmetric along the tunnel’s vertical axis, and the quantity values are larger. For example, the range of resultant displacements of the tunnel periphery is 4.7-10.9 cm, and there are 148 plastic elements in the surrounding rock at the computation end. However, when the rock is an isotropic one, the distributions of the displacements, plastic zones, σ_x and σ_y contours in the surrounding rock are axisymmetric along the tunnel’s vertical axis, and the quantity values are smaller. For example, the range of resultant displacements of the tunnel periphery is 2.3- 6.7 cm, and there are 78 plastic elements in the surrounding rock at the computation end. There are distinct differences between the two cases, and the stability of the tunnel in transversely isotropic rock is relatively low.

Table 1 Mechanical parameters of rock

Table 2 Displacements at some points on opening boundary (Unit: cm)

Fig.8 Calculation case for transversely isotropic rock: (a) Displacement vectors; (b) Plastic zones; (c) Contours of σ_x (10^-1 MPa);    (d) Contours of σ_y (10^-1 MPa)

Fig.9 Calculation case for isotropic rock: (a) Displacement vectors; (b) Plastic zones; (c) Contours of σ_x (10^-1 MPa); (d) Contours of σ_y (10^-1 MPa)

So, it can be known that when a layered rock has very different mechanical properties in the directions parallel to and perpendicular to the foliated plane separately, and its occurrence shows certain geometrical relationship with the excavation factors (for example, the tunnel’s axial line), the dynamic state of the surrounding rock will be more complicated compared with the case of isotropic rock. So if a layered rock is met in the practice of rock mechanics and engineering, it is necessary to consider the anisotropy of mechanical behavior of the rock.

5 Conclusions

(1) The empirical expression of C and  for layered rock is proposed and compared with a test result made by former researchers. The results show that both have a good agreement, indicating that the empirical expression is reliable.

(2) The principles of 3D elasto-plastic FEM analysis are described, in which a layered rock is taken as a transversely isotropic one and the suggested empirical expression of C and  is used. Among them the important thing is to adopt a search-trial method to find orientation angle ρ of shear failure plane (or weakest shear plane) with respect to the major principal stress σ₁ as well as C and

(3) The calculation examples show that when a layered rock possesses obvious anisotropic mechanical behavior and its occurrence presents certain geometrical relationship with the excavation factors (for example, a tunnel’s axial line), the dynamic state of surrounding rock is more completed compared with the case of isotropic rock. So, this feature should be considered in numerical analysis and construction of a project.

References

[1] SERRANO A, OLALLA C. Ultimate bearing capacity of an anisotropic discontinuous rock mass. Part I: Basic modes of failure [J]. Int J Rock Mech Min Sci, 1998, 35(3): 301-324.

[2] LEE Y K, PIETRUSZCZAK S. Application of critical plane approach to the prediction of strength anisotropy in transversely isotropic rock masses [J]. Int J Rock Mech Min Sci, 2008, 45(3): 513-523.

[3] SAROGLOU H, TSIAMBAOS G. A modified Hoek-Brown failure criterion for anisotropic intact rocks [J]. Int J Rock Mech Min Sci, 2008, 45(3): 223-234.

[4] XU Wei-ya, ZHANG Gui-ke. Study on orthotropic equivalent strength parameters of jointed rock mass [J]. Chinese Journal of Geotechnical Engineering, 2007, (6): 806-810. (in Chinese)

[5] YANG Qiang, CHEN Xin, ZHOU Wei-yuan. Anisotropic yield criterion for jointed rock masses based on a two-order damage tensor [J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(8): 1275-1282. (in Chinese)

[6] YU Tian-tang. Modeling of inherent anisotropy for geotechnical material [J]. Chinese Journal of Rock Mechanics and Engineering, 2004, 23(10): 1604-1607. (in Chinese)

[7] COLAK K, UNLU T. Effect of transverse anisotropy on the Hoek-Brown strength parameter ‘m_i’ for intact rocks [J]. Int J Rock Mech Min Sci, 2004, 41(3): 1045-1052.

[8] LYDZBA D, PIETRUSZCZAK S, SHAO J F. On anisotropy of stratified rocks: Homogenization and fabric tensor approach [J]. Computers and Geotechnics, 2003, 30(3): 289-302.

[9] PIETRUSZCZAK S, LYDZBA D, SHAO J F. Modelling of inherent anisotropy in sedimentary rocks [J]. Int J Solids Struct, 2002, 39(3): 637-648.

[10] USHAKSARAEI R, PIETRUSZCZAK S. Failure criterion for structural masonry based on critical plane approach [J]. J Eng Mech, 2002, 128(7): 769-778.

[11] PIETRUSZCZAK S, MROZ Z. On failure criteria for anisotropic cohesive-frictional materials [J]. Int J Numer Anal Meth Geomech, 2001, 25(3): 509-524.

[12] SMITH B, CHEATHAM B. Anisotropic compacting yield condition applied to porous limestone [J]. Int J Rock Mech Min Sci, 1980, 17(3): 159-165.

[13] NOVA R. The failure of transversely isotropic rocks in triaxial compression [J]. Int J Rock Mech Min Sci, 1980, 17(6): 325-332.

[14] JAEGER J C. Shear failure of anisotropic rock [J]. Geol Mag, 1960, 97: 65-72.

[15] ZHANG Yu-jun, TANG Yi-xing. 2-D FEM analysis for an underground opening considering the strength-anisotropy of the layered rock mass [J]. Chinese Journal of Geotechnical Engineering, 1999, 21(3): 307-310. (in Chinese)

[16] ATTEWELL B, SANDFORD R. Intrinsic shear strength of a brittle, anisotropic rock (I): Experimental and mechanical interpretation [J]. Int J Rock Mech Min Sci, 1974, 11(11): 423-430.

[17] ZHOU Wei-yuan. Higher rock mechanics [M]. Beijing: Water Conservancy and Power Press, 1990: 116-119. (in Chinese)

[18] OWEN J, HINTON E. Finite elements in plasticity: Theory and practice [M]. West Cross, Swansea: Pineridge Press Limited, 1980: 215-231.

[19] BRADY G, BROWN T. Rock mechanics for underground mining [M]. Beijing: Coal Industry Press, 1990: 83-84. (in Chinese)

[20] YU Xue-fu, ZHENG Ying-ren, LIU Huai-heng, FANG Zheng-chang. Stability analysis for surrounding rock masses of underground projects [M]. Beijing: Coal Industry Press, 1983: 505-508. (in Chinese)

(Edited by CHEN Wei-ping)

Corresponding author: ZHANG Yu-jun, PhD; Tel: +86-27-87198482; E-mail: yjzhang@whrsm.ac.cn