- Email: [email protected]
An elasto-plastic damage constitutive model with double yield surfaces for saturated soft rock
ARTICLE IN PRESS International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
Contents lists available at ScienceDirect
International Journal of Rock Mechanics & Mining Sciences journal homepage: www.elsevier.com/locate/ijrmms
An elasto-plastic damage constitutive model with double yield surfaces for saturated soft rock C.Y. Zhou , F.X. Zhu School of Engineering, Geotechnical Engineering & Information Technology, Sun Yatsen University, No. 135 Xingang Xi Lu, Guangzhou 510275, China
a r t i c l e in fo
abstract
Article history: Received 4 May 2009 Received in revised form 2 December 2009 Accepted 7 January 2010 Available online 16 February 2010
An elasto-plastic damage constitutive model with double yield surfaces is developed based on irreversible thermodynamics theory and damage mechanics theory. Two kinds of plastic deformation mechanisms, including plastic friction and plastic pore deformation mechanisms, are considered. The plastic friction yield criterion is established by a parabolic open function that incorporates the volumetric deformation effect, and the motion of yield function in stress space is governed by its center position and rotation hardening rule. Meanwhile, a plastic pore yield criterion is established by adopting Gurson’s criterion using the proposed friction yield model to determine the matrix deformation of porous materials. A damage variable is defined to describe the development of various microscopic defects. Comparisons between numerical predictions and experimental data of triaxial compression tests are presented with various confining pressure conditions to verify the rationality of the proposed model. & 2010 Elsevier Ltd. All rights reserved.
Keywords: Plastic deformation mechanism Elasto-plastic damage constitutive model Double yield surfaces Saturated soft rock
1. Introduction Some serious hazards such as landslides, collapses and debris flows caused by soft rocks have been encountered frequently in the soft rock distributed areas due to the engineering construction activities. A series of experimental studies including mechanical behavior investigations and microstructure analysis of various kinds of soft rocks has been performed [1–12]. The results have shown that the mechanical behavior of the soft rocks after saturation is very complex and is characterized by various features such as elasto-plasticity, strain hardening, and performance deterioration. Several constitutive models such as the hyperbolic model, statistical damage model, elasto-plastic models with modified Mises or Drucker–Prager yield surface, modified cap model and strain softening model with unified strength law have been proposed by various researchers to show the mechanical behavior of this material [13–17]. Although all the models mentioned above can successfully predict the deformation features of some soft rocks in different ways, the following issues have not been fully resolved: (1) some of the models have been established by fitting the experimental data. The physical meaning of the model parameters is not clear and the theoretical background of the model is confusing; (2) most of them have been satisfied with the assumptions of the Generalized Hooke’s Law or the adaptation of classic yield criterion used in soil mechanics.
Corresponding author. Tel./fax: + 86 20 84036771.
E-mail address: [email protected] (C.Y. Zhou). 1365-1609/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijrmms.2010.01.002
These models provide a convenient conceptual tool for the evaluation of rock deformation. However, a general constitutive model of mechanical behavior for the soft rocks in their saturated state has not been proposed. Therefore, it is essential to develop a general constitutive model in order to have a good description of the complex features of saturated soft rocks. Two kinds of plastic deformation mechanisms can be identified based on the mineral compositions and microstructure features of soft rocks obtained from the experimental data: (1) due to the water reduced dilation of the inner pores, saturated soft rocks become high porosity materials, which behave like porous viscous material. Thus, the plastic deformation may be partly generated by inelastic pore collapse by breaking the contact forces between grains; (2) soft rocks consist of granular and powder ingredients, and these can be considered as inner friction material. This characteristic is not obvious in the initial stage of compressive deformation for the materials of saturated soft rocks. But with continuous pore collapse, the contact area between grains will be increased, which may cause the enhancement of inner friction. Consequently, another plastic deformation may be generated by mutual sliding between the framework grain interfaces. Therefore, it can be assumed that the plastic flow of soft rocks is a combination of plastic friction mechanism and plastic pore deformation mechanism. The plastic friction mechanism is usually described by the Mohr–Coulomb criterion and Drucker–Prager criterion is based on the assumption of neglecting plastic volumetric deformation and it is the typical plastic shearing mechanism for geological materials [18–21]. It is obvious that the mechanism of volumetric
ARTICLE IN PRESS 386
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
consolidation is related to the microstructure of material. For granular and powder materials, the volumetric consolidation may be caused by the rearrangement of grains after reducing of pore space. For the study of soft rocks in this paper, due to the nonhomogeneous distribution of clay minerals and framework grains, the stress distribution is not uniform in the materials under loading. So the plastic deformation related to grains rearrangement may contain two parts: volumetric deformation and shearing deformation, both induce energy dissipation. For the plastic pore deformation mechanism, this porous material is usually assumed to be composed of an equivalent solid matrix and a connected porosity. But the mechanism has been described by some classic shearing yield criteria [15,22,23], which are obviously unsuitable. Recently, some special yield surfaces [24–26] have been established for some geomaterials based on Gurson’s criterion, which was a general macroscopic yield criterion for porous materials with homogeneous solid matrix [27–29]. In this criterion, the plastic behavior of the solid matrix is described by the von-Mises perfect plastic criterion. The micromechanics feature under this criterion has been studied by Perrin and Leblond [30] and Leblond and Perrin [31], which can be obtained by using a rigorous homogenization technique and is explicitly related to the yield stress of the equivalent solid matrix and porosity. Thus, Gurson’s criterion can provide a beneficial inspiration for this study on porous saturated soft rocks. In addition, many authors have already noted that when rock materials are subjected to external loading, the elasto-plastic deformation is often accompanied by internal damage due to the development of various microscopic defects [23,32–34]. This determines another deformation feature for soft rocks. Furthermore, such damage will have influence on the elastic and plastic properties of the material; its evolution depends largely on the state of stress and strain [35–37]. Thus, in order to develop an elasto-plastic constitutive model to predict the inelastic behavior of these materials, a systematic framework for describing the coupled process of elasto-plastic deformation and material damage is required. One of the most promising approach for the systematic formulation of inelastic constitutive equations has been provided by the irreversible thermodynamics theory and this approach has been successfully applied to solving the problems of rock damage [33,34,38,39]. This paper considers the constitutive and damage evolution modelling of saturated soft rocks under the framework of irreversible thermodynamics theory and damage mechanics theory. Firstly, the functions of double yield surfaces are conducted to describe the plastic friction mechanism and the plastic pore deformation mechanism. After that, an elasto-plastic damage constitutive model for saturated soft rock is proposed based on the non-associated flow rule. Finally, the results obtained from the prediction of the model are compared with the laboratory results of the saturated soft rock samples.
2. Plastic friction mechanism of soft rocks In this section, volumetric deformation effect is accounted for the analysis of plastic deformation, and the plastic friction yield criterion is proposed for soft rock materials in the thermodynamics theoretical framework for constructing elastic/plastic constitutive models established by Collins. This is an improvement of classical shearing criteria. 2.1. Plastic friction yield criterion As mentioned above, for soft rock subjected to external loading, the strain increment contains two parts: volumetric
strain increment dev and shearing strain increment des. And these two strain tensors comprise elastic deformation and plastic deformation, respectively: dev ¼ deev þdepv ;
des ¼ dees þ deps
ð1Þ
When the deformation is isothermal, the incremental work done by the applied stresses is the sum of the increment in the Helmholtz free energy function, c, and the increment in the dissipation function, f. These functions are defined per unit volume. The equation can be expressed as follows [40,41]: p0 dev þ q des ¼ dc þ df
ð2Þ
df Z 0
This inequality is the statement of the second law of thermodynamics, which is appropriate for isothermal deformations, and it is a strict ‘greater than’ whenever any irreversible plastic deformations occur. The variables p0 and q denote the stress invariants: p0 ¼ sm rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 q ¼ seq ¼ ss 2 d d
ð3Þ
where seq is the effective stress, sd is the deviatoric stress tensor, and sm = skk/3 is hydrostatic stress. According to the results obtained by Collins, the increment in the free energy and dissipation function can be written as dc ¼ ð@c=@ev Þdev þ ð@c=@es Þdes þ ð@c=@epv Þdev þ ð@c=@eps Þdev x df ¼ ½@ðdfÞ=@ðdepv Þdepv þ ½@ðdfÞ=@ðdeps Þdeps ¼ p0 depv þ t deps
ð4Þ
where p0 ¼ @ðdfÞ=@ðdepv Þ is the dissipative effective pressure and t ¼ @ðdfÞ=@ðdeps Þ is the dissipative shear invariant. For decoupled materials, the free energy function c can be expressed as the sum of a function of only the elastic strains, plus a function of only the plastic strains:
c ¼ c1 ðeev ; ees Þ þ c2 ðeev ; ees Þ
ð5Þ
The total work increment associated with the effective stress can be written as the sum of the elastic and plastic components: dwe ¼ p0 deev þq dees ¼ dc1 dwp ¼ p0 depv þ q deps ¼ dc2 df ¼ ðp0 þ r0 Þdepv þ ðt þ zÞde
ð6Þ
where r0 ¼ @c2 =@epv ; x ¼ @c2 =@eps is the ‘‘shift stress’’ [40,42]. Hence, by using (6), the following relationship can be obtained: p0 ¼ p0 þ r0 ; q ¼ t þ x and (
sm r0 ¼ p0 seq z ¼ t
ð7Þ
A dissipation increment function including the combination dissipation of plastic volumetric and shearing deformation is [40] df ¼ ½Aðdepv Þ2 þ Bðdeps Þ2 1=2
ð8Þ
where A denotes the dissipation generated by the unit plastic volumetric strain and B denotes the dissipation generated by the unit plastic shearing strain. To simplify this model, Collins made the additional assumption that the second part of the free energy function, c2, is dependent on the volumetric part of the plastic strain, but not on its shear component. As a consequence, the shift stress only has an isotropic component, r0 , and the other shift stress, x, satisfies x ¼ @c2 =@eps ¼ 0. Hence the dissipative stresses are
p0 ¼ @ðdfÞ=@ðdepv Þ ¼ A depv =df;
t ¼ @ðdfÞ=@ðdeps Þ ¼ B deps =df
ð9Þ
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
387
τ CSL
B
O
NCL
A
π'
Fig. 1. Yield locus in dissipative stress space (from Ref. [40]).
Combining Eqs. (8) and (9), we can obtain the yield locus in dissipative stress space as follows (Fig. 1):
p0 2
t2
¼1 ð10Þ B By introducing Eq. (7) in (10), the yield locus in true stress space is obtained: A
þ
ðsm r0 Þ2 ðseq Þ2 þ ¼1 A B
ð11Þ
where the shift stress r0 leads to the integral translation of the yield surface, and the parameters A and B have the dimensions of stress and hence must be functions connected with three defining variables, {sm, seq, r0 } [43], i.e., A= f{sm, seq, r0 } and B =f2{sm, seq, r0 }. This yield function is capable of describing a wide range of situations by adjusting the form of A and B. For example, Collins and Kelly [40] and Collins and Houlsby [44] proved that the classic linear criterion, such as the Drucker–Prager yield criterion, can be recovered by choosing B ¼ M2 s2m , without considering plastic volumetric strain, where M denotes the slope of the critical state line. The problem in this section is how to describe the plastic deformation generated by inner friction effects of soft rocks, and the proposed model is similar to the classic shear friction model to some extent except for the consideration of the energy dissipation generated by plastic volumetric deformation. Two principles should be satisfied when choosing the function forms of A and B: (1) the yield surface model should accord with the function form of classic shearing yield surface, that is a parabolic function; (2) if the power exponent in the yield function is changed to ‘1’, this function will degrade to the strength envelop curve for rocks. On this basis, the functions A, B, and r0 are expressed as the following specific forms: A ¼ r0 sm r0 2 B ¼ ðM r0 Þ2 p r0 ¼ r 2
Fig. 2. The friction yield surface considering rotation hardening rule.
Assuming the component of Eq. (14) equals ‘1’, the yield function transforms to one type of Drucker–Prager yield criterion satisfying postulate (2), which is given by
seq þ
6sin j ðsm Cctg jÞ ¼ 0 3sin j
2.2. Hardening rule The yield criterion proposed above is based on the assumption by Collins, that the shift stress only has an isotropic component, r0 , ignoring the effect of shift stress, x, which reflects the rotation regulation of yield surface. Therefore, it is necessary to introduce a hardening function, Z, to describe the rotation behavior of soft rocks from initial plastic yield to asymptotic failure in some extent, as illustrated in Fig. 2. Hence the plastic friction yield criterion can be rewritten as follows: F1 ¼ s2eq þ 12kZM 2 pr ðsm pr Þ ¼ 0
where pr =C cot j denotes the intercept of the yield surface locus with the sm-axis, j is the internal friction angle, C is the cohesion, and M=6sin j/(3 sin j). Substituting r0 , A, and B into Eq.(11), the following equation is obtained: p 2 sm r ðseq Þ2 2 ð13Þ ¼1 p 2 þ pr Mpr 2 sm r 2 2 2 Simplifying Eq. (12), the yield surface has been transformed to a parabola open function from an ellipse close function. That is, Mpr 2sm 1=2 2 seq ¼ ð14Þ 2 pr
ð16Þ
where the parameter k is introduced to reflect the influence of initial loading condition on the hardening function. Based on the hardening model proposed for hard argillites [13,45], a kind of mudstone, the hardening function Z can be expressed in a specific form:
Z ¼ Z0 þ ðZm Z0 Þ ð12Þ
ð15Þ
eps Rþ eps
ð17Þ
pffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 ð p1 p2 Þ2 þ ð p2 p3 Þ2 þð p3 p1 Þ2 , 3
where eps ¼ epeq ¼ e e e e e e Z0 is the threshold value for initial plastic yield, Zm is the parameter to determine the final rotation angle of yield surface, which corresponds to failure surface, and ‘1’ is usually selected for Zm, and R denotes the rotation velocity of yield surface.
3. Plastic pore deformation mechanism of soft rocks According to the previous idealized representation, the saturated soft rock material is considered as a typical porous medium with equivalent homogeneous solid matrix and connected macroscopic porosity. A pore deformation mechanism is studied in this section and a yield function is established based on Gurson’s criterion.
ARTICLE IN PRESS 388
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
3.1. Plastic pore deformation yield criterion Gurson’s criterion has been widely applied to some fields due to its clear physical meaning and simple form. Nevertheless, some limitations still exist, as this model only considers the pore evolution in an infinite rigid-plastic matrix satisfying the Mises criteria. Some researchers have introduced other yield criteria into Gurson’s physical model and established a yield criterion with different matrix materials based on the assumption that the pore materials exhibit symmetric responses under compression and tension. For example, Jeong [25] and Joong and Pan [46] used the Mohr–Coulomb criterion to describe the plastic deformation of matrix for pressure-sensitive materials, Zhao and Tao [47] and Yang et al. [48] proposed a plastic damage yield surface for geomaterials with the assumption that the matrix satisfies the Drucker–Prager yield model. On this basis, in order to adapt the models for porous soft rock, we reconsider Gurson’s model from the following aspects. (1) The yield criterion of soft rock’s matrix: For the Mises matrix, the effective stress satisfies seq ¼ s . The classic Gurson’s model satisfying the Mises criteria is as follows: S2 3Sm 1f 2 ¼ 0 FðS; s ; f Þ ¼ eq þ 2f cosh ð18Þ 2 2s s The stress term s in Gurson’s model actually equals to the equivalent stress seq of the matrix. For soft rocks, it is assumed that the matrix material obeys the friction yield criteria proposed in the upper section. So, the equivalent stress in matrix can be expressed as Mpr 2sm 1=2 kZ0 2 seq ¼ ð19Þ 2 pr (2) The relationship between the matrix’s spherical stress and the macroscopic spherical stress: This relationship has been studied by some scholars and different results are obtained, but none of them was proved perfectly. Nevertheless, one of the conclusions is agreeable due to the following proving process. For a spherical voxel of porous media, which bears a macroscopic stress (Sm) in the location where the radius r equals b, Timoshenko and Goodier [49] give a linear-elastic analytical solution of radial stress (srr) and circumferential stress (syy,sjj):
srr ¼
b3 ðr 3 a3 Þ Sm ; r 3 ðb3 a3 Þ
syy ¼ sjj ¼
b3 ð2r 3 þa3 Þ Sm 2r 3 ðb3 a3 Þ
ð20Þ
Assuming that f= (a/b)3 denotes the porosity, Eq. (20) can be written as ðr 3 a3 Þ srr ¼ 3 Sm ; r ð1f Þ
ð2r 3 þa3 Þ syy ¼ sjj ¼ 3 Sm 2r ð1f Þ
ð21Þ
Therefore, 1 3
sm ¼ ðsrr þ syy þ sjj Þ ¼
1 Sm ð1f Þ
(2) Considering the single action of hydrostatic stress, the analytical solution of hydrostatic yield stress can be calculated by the yield function, and when the matrix is assumed to satisfy the Mises yield criteria, the result will be equal to ð2=3Þs log f , which is the hydrostatic yield stress solution of the classic Gurson model. (3) When the porosity f equals 0, the macroscopic yield surface will reduce to the yield function of the matrix. On this basis, the yield surface function for describing the plastic pore deformation of soft rocks is proposed as follows: F2 ðS; seq ; f Þ ¼ where seq ¼
S2eq
s2eq
Mpr 2
þ 2f cosh
3Sm 1f 2 ¼ 0 2seq i1=2
ð24Þ
h Z0 Sm 2kZ0 2k ð1f Þpr
Note that this yield surface contains three features: (1) It is a convex and smooth surface. (2) It has clear geometric meaning and definite physical connotation to choose the porosity as an inner variable. Two advantages are presented: reflecting the irreversible change of inner state and easily establishing the connection with macroscopic mechanical behavior. (3) Assuming the porous material presents symmetric responses under compression and tension and ignoring pore nucleation, when the pore volume percentage decreases, the yield surface will dilate gradually, making the plastic pore deformation weak. In addition, it is assumed that the changing of pore deformation yield surface is related to the plastic hardening of matrix material, according to the study of plastic deformation of porous rocks [24], the rotation hardening function Z0 is chosen to describe the hardening rule of the model:
Z0 ¼ Z0 þ ðZm Z0 Þ
eps R þ eps
ð25Þ
where eps is the plastic shearing strain of the matrix, and equals the plastic equivalent strain epeq , and other parameters are defined as above. Moreover, the relationship between the matrix’s plastic deformation and macroscopic plastic deformation of soft rocks should be established. Based on the definition of macroscopic deformation rate by Bishop and Hill, Gurson [28] established a dissipation reciprocal relationship using the principle of virtual work: Z Z p sij e_ pij dV ¼ sij e_ pij dV ð26Þ V Sij E_ ij ¼ V
Vmatrix
The equivalent stress (strain) can be defined as a linear stress (strain) state substituting for a complex stress (strain) state with the same effect. So, Eq. (26) can be rewritten as Z p V Sij E_ ij ¼ sij e_ pij dV ¼ Vmatrix seq e_ peq ð27Þ Vmatrix
ð22Þ
Hence, the expression of the equivalent stress of matrix can be obtained: 1=2 Mpr 2kZ0 Sm 2kZ0 seq ¼ ð23Þ 2 ð1f Þpr Moreover, according to the derivation of the yield function with Gurson’s form introduced in [25,46], several requirements should be satisfied: (1) The macroscopic yield surface will be reduced to Gurson’s yield surface, when the matrix obeys the Mises yield criteria.
or V Sij dEpij ¼ Vmatrix seq depeq
ð28Þ
Hence, depeq ¼
1 S dEp ð1f Þseq ij ij
ð29Þ
namely, deps ¼
1 S dEp ð1f Þseq ij ij
ð30Þ
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
3.2. Pore evolution model
where
During the deformation of soft rock, pore evolution is determined by the volume change of voids in the material. It is assumed that the pore evolution is only determined by plastic pore deformation mechanism, not taking account of the pore nucleation. Therefore,
s0eq ¼
V ¼ Vvoid þ Vmatrix
ð31Þ
where Vvoid is total pore volume, Vmatrix is the volume of matrix, and V is the macroscopic volume of soft rock. Using f= Vvoid/V, the following equation can be obtained: V dVvoid Vvoid dV dVdVmatrix Vvoid dV df ¼ d void ¼ ¼ V V V V V V V dV dVmatrix Vmatrix Vvoid dV ð32Þ ¼ Vmatrix V V V V
1=2 Mpr 2kZ0 sm 2kZ0 2 ð1f Þpr
Z ¼ Z0 þ ðZm Z0 Þ
eps Rþ eps
Z0 ¼ Z0 þ ðZm Z0 Þ de p2 ¼
e p2 R þ e p2
1 s dep2 ð1f Þs0eq ij ij
389
ð43Þ
ð44Þ
ð45Þ
ð46Þ
4. Elasto-plastic damage constitutive model of soft rocks 4.1. Constitutive theory of irreversible thermodynamics
Using 1 f= Vmatrix/V, Eq. (32) can be written as dV dV dV dV dVmatrix df ¼ ð1f Þ matrix f ¼ ð1f Þ Vmatrix Vmatrix V V V
As a thermodynamic potential function, the Helmholtz energy ð33Þ dEpv ,
where dV/V is the macroscopic volumetric strain increment and dVmatrix/Vmatrix denotes volumetric strain increment depv of the matrix. Hence, Eq. (33) can be expressed as df ¼ ð1f ÞðdEpv depv Þ
ð34Þ
As the yield criteria for the matrix material can be transformed to one of the Drucker–Prager yield criteria under some conditions: FDP ¼ seq þ
6sin j ðsm C tan jÞ ¼ 0 3sin j
ð35Þ
According to the rule of orthogonal flow, the relative plastic strain increment can be expressed as depv ¼ dl
@FDP ; @sm
deps ¼ dl
@FDP @seq
Hence, we obtain @FDP 6sin j depv ¼ =@seq deps ¼ dep @sm 3sin j s
ð37Þ
Inserting Eq. (38) into Eq. (34), the pore evolution equation is given by p 6sin j Sij dEij df ¼ ð1f ÞdEpv ð39Þ 3sin j seq Finally, in order to keep the variables of the above models consistent with each other, s is chosen to denote the macroscopic stress, ep2 to denote the macroscopic plastic strain caused by pore deformation, e p2 to denote the plastic shearing strain of matrix, and s0eq to denote the equivalent stress of matrix. So, the above models can be rewritten as
p2 6sin j sij deij df ¼ ð1f Þdep2 v 0 3sin j seq
c ¼ cðee ; T; vk Þ
ð47Þ
From which it follows that
c_ ¼
@c e @c _ @c :Tþ : e_ þ : v_ k @ee @T @vk
ð48Þ
In the case of infinitesimal deformation, the total strain can be expressed as the sum of the elastic and the plastic strain ð49Þ
By differentiating Eqs. (48) and (49) with respect to time and substituting the resulting expression into the Clausius–Duhem inequality, the following elastic constitutive equation as a result of the non-negative requirement of entropy production can be obtained: @c @c ð50Þ ; S¼ @ee @T The thermodynamic conjugate forces Ak to vk can be defined as follows:
s¼r
Ak ¼ r
@c ðk ¼ 1; 2; . . . ; nÞ @vk
ð51Þ
Furthermore, the mechanical flux vector J and their thermodynamic conjugate force vector X can be defined as follows: n T g oT ð52Þ J ¼ r e_ p ; v_ k ; q ; X ¼ s; Ak ; T The rate of entropy production can be expressed as the scalar product of X and J as follows [38]:
r F ¼ XJ Z0
ð53Þ
ð40Þ
If a component Jm of J is assumed to be a functional of the corresponding thermodynamic force only, the existence of a dissipation potential may be proved
ð41Þ
and Jm is obtained as follows:
F ¼ FðXm ; vk ; TÞ
!
s2eq 3sm F2 ¼ 0 2 þ 2f cosh 1f 2 ¼ 0 2s0eq s eq
internal variables vk(k= 1,2, y, n), where ee and T are conventional internal variables, and vk denotes the changes in density and configuration of dislocations, and the development of microscopic cavities and other imperfections. By using these internal variables, the Helmholtz free energy of the materials can be written as [39]
e ¼ ee þ ep ð36Þ
Substituting Eq. (30) into Eq. (37), the relationship between matrix’s plastic deformation and macroscopic plastic deformation of soft rock is established as Sij dEpij 6sin j depv ¼ ð38Þ 3sin j ð1f Þseq
F1 ¼ s2eq þ 12kZM 2 pr ðsm pr Þ ¼ 0
c is related to elastic strain ee, absolute temperature T and other
Jm ¼ ð42Þ
@F @Xm
ð54Þ
ð55Þ
where the quantities before the semicolon (;) denote variables while those after (;) are parameters.
ARTICLE IN PRESS 390
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
The specific forms of Eq. (55) can be written as @F @F @F e_ p ¼ ; v_ k ¼ ; q¼ @s @Ak @ðg=TÞ
ð56Þ
These equations imply that the evolution equations of internal state variables can be derived from the dissipation potential F. It can be seen that the first item of Eq. (56) represents the inelastic constitutive equations, which is formally similar to the flow rule of plastic mechanics. The above thermodynamic theory provides a systematic framework for establishing constitutive and damage evolution equations. 4.2. Damage variable and damage evolution model When saturated soft rock materials are subjected to external loading, the elasto-plastic deformation is often accompanied by internal damage due to the development of various microscopic defects such as micropores and microcracks. Therefore, it is necessary to define a unified damage variable to describe this performance deterioration. According to Lemaitre [36], a scalar damage variable D is chosen, just considering the isotropic damage condition, and the effective stress tensor can be given by
sij
s~ ij ¼
ð57Þ 1D The damage influences of elastic compliance and plastic yield condition are taken into account. Hence the effective elastic compliance matrix is expressed as 1
1 C~ ¼
ð1DÞ2
C 1 2
1 6 6 v0 6 6 v0 1 6 ¼ 6 2 0 ð1DÞ E0 6 6 6 0 4 0
3
v0
v0
0
0
0
1
v0
0
0
0
v0
1
0
0
0
0
0
2ð1 þ v0 Þ
0
0
0 0
0 0
0 0
2ð1þ v0 Þ 0
0 2ð1 þ v0 Þ
7 7 7 7 7 7 7 7 7 5
ð58Þ By introducing the damage variable into the plastic friction yield function, the effective yield criterion can be obtained: F~ 1 ¼ ðs~ eq Þ2 þ 12kZM 2 pr ðs~ m pr Þ ¼ 0
ð59Þ
where s~ m ¼ sm =ð1DÞ, s~ eq ¼ seq =ð1DÞ. For the plastic pore deformation yield condition, it is only considered the damage effect of matrix induced by the initiation and development of microcracks. So the effective yield criterion can be written as " # s2eq 3ðsm Þ 1f 2 ¼ 0 ð60aÞ F~ 2 ¼ 0 2 þ 2f cosh 2s0eq s eq p2 6sin j sij deij df ¼ ð1f Þdep2 v 3sin j s0eq 1=2 Mpr 2kZ0 s~ m 2kZ0 2 ð1f Þpr
e p2 Z0 ¼ Z0 þ ðZm Z0 Þ R þ e p2 de
p2
1 ¼ s dep2 ð1f Þs0eq ij ij
s~ m ¼
sm 1D
Vk ¼ fvp ; D; bg
ð60bÞ
ð60cÞ
Ak ¼ fRp ; Y; Bg
ð62Þ
Hence, the Helmholtz free energy function per unit volume of material may be written as
rCðee ; D; vp ; bÞ ¼ rCe ðeeij ; DÞ þ rCp ðvp Þ þ rCd ðbÞ e
ð63Þ
e ij ; DÞ
denotes the elastic strain energy incorporating where rC ðe the effects of the damage variable D on the elastic response of the material, whereas rCp ðvp Þ and rCd(b) are the free energy induced by isotropic hardening and damage development, respectively. The term rCd(b) can be expressed as a simple form with internal variable b [35]: d
pc ðbÞ ¼ 12Cd b
2
ð64Þ
Combining Eqs. (51), (63), and (64), the thermodynamic conjugate force corresponding to b is @c ¼ Cd b @b
ð65Þ
where Cd is a material constant. The other thermodynamic conjugate force Y corresponding to damage variable D can be written as a form of strain energy release rate [50,51], so under the isotropic damage condition, Y may be expressed as Y¼
ð60dÞ
ð61Þ
The current damage state of the material is given by the damage variable D and the internal variable b has been introduced to prescribe the effects of current damage state on the further development of damage, which is similar to the hardening variable of plastic deformation. The thermodynamic conjugate forces corresponding to the internal variables are expressed as
B¼r
where
s0eq ¼
The method adopted for describing damage evolution regulation of rocks mainly concentrates in two aspects: (1) From the angle of phenomenal mechanics and based on statistical distribution theory, the damage variable is assumed to obey the regulation of Weibull distribution or average distribution, and the damage evolution function can be conducted [33,34]. (2) According to the continuum mechanics theory [38], some internal variables are chosen to deal with the damage problem. The former can describe the damage problem in some extent, but the model parameters are usually determined by fitting the experimental data, and their physical meaning are not clear. The latter is based on the irreversible thermodynamic theory, and the energy dissipation generated by plastic deformation or crack development during the damage evolution process is taken into account. The damage is seen as an internal variable in the thermodynamic system and a dissipation potential function is introduced to describe the damage surface. With this method, a systematic damage constitutive framework can be established. Therefore, the latter method is adopted here to analyze the damage evolution law of soft rock. Considering an isothermal condition, it will now employ a scalar isotropic hardening variable vp for plastic deformation, a scalar damage variable D to describe the anisotropic damage states, and another scalar variable b prescribing the current state of damage:
ðs~ Þ2 E0 ð1DÞ
ð66Þ
where ð60eÞ
ð60fÞ
s~ ¼ s~ eq 23ð1þ v0 Þ þ3ð12v0 Þðs~ m =s~ eq Þ
1=2
ð67Þ
The damage evolution equation can be obtained with a damage potential function. These experimental studies have revealed that a limited area, where the damage does not develop
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
and the material presents the elastic behavior, exists in the stress(strain) space [52,53]. Similar to plastic yield surface, some researchers proposed a damage surface Fd = 0 in the thermodynamic conjugate force space and established the relative damage criteria and dissipation potential function. However, there are no systematic results available for the determination of this damage surface so far. In this paper, based on the assumption of the existence of a damage surface and initial damage threshold, one of the simplest convex surfaces for the damage dissipation potential may be given as follows: Fd ðY; B; DÞ ¼ Yd BY0 ¼ 0
ð68aÞ
pffiffiffiffiffiffi Yd ¼ YY
ð68bÞ
where Y0 is the value of Yd corresponding to the initial damage threshold. Substituting Eq. (68) into (56), the incremental evolution equations of D and b can be obtained: dD ¼ dld
@Fd @Y
ð69Þ
db ¼ dld
@Fd @ðBÞ
ð70Þ
where ld Z0 denotes the damage consistence parameter. The loading and unloading conditions may be expressed as follows: (1) If Fd o0, the damage does not begin. Namely the multiplier satisfies dld = 0. (2) If Fd = 0, the new damage occurs. Namely the multiplier dld 40 and it can be obtained by using the damage consistence condition : @Fd dY @Y ð71Þ dld ¼ @F dB d @B db or dld ¼
@Fd =@Y dY Cd ð@Fd =@BÞ
ð72Þ
Based on the irreversible thermodynamic theory, one assumption has been made that the plastic deformation and damage evolution of soft rock are mutually independent considering isothermal condition. So, the Helmholtz free energy of the materials can be written as
rC ¼ rCe ðeeij ; DÞ þ rCp ðvp Þ þ rCd ðbÞ ¼ 12eeij C~ ijkl eekl þ rCp ðvp Þ þ rCd ðbÞ
ð73Þ
where C~ ijkl denotes effective elastic tensor and g(vp) is the potential function related to plastic dissipation. Substituting Eq. (73) into (50), the following equation can be obtained:
sij ¼ r
@C ¼ C~ ijkl eekl ¼ C~ ijkl ðekl epkl Þ @eeij
ð74Þ
The incremental form is dsij ¼ dC~ ijkl ðekl epkl Þ þ C~ ijkl ðdekl depkl Þ or 1 1 deij ¼ C~ ijkl dskl C~ ijkl dC~ ijkl ðeij epij Þ þ depij
1 1 ¼ C~ ijkl dskl þ dC~ ijkl skl þ depij 1
1
¼ C~ ijkl dskl þ
@C~ ijkl dD skl þ depij @D
1
ð75Þ
ð76Þ 1
where C~ ijkl dskl denotes the elastic strain increment, ð@C~ ijkl =@DÞdD skl denotes the damage strain increment and depij denotes the plastic strain increment. That means the total elasto-plastic strain increment of soft rock consists of three parts: elastic strain increment, damage strain increment, and plastic strain increment. Based on the study in the previous section, the plastic part is divided into a plastic friction component dep1 ij and a plastic shear : component dep2 ij p2 depij ¼ dep1 ij þ deij
ð77Þ
A potential function with a similar form of other mudstones [13] is chosen to describe the friction deformation of soft rock for the reason of lacking experimental data about plastic yield characteristics: sm pr ð78Þ Q1 ¼ seq þ aðsm pr Þln I0 where the variable I0 o0 denotes the intersection point between the potential surface and the axis sm, and the parameter a is the slope of the dividing line between the compression region and the dilation region. By using the potential function, the stress space is divided into two zones, corresponding to plastic compressibility and dilation, respectively. The boundary between these two zones is defined by the condition qQ1/qsm =0. In this model, the dividing line function can be expressed as [24] fb ¼ seq aðsm pr Þ ¼ 0
ð79Þ
Combining these two boundary functions, the expression of I0 can be obtained. Hence, considering the damage effect and nonassociated flow rule, the effective plastic potential function and dividing line can be expressed as follows: s~ m pr ð80Þ Q~ 1 ¼ s~ eq þ aðs~ m pr Þln I0 fb ¼ s~ eq aðs~ m pr Þ ¼ 0
4.3. Elasto-plastic damage constitutive law of soft rocks
391
ð81Þ
where s~ m ¼ sm =ð1DÞ, s~ eq ¼ seq =ð1DÞ. For plastic pore deformation, the effective plastic potential function is given formally similar to the yield criterion: " # s2eq 3sm ¼0 ð82Þ Q~ 2 ¼ 0 2 þ 2f cosh 2s0eq s eq where the parameters are the same as above. According to the non-associated rule of orthogonal flow, the following two equations can be obtained: dep1 ij ¼
dl1 @Q~ 1 @sij
ð83Þ
dep2 ¼ ij
dl2 @Q~ 2 @sij
ð84Þ
The plastic multiplier dl can be determined by the consistency principle dF~ 1 ¼
@F~ 1 @F~ 1 @F~ 1 @F~ 1 @F~ 1 dD þ dD dsij þ dZ ¼ dsij þ @sij @D @Z @sij @D @F~ 1 @Z p @F~ 1 @F~ 1 @F~ 1 @Z @Q~ 1 þ dD þ de ¼ dsij þ dl1 ¼0 @Z @eps s @sij @D @Z @eps @seq ð85Þ
@F~ 2 @F~ 2 @F~ 2 @F~ 2 @F~ 2 dD þ df ¼ dF~ 2 ¼ dsij þ dZ0 þ ds @sij @D @Z0 @f @sij ij
ARTICLE IN PRESS 392
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
þ
¼
@F~ 2 @D
dD þ
" p2 # @F~ 2 @Z0 @F~ 2 6sin j sij deij p2 p2 ð1f Þd d e þ e v @Z0 @e p2 @f 3sin j s0eq
@F~ 2 @F~ 2 @F~ 2 @Z0 1 @Q~ dD þ dsij þ sij dl2 2 0 0 p2 @sij @D @Z @e ð1f Þseq @sij " # @F~ 2 @Q~ 2 6sin j sij @Q~ 2 ¼0 ð1f Þdl2 þ dl2 @f @seq 3sin j s0eq @sij
ð86Þ
So,
dl1 ¼
@F~ 1 ds þ @F~ 1 dD ij @sij @D
dl2 ¼
ð87Þ
H1 @F~ 2 ds þ @F~ 2 dD ij @sij @D
ð88Þ
H2
where H1 ¼
@F~ 1 @Z @Q~ 1 @Z @eps @seq
ð89aÞ
@F~ 2 @Z0 1 @Q~ @F~ sij 2 2 0 0 p2 @Z @e ð1f Þseq @sij @f " # ~ ~ s @Q 2 6sin j ij @Q 2 ð1f Þ @seq 3sin j s0eq @sij
H2 ¼
ð89bÞ
Substituting dl1 and dl2 into Eqs. (83) and (84), two types of plastic deformation increment are obtained. During a general loading history, the plastic deformation of soft rock is determined by two plastic deformation mechanisms together (see Fig. 3). The present stress satisfies the compression state. In Fig. 3, the plane seq sm is divided into three areas by failure surface Fmax, friction yield surface F1 and pore deformation yield surface F2, including the elastic area A0, the plastic area A1, which is only relative to pore deformation mechanism and the plastic area A2, which is relative to two deformation mechanism simultaneously. At the initial stage of loading, the stress point falls into A0. With the increase of seq, the stress point reaches A1, and it represents the plastic pore deformation. And with continuous pore collapse, the contact area between grains will be increased. It may cause the enhancement of inner friction. So the friction mechanism is activated and plays a role in the pore deformation mechanism simultaneously, and the stress point reaches the A2 area; then more ligament fractures of skeleton occur with the further increase of seq and the deformation state of soft rock is basically determined by the transfixion of a fracture plane and development of a shear zone. Now the plastic friction mechanism plays a more important role than pore deformation mechanism, and the failure surface is reached finally.
eq Fmax: Failure Surface
O
m
Fig. 3. The plastic deformation determined by two plastic deformation mechanisms.
Three distinct constitutive domains can be identified as follows: (1) If F~ 1 o0; F~ 2 o 0, the applied stress state is fully inside the elastic domain or leads to an elastic unloading. No plastic flow occurs and therefore dl1 =dl2 = 0, namely depij ¼ 0. (2) If F~ 2 ¼ 0; dF~ 2 ¼ 0, but F~ 1 o 0 or F~ 1 ¼ 0; dF~ 1 o0, then the plastic pore deformation mechanism is activated while the plastic friction yield surface is not reached. The plastic multiplier dl2 40 is determined by Eq. (88) and dl =0. Namely, the plastic strain increment can be expressed as . depij ¼ dep2 ij ~ (3) If F 1 ¼ 0; dF~ 1 ¼ 0 and F~ 2 ¼ 0; dF~ 2 ¼ 0, the two plastic deformation mechanisms are both activated. The plastic multipliers dl1 40 and dl1 40 can be determined by Eqs. (87) and (88) together. The plastic strain increment can be expressed p2 as depij ¼ dep1 ij þdeij .
5. Verification of the proposed model The proposed model involves several parameters, which can be determined step by step from direct interpretation of triaxial compression test results of soft rocks and other similar mudstones. The two elastic parameters, E0 and v0, are obtained from the initial linear part of the stress–strain curve obtained from a conventional triaxial compression test. The values of E0 may depend on the confining pressure. This phenomenon is considered in this study. The cohesion C and inner friction angle j are obtained by drawing the locus of peak stresses in the p–q plane. The hardening parameters, Z0,Zm, R are determined by the study of other similar mudstones [13], and the parameter k is influenced intensively by the initial confining pressure. After fitting the test data, k is found to satisfy the logistic growth function: k =23.03 21.53/ (1+(sm0 /0.76)4.84), where sm0 is the initial confining pressure. The value of parameter a can be obtained by identifying the stress point corresponding to depv ¼ 0 in the stress–strain curve. By defining the locus of stresses at the onset of elastic strain as the initial damage stress point, and substituting this value into Eqs. (66) and (67), the initial damage threshold can be obtained. The value of material constant Cd can be determined by the experimental data of other elasto-plastic materials [35]. To verify the capacity of the proposed model for describing the mechanical behavior of soft rock, numerical simulations of several groups of triaxial compression tests by using the developed model are presented. These tests are performed on saturated silty mudstone, a kind of soft rock from the southern coastal region of China. The following values of parameters are used in Table 1. Note that f0 is the porosity of soft rock bearing no loading, ignoring the porosity variation during the process of applying the confining pressure in triaxial tests. Figs. 4–7 present the simulations of triaxial compression tests at four levels of confining pressures: 0, 0.5, 1.0, and 3.0 MPa. For convenience, compression is assumed to be positive. From Figs. 4–7, it can be seen that good agreement between the numerical simulations and experimental data is obtained. The mechanical behaviors of soft rock, such as plastic deformation, material damage, and transition from plastic compressibility to dilation, are correctly predicted by the proposed model. It can also be seen that the mechanical behavior is sensitive to confining pressure. At low confining pressures, say 0, 0.5, and 1 MPa, a brittle failure is observed, with a small softening regime and the material failure occurs with a small value of the axial
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
393
Table 1 Values of parameters in the model. E0 (MPa) Confining pressure 0 MPa
0.5 MPa
1 MPa
3 MPa
610
1280
1740
1765
v0
C (MPa)
j (deg)
Z0
Zm
R
a
f0 (%)
Cd
0.35
1.61
30.16
0.017
1
0.00008
1.41
15.1
5
Fig. 4. Simulation of a triaxial compression test on saturated soft rock with 0 MPa confining pressure.
Fig. 5. Simulation of a triaxial compression test on saturated soft rock with 0.5 MPa confining pressure.
strain. This means with the increase of seq, the friction mechanism is activated instead of the pore deformation mechanism in a short time. However, at a slightly higher confining pressure, say 3 MPa, the axial strain corresponding to the failure state is much higher, which means that the plastic friction yield surface is reached at a higher axial stress and the pore deformation mechanism plays a more important role in the whole process.
Fig. 6. Simulation of a triaxial compression test on saturated soft rock with 1 MPa confining pressure.
Fig. 7. Simulation of a triaxial compression test on saturated soft rock with 3 MPa confining pressure.
A falling phenomenon is shown at the last stage of the stress– strain curve. This means the transfixion of the fracture plane and development of the shear zone in soft rock finally causes the failure of material, and this strain softening behavior generally implies the formation of strain localization mode. Suitable regularization methods like non-local formulation and secondgradient modelling should be employed. However, this feature is beyond the scope of study of the present work and an asymptotic
ARTICLE IN PRESS 394
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
at various confining pressure conditions. It has been shown that the proposed model correctly describes the main features of the mechanical behavior of saturated soft rock. In addition, the related results of other mudstone have been used in the determination of model parameters. The comparisons between the simulations and experimental data presented in this paper are only represent the first phase of validation of the model. Further advanced tests are still needed to adjust these parameters and check the general validity of the proposed model.
Acknowledgments
Fig. 8. Damage evolution curves at different confining pressures (0, 1, and 3 MPa).
perfect plastic regime is obtained to reflect the post-peak failure behavior of soft rock. The results of numerical simulations are appreciably greater than experimental data, mainly occurring after the yield locus. It is because the pore water pressure is not taken into account in the present model, inducing the effective stress of soft rock for numerical calculations greater than the actual situation. Nevertheless, the differentiation is still in the range of allowable error. In Fig. 8, the damage evolution curves at different confining pressures (0, 1, and 3 Mpa) are presented. The amount of elastic damage only happened at a small percentage of total structural damage of soft rock, and most parts of the damage took place during the plastic deformation phase, which mainly occurred after the activation of plastic friction mechanism. In this phase, the development of various microscopic defects such as micropores and microcracks can all be considered to be the performance of damage evolution. Moreover, by comparing the curves at different levels of confining pressures, it can also be found that the amount of damage at higher confining pressure is more than the one at lower confining pressure. That means the increase in confining pressure can strengthen the structure and reduce the amount of total damage of frozen soil.
6. Conclusion In this paper, an elasto-plastic damage constitutive model based on irreversible thermodynamics theory and damage mechanics theory is developed for saturated soft rock. Two plastic deformation mechanisms are taken into account: plastic friction and plastic pore deformation. By incorporating the volumetric deformation effect and considering the motion of yield function in stress space governed by its center position and rotation hardening rule, the plastic friction yield criterion is proposed. Meanwhile the plastic pore yield criterion is established as an adaptation of Gurson’s criterion using the proposed friction yield model to determine the matrix deformation of porous materials. Moreover, a damage variable describing the development of various microscopic defects is defined. Based on this variable, an elasto-plastic damage constitutive model for saturated soft rock is proposed for taking into account the two plastic mechanisms. The non-associated flow rule is adopted in the model, and the plastic strain increment-type expression is obtained. Numerical simulations have been compared with experimental data from triaxial compression tests that were performed on the saturated samples
The research is partly supported by the Natural Science Foundation of China (NSFC) through Grant no. 40672194, Guangdong Natural Science Foundation through Grant no. 06104932 and Doctorial Discipline Foundation through Grant no. 20060558060. The authors gratefully acknowledge the contribution of Engineer You-Rou Wang in completing the experiment. References [1] Mondo NH, Bjorlykke K, Jahren J. Experimental mechanical compaction of clay mineral aggregates—changes in physical properties of mudstones during burial. Mar. Petrol Geol. 2007;24:289–311. [2] Koncagul EC. Comparison of uniaxial compressive strength test and slake durability index values of shale samples from the breathtt formation, Kentucky. PhD thesis, University of Missouri, Rolla, 1998. [3] Botts ME. The effects of slaking on the engineering behavior of clay shale. PhD thesis, University of Colorado, Boulder, 1986. [4] Yamaguchi H, Yoshida K, Kuroshima I, Fukuda M. Slaking and shear properties of mudstone. In: Rock mechanics and power plants. Rotterdam: Balkema; 1988. p. 133–44. [5] Yoshinaka R, Tran TV, Osada M. Non-linear, stress and strain-dependent behavior of soft rocks under cyclic triaxial conditions. Int J Rock Mech Min Sci 1998;35:941–55. [6] Ye GL, Naito K, Sawada K. Experimental study on soft sedimentary rock under plane strain compression and creep tests. In: Proceedings of the third Asian rock mech symp, 2004. p. 865–70. [7] Singh TN, Verma AK. Prediction of creep characteristic of rock under varying environment. Environ Geol 2005:559–68. [8] Eseme E, Littke R, Krooss BM. Factors controlling the thermo-mechanical deformation of oil shales: implications for compaction of mudstones and exploitation. Mar Petrol Geol 2006;23:715–34. [9] Corkum AG, Martin CD. The mechanical behaviour of weak mudstone (Opalinus Clay) at low stresses. Int J Rock Mech Min Sci 2007;44:196–209. [10] Szymakowski J. Shear testing of soft rock masses. Geotech Test J 2007; 30:133–40. [11] Billiotte J, Yang D, Su K. Experimental study on gas permeability of mudstones. Phys Chem Earth 2008;33:231–6. [12] Kodaka T, Oka F, Kitahara H, Ohta H, Otani J. Observation of compaction bands under triaxial conditions for diatomaceous mudstone. In: Hyodo M, Murata H, Nakata Y, editors. Proc int symp geomech geotech particulate media, Ube Yamaguchi, Japan, 12–14 September 2006. Rotterdam: Balkema; 2006. p. 69–75. [13] Shao JF, Zhu QZ, Su K. Modeling of creep in rock materials in terms of material degradation. Comp Geotech 2003:549–55. [14] Li Hangzhou, Liao Hongjian, Sheng Qian. Study on statistical damage constitutive model of soft rock based on unified strength theory. Chin J Rock Mech Eng 2006;25:1331–6. [15] Dashnor H, Albert G, Francoise H, Christophe A. Saturated and unsaturated behaviour modelling of Meuse–Haute/Marne argillite. Int J Plast 2007;23: 733–66. [16] Adachi T, Oka F. An elasto-plastic constitutive model for soft rock with strain softening. Int J Numer Anal Meth Geonmech 1995;19:233–47. [17] Jun Yinglian, Tai Quanzhou, Yuan Hua, Guo Liangdai. Reliability analysis of railway tunnel construction process within soft and weak rock mass using the unified elasto-plastic strength criterion. Key Eng Mater 2007:2521–4. [18] Lade PV. Modelling the strengths of engineering materials in three dimensions. Mech Cohes-Frictional Mater 1997;2:339–56. [19] Desai CS. Mechanics of materials and interfaces: the disturbed state concept. Boca Raton: CRC Press; 2001. [20] Peric D, Ayari MA. On the analytical solutions for the three-invariant Camclay model. Int J Plast 2002;18:1061–82. [21] Schroeder C. Du coccolithe au re´servoir pe´trolier: approche phe´nome´nologique du comportement me´canique de la craie en vue de sa mode´lisation a´ diffe´rentes e´chelles. PhD thesis, Univ Lie´ge, France, 2003.
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
[22] Zhang H, Zhou L, Huang H. Strain localization analysis of anisotropic saturated and partially saturated porous media. Rock Soil Mech 2004;25: 675–80. [23] Dai Yonghao, Chen Weizhong, Wu Guojun, Zhou Xide. Study on elastoplastic damage model of unsaturated rock mass and its application. Chin J Rock Mech Eng 2008;27:728–35. [24] Xie SY, Shao JF. Elastoplastic deformation of a porous rock and water interaction. Int J Plast 2006:2195–225. [25] Jeong HY. A new yield function and a hydrostatic stress-controlled void nucleation model for porous solids with pressure-senstitive matrices. Int J Solids Struct 2002;39:1385–403. [26] Gurson AL. Porous rigid-plastic materials containing rigid inclusions-yield function: plastic potential and void nucleation. In: Taplin DMR, editor. Proc int conf fract. Oxford: Pergamon; 1977. p. 357–64. [27] Gurson AL. Continuum theory of ductile rupture by void nucleation and growth: part I—yield criteria and flow rules for porous ductile media. J Eng Mater Technol 1977;99:2–15. [28] Gurson AL. Plastic flow and fracture behavior of ductile materials incorporating void nucleation, growth and interaction. PhD thesis, Brown Univ, Providence, RI, 1975. [29] Tvergaard V. Material failure by void growth to coalescence. Adv Appl Mech 1990;27:83–147. [30] Perrin G, Leblond JB. Accelerated void growth in porous ductile solids containing two populations of cavities. Int J Plast 2000;16:91–120. [31] Leblond JB, Perrin G. Introduction a´ la mcanique de la rupture ductile des me´taux. Lecture notes. Paris: Ecole Polytechnique; 1996 [in French]. [32] Fu Zhiliang. Nonlinear creep damage model for soft rock under uniaxial compression. J Coal Sci Eng 2006;12:37–9. [33] Cao Wengui, Zhao Minghua, Liu Chengxue. Study on the model and its modifying method for rock softening and damage based on Weibull random distribution. Chin J Rock Mech Eng 2004;23:3226–31. [34] Xu Weiya, Wei Lide. Study on statistical damage constitutive model for rock. Chin J Rock Mech Eng 2002;21:787–91. [35] Hayakawa K, Murakami S, Liu Y. An irreversible thermodynamics theory for elastic–plastic-damage materials. Eur J Mech A Solids 1998;17:13–32. [36] Lemaitre J. A course of damage mechanics, 2nd ed. Berlin: Springer; 1992. [37] Murakami S. Anisotropic damage in metals. In: Boehler JP, editor. Failure criteria of structured media. Rotterdam: Balkema; 1993. p. 99–119.
395
[38] Krajcinovic D, Lemaitre J. Continumn damage mechanics—theory applications. Vienna: Springer; 1986. [39] Lemaitre J, Chaboche JL, Shrivastava B. Mechanics of solid materials. Cambridge: Cambridge University Press; 1990. [40] Collins IF, Kelly PA. A thermomechanical analysis of a family of soil models. Geotechnique 2002;52:507–18. [41] Collins IF, Muhunthan B. On the relationship between stress—dilatancy, anisotropy, and plastic dissipation for granular materials. Geotechnique 2003;53:611–8. [42] Collins IF, Hilder T. A theoretical framework for constructing elastic/plastic constitutive models of triaxial tests. Int J Numer Anal Meth Geomech 2002; 26:1313–47. [43] Shen Zhujiang. Summary on the failure criteria and yield functions. Chin J Geotech Eng 1995;17:1–8. [44] Collins IF, Houlsby GT. Application of thermomechanical principles to the modelling of geotechnical materials. Proc R Soc London 1997: 1975–2001. [45] Chiarelli AS. Experimental investigation and constitutive modeling of coupled elastoplastic damage in hard claystones. PhD thesis, Univ Lille, France, 2000. [46] Jeong HY, Pan J. A macroscopic constitutive law for porous solids with pressure-sensitive matrices and its implications to plastic flow localization. Int J Solids Struct 1995;32:3669–91. [47] Zhao Jisheng, Tao Xiaxin. A constitutive model of the porous media based on micro-mechanics damage. J Xian Min Inst 2002;17:73–8. [48] Yang qiang, Chen xin, Zhou weiyuan. Elasto-plastic damage model for geomaterials and strain localization analyses. Chin J Rock Mech Eng 2004; 23:3577–83. [49] Timoshenko SP, Goodier JN. Theory of elasticity. New York: McGraw-Hill; 1951. [50] George ZV, Deliktas B. A coupling anisotropic damage model for the inelastic response of composite materials. Comput Meth Appl Mech Eng 2000;183: 159–99. [51] Abu Al-Rub RK, Voyiadjis GZ. On the coupling of anisotropic damage and plasticity models for ductile materials. Int J Solids Struct 2003;40:2611–43. [52] Holcomb DJ, Costin LS. Detecting damage surface in brittle materials using acoustic emission. J Appl Mech 1986;53:536–44. [53] Holcomb DJ. Discrete memory in rock: a review. J Rheol 1994;28:725–58.
Contents lists available at ScienceDirect
International Journal of Rock Mechanics & Mining Sciences journal homepage: www.elsevier.com/locate/ijrmms
An elasto-plastic damage constitutive model with double yield surfaces for saturated soft rock C.Y. Zhou , F.X. Zhu School of Engineering, Geotechnical Engineering & Information Technology, Sun Yatsen University, No. 135 Xingang Xi Lu, Guangzhou 510275, China
a r t i c l e in fo
abstract
Article history: Received 4 May 2009 Received in revised form 2 December 2009 Accepted 7 January 2010 Available online 16 February 2010
An elasto-plastic damage constitutive model with double yield surfaces is developed based on irreversible thermodynamics theory and damage mechanics theory. Two kinds of plastic deformation mechanisms, including plastic friction and plastic pore deformation mechanisms, are considered. The plastic friction yield criterion is established by a parabolic open function that incorporates the volumetric deformation effect, and the motion of yield function in stress space is governed by its center position and rotation hardening rule. Meanwhile, a plastic pore yield criterion is established by adopting Gurson’s criterion using the proposed friction yield model to determine the matrix deformation of porous materials. A damage variable is defined to describe the development of various microscopic defects. Comparisons between numerical predictions and experimental data of triaxial compression tests are presented with various confining pressure conditions to verify the rationality of the proposed model. & 2010 Elsevier Ltd. All rights reserved.
Keywords: Plastic deformation mechanism Elasto-plastic damage constitutive model Double yield surfaces Saturated soft rock
1. Introduction Some serious hazards such as landslides, collapses and debris flows caused by soft rocks have been encountered frequently in the soft rock distributed areas due to the engineering construction activities. A series of experimental studies including mechanical behavior investigations and microstructure analysis of various kinds of soft rocks has been performed [1–12]. The results have shown that the mechanical behavior of the soft rocks after saturation is very complex and is characterized by various features such as elasto-plasticity, strain hardening, and performance deterioration. Several constitutive models such as the hyperbolic model, statistical damage model, elasto-plastic models with modified Mises or Drucker–Prager yield surface, modified cap model and strain softening model with unified strength law have been proposed by various researchers to show the mechanical behavior of this material [13–17]. Although all the models mentioned above can successfully predict the deformation features of some soft rocks in different ways, the following issues have not been fully resolved: (1) some of the models have been established by fitting the experimental data. The physical meaning of the model parameters is not clear and the theoretical background of the model is confusing; (2) most of them have been satisfied with the assumptions of the Generalized Hooke’s Law or the adaptation of classic yield criterion used in soil mechanics.
Corresponding author. Tel./fax: + 86 20 84036771.
E-mail address: [email protected] (C.Y. Zhou). 1365-1609/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijrmms.2010.01.002
These models provide a convenient conceptual tool for the evaluation of rock deformation. However, a general constitutive model of mechanical behavior for the soft rocks in their saturated state has not been proposed. Therefore, it is essential to develop a general constitutive model in order to have a good description of the complex features of saturated soft rocks. Two kinds of plastic deformation mechanisms can be identified based on the mineral compositions and microstructure features of soft rocks obtained from the experimental data: (1) due to the water reduced dilation of the inner pores, saturated soft rocks become high porosity materials, which behave like porous viscous material. Thus, the plastic deformation may be partly generated by inelastic pore collapse by breaking the contact forces between grains; (2) soft rocks consist of granular and powder ingredients, and these can be considered as inner friction material. This characteristic is not obvious in the initial stage of compressive deformation for the materials of saturated soft rocks. But with continuous pore collapse, the contact area between grains will be increased, which may cause the enhancement of inner friction. Consequently, another plastic deformation may be generated by mutual sliding between the framework grain interfaces. Therefore, it can be assumed that the plastic flow of soft rocks is a combination of plastic friction mechanism and plastic pore deformation mechanism. The plastic friction mechanism is usually described by the Mohr–Coulomb criterion and Drucker–Prager criterion is based on the assumption of neglecting plastic volumetric deformation and it is the typical plastic shearing mechanism for geological materials [18–21]. It is obvious that the mechanism of volumetric
ARTICLE IN PRESS 386
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
consolidation is related to the microstructure of material. For granular and powder materials, the volumetric consolidation may be caused by the rearrangement of grains after reducing of pore space. For the study of soft rocks in this paper, due to the nonhomogeneous distribution of clay minerals and framework grains, the stress distribution is not uniform in the materials under loading. So the plastic deformation related to grains rearrangement may contain two parts: volumetric deformation and shearing deformation, both induce energy dissipation. For the plastic pore deformation mechanism, this porous material is usually assumed to be composed of an equivalent solid matrix and a connected porosity. But the mechanism has been described by some classic shearing yield criteria [15,22,23], which are obviously unsuitable. Recently, some special yield surfaces [24–26] have been established for some geomaterials based on Gurson’s criterion, which was a general macroscopic yield criterion for porous materials with homogeneous solid matrix [27–29]. In this criterion, the plastic behavior of the solid matrix is described by the von-Mises perfect plastic criterion. The micromechanics feature under this criterion has been studied by Perrin and Leblond [30] and Leblond and Perrin [31], which can be obtained by using a rigorous homogenization technique and is explicitly related to the yield stress of the equivalent solid matrix and porosity. Thus, Gurson’s criterion can provide a beneficial inspiration for this study on porous saturated soft rocks. In addition, many authors have already noted that when rock materials are subjected to external loading, the elasto-plastic deformation is often accompanied by internal damage due to the development of various microscopic defects [23,32–34]. This determines another deformation feature for soft rocks. Furthermore, such damage will have influence on the elastic and plastic properties of the material; its evolution depends largely on the state of stress and strain [35–37]. Thus, in order to develop an elasto-plastic constitutive model to predict the inelastic behavior of these materials, a systematic framework for describing the coupled process of elasto-plastic deformation and material damage is required. One of the most promising approach for the systematic formulation of inelastic constitutive equations has been provided by the irreversible thermodynamics theory and this approach has been successfully applied to solving the problems of rock damage [33,34,38,39]. This paper considers the constitutive and damage evolution modelling of saturated soft rocks under the framework of irreversible thermodynamics theory and damage mechanics theory. Firstly, the functions of double yield surfaces are conducted to describe the plastic friction mechanism and the plastic pore deformation mechanism. After that, an elasto-plastic damage constitutive model for saturated soft rock is proposed based on the non-associated flow rule. Finally, the results obtained from the prediction of the model are compared with the laboratory results of the saturated soft rock samples.
2. Plastic friction mechanism of soft rocks In this section, volumetric deformation effect is accounted for the analysis of plastic deformation, and the plastic friction yield criterion is proposed for soft rock materials in the thermodynamics theoretical framework for constructing elastic/plastic constitutive models established by Collins. This is an improvement of classical shearing criteria. 2.1. Plastic friction yield criterion As mentioned above, for soft rock subjected to external loading, the strain increment contains two parts: volumetric
strain increment dev and shearing strain increment des. And these two strain tensors comprise elastic deformation and plastic deformation, respectively: dev ¼ deev þdepv ;
des ¼ dees þ deps
ð1Þ
When the deformation is isothermal, the incremental work done by the applied stresses is the sum of the increment in the Helmholtz free energy function, c, and the increment in the dissipation function, f. These functions are defined per unit volume. The equation can be expressed as follows [40,41]: p0 dev þ q des ¼ dc þ df
ð2Þ
df Z 0
This inequality is the statement of the second law of thermodynamics, which is appropriate for isothermal deformations, and it is a strict ‘greater than’ whenever any irreversible plastic deformations occur. The variables p0 and q denote the stress invariants: p0 ¼ sm rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 q ¼ seq ¼ ss 2 d d
ð3Þ
where seq is the effective stress, sd is the deviatoric stress tensor, and sm = skk/3 is hydrostatic stress. According to the results obtained by Collins, the increment in the free energy and dissipation function can be written as dc ¼ ð@c=@ev Þdev þ ð@c=@es Þdes þ ð@c=@epv Þdev þ ð@c=@eps Þdev x df ¼ ½@ðdfÞ=@ðdepv Þdepv þ ½@ðdfÞ=@ðdeps Þdeps ¼ p0 depv þ t deps
ð4Þ
where p0 ¼ @ðdfÞ=@ðdepv Þ is the dissipative effective pressure and t ¼ @ðdfÞ=@ðdeps Þ is the dissipative shear invariant. For decoupled materials, the free energy function c can be expressed as the sum of a function of only the elastic strains, plus a function of only the plastic strains:
c ¼ c1 ðeev ; ees Þ þ c2 ðeev ; ees Þ
ð5Þ
The total work increment associated with the effective stress can be written as the sum of the elastic and plastic components: dwe ¼ p0 deev þq dees ¼ dc1 dwp ¼ p0 depv þ q deps ¼ dc2 df ¼ ðp0 þ r0 Þdepv þ ðt þ zÞde
ð6Þ
where r0 ¼ @c2 =@epv ; x ¼ @c2 =@eps is the ‘‘shift stress’’ [40,42]. Hence, by using (6), the following relationship can be obtained: p0 ¼ p0 þ r0 ; q ¼ t þ x and (
sm r0 ¼ p0 seq z ¼ t
ð7Þ
A dissipation increment function including the combination dissipation of plastic volumetric and shearing deformation is [40] df ¼ ½Aðdepv Þ2 þ Bðdeps Þ2 1=2
ð8Þ
where A denotes the dissipation generated by the unit plastic volumetric strain and B denotes the dissipation generated by the unit plastic shearing strain. To simplify this model, Collins made the additional assumption that the second part of the free energy function, c2, is dependent on the volumetric part of the plastic strain, but not on its shear component. As a consequence, the shift stress only has an isotropic component, r0 , and the other shift stress, x, satisfies x ¼ @c2 =@eps ¼ 0. Hence the dissipative stresses are
p0 ¼ @ðdfÞ=@ðdepv Þ ¼ A depv =df;
t ¼ @ðdfÞ=@ðdeps Þ ¼ B deps =df
ð9Þ
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
387
τ CSL
B
O
NCL
A
π'
Fig. 1. Yield locus in dissipative stress space (from Ref. [40]).
Combining Eqs. (8) and (9), we can obtain the yield locus in dissipative stress space as follows (Fig. 1):
p0 2
t2
¼1 ð10Þ B By introducing Eq. (7) in (10), the yield locus in true stress space is obtained: A
þ
ðsm r0 Þ2 ðseq Þ2 þ ¼1 A B
ð11Þ
where the shift stress r0 leads to the integral translation of the yield surface, and the parameters A and B have the dimensions of stress and hence must be functions connected with three defining variables, {sm, seq, r0 } [43], i.e., A= f{sm, seq, r0 } and B =f2{sm, seq, r0 }. This yield function is capable of describing a wide range of situations by adjusting the form of A and B. For example, Collins and Kelly [40] and Collins and Houlsby [44] proved that the classic linear criterion, such as the Drucker–Prager yield criterion, can be recovered by choosing B ¼ M2 s2m , without considering plastic volumetric strain, where M denotes the slope of the critical state line. The problem in this section is how to describe the plastic deformation generated by inner friction effects of soft rocks, and the proposed model is similar to the classic shear friction model to some extent except for the consideration of the energy dissipation generated by plastic volumetric deformation. Two principles should be satisfied when choosing the function forms of A and B: (1) the yield surface model should accord with the function form of classic shearing yield surface, that is a parabolic function; (2) if the power exponent in the yield function is changed to ‘1’, this function will degrade to the strength envelop curve for rocks. On this basis, the functions A, B, and r0 are expressed as the following specific forms: A ¼ r0 sm r0 2 B ¼ ðM r0 Þ2 p r0 ¼ r 2
Fig. 2. The friction yield surface considering rotation hardening rule.
Assuming the component of Eq. (14) equals ‘1’, the yield function transforms to one type of Drucker–Prager yield criterion satisfying postulate (2), which is given by
seq þ
6sin j ðsm Cctg jÞ ¼ 0 3sin j
2.2. Hardening rule The yield criterion proposed above is based on the assumption by Collins, that the shift stress only has an isotropic component, r0 , ignoring the effect of shift stress, x, which reflects the rotation regulation of yield surface. Therefore, it is necessary to introduce a hardening function, Z, to describe the rotation behavior of soft rocks from initial plastic yield to asymptotic failure in some extent, as illustrated in Fig. 2. Hence the plastic friction yield criterion can be rewritten as follows: F1 ¼ s2eq þ 12kZM 2 pr ðsm pr Þ ¼ 0
where pr =C cot j denotes the intercept of the yield surface locus with the sm-axis, j is the internal friction angle, C is the cohesion, and M=6sin j/(3 sin j). Substituting r0 , A, and B into Eq.(11), the following equation is obtained: p 2 sm r ðseq Þ2 2 ð13Þ ¼1 p 2 þ pr Mpr 2 sm r 2 2 2 Simplifying Eq. (12), the yield surface has been transformed to a parabola open function from an ellipse close function. That is, Mpr 2sm 1=2 2 seq ¼ ð14Þ 2 pr
ð16Þ
where the parameter k is introduced to reflect the influence of initial loading condition on the hardening function. Based on the hardening model proposed for hard argillites [13,45], a kind of mudstone, the hardening function Z can be expressed in a specific form:
Z ¼ Z0 þ ðZm Z0 Þ ð12Þ
ð15Þ
eps Rþ eps
ð17Þ
pffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 ð p1 p2 Þ2 þ ð p2 p3 Þ2 þð p3 p1 Þ2 , 3
where eps ¼ epeq ¼ e e e e e e Z0 is the threshold value for initial plastic yield, Zm is the parameter to determine the final rotation angle of yield surface, which corresponds to failure surface, and ‘1’ is usually selected for Zm, and R denotes the rotation velocity of yield surface.
3. Plastic pore deformation mechanism of soft rocks According to the previous idealized representation, the saturated soft rock material is considered as a typical porous medium with equivalent homogeneous solid matrix and connected macroscopic porosity. A pore deformation mechanism is studied in this section and a yield function is established based on Gurson’s criterion.
ARTICLE IN PRESS 388
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
3.1. Plastic pore deformation yield criterion Gurson’s criterion has been widely applied to some fields due to its clear physical meaning and simple form. Nevertheless, some limitations still exist, as this model only considers the pore evolution in an infinite rigid-plastic matrix satisfying the Mises criteria. Some researchers have introduced other yield criteria into Gurson’s physical model and established a yield criterion with different matrix materials based on the assumption that the pore materials exhibit symmetric responses under compression and tension. For example, Jeong [25] and Joong and Pan [46] used the Mohr–Coulomb criterion to describe the plastic deformation of matrix for pressure-sensitive materials, Zhao and Tao [47] and Yang et al. [48] proposed a plastic damage yield surface for geomaterials with the assumption that the matrix satisfies the Drucker–Prager yield model. On this basis, in order to adapt the models for porous soft rock, we reconsider Gurson’s model from the following aspects. (1) The yield criterion of soft rock’s matrix: For the Mises matrix, the effective stress satisfies seq ¼ s . The classic Gurson’s model satisfying the Mises criteria is as follows: S2 3Sm 1f 2 ¼ 0 FðS; s ; f Þ ¼ eq þ 2f cosh ð18Þ 2 2s s The stress term s in Gurson’s model actually equals to the equivalent stress seq of the matrix. For soft rocks, it is assumed that the matrix material obeys the friction yield criteria proposed in the upper section. So, the equivalent stress in matrix can be expressed as Mpr 2sm 1=2 kZ0 2 seq ¼ ð19Þ 2 pr (2) The relationship between the matrix’s spherical stress and the macroscopic spherical stress: This relationship has been studied by some scholars and different results are obtained, but none of them was proved perfectly. Nevertheless, one of the conclusions is agreeable due to the following proving process. For a spherical voxel of porous media, which bears a macroscopic stress (Sm) in the location where the radius r equals b, Timoshenko and Goodier [49] give a linear-elastic analytical solution of radial stress (srr) and circumferential stress (syy,sjj):
srr ¼
b3 ðr 3 a3 Þ Sm ; r 3 ðb3 a3 Þ
syy ¼ sjj ¼
b3 ð2r 3 þa3 Þ Sm 2r 3 ðb3 a3 Þ
ð20Þ
Assuming that f= (a/b)3 denotes the porosity, Eq. (20) can be written as ðr 3 a3 Þ srr ¼ 3 Sm ; r ð1f Þ
ð2r 3 þa3 Þ syy ¼ sjj ¼ 3 Sm 2r ð1f Þ
ð21Þ
Therefore, 1 3
sm ¼ ðsrr þ syy þ sjj Þ ¼
1 Sm ð1f Þ
(2) Considering the single action of hydrostatic stress, the analytical solution of hydrostatic yield stress can be calculated by the yield function, and when the matrix is assumed to satisfy the Mises yield criteria, the result will be equal to ð2=3Þs log f , which is the hydrostatic yield stress solution of the classic Gurson model. (3) When the porosity f equals 0, the macroscopic yield surface will reduce to the yield function of the matrix. On this basis, the yield surface function for describing the plastic pore deformation of soft rocks is proposed as follows: F2 ðS; seq ; f Þ ¼ where seq ¼
S2eq
s2eq
Mpr 2
þ 2f cosh
3Sm 1f 2 ¼ 0 2seq i1=2
ð24Þ
h Z0 Sm 2kZ0 2k ð1f Þpr
Note that this yield surface contains three features: (1) It is a convex and smooth surface. (2) It has clear geometric meaning and definite physical connotation to choose the porosity as an inner variable. Two advantages are presented: reflecting the irreversible change of inner state and easily establishing the connection with macroscopic mechanical behavior. (3) Assuming the porous material presents symmetric responses under compression and tension and ignoring pore nucleation, when the pore volume percentage decreases, the yield surface will dilate gradually, making the plastic pore deformation weak. In addition, it is assumed that the changing of pore deformation yield surface is related to the plastic hardening of matrix material, according to the study of plastic deformation of porous rocks [24], the rotation hardening function Z0 is chosen to describe the hardening rule of the model:
Z0 ¼ Z0 þ ðZm Z0 Þ
eps R þ eps
ð25Þ
where eps is the plastic shearing strain of the matrix, and equals the plastic equivalent strain epeq , and other parameters are defined as above. Moreover, the relationship between the matrix’s plastic deformation and macroscopic plastic deformation of soft rocks should be established. Based on the definition of macroscopic deformation rate by Bishop and Hill, Gurson [28] established a dissipation reciprocal relationship using the principle of virtual work: Z Z p sij e_ pij dV ¼ sij e_ pij dV ð26Þ V Sij E_ ij ¼ V
Vmatrix
The equivalent stress (strain) can be defined as a linear stress (strain) state substituting for a complex stress (strain) state with the same effect. So, Eq. (26) can be rewritten as Z p V Sij E_ ij ¼ sij e_ pij dV ¼ Vmatrix seq e_ peq ð27Þ Vmatrix
ð22Þ
Hence, the expression of the equivalent stress of matrix can be obtained: 1=2 Mpr 2kZ0 Sm 2kZ0 seq ¼ ð23Þ 2 ð1f Þpr Moreover, according to the derivation of the yield function with Gurson’s form introduced in [25,46], several requirements should be satisfied: (1) The macroscopic yield surface will be reduced to Gurson’s yield surface, when the matrix obeys the Mises yield criteria.
or V Sij dEpij ¼ Vmatrix seq depeq
ð28Þ
Hence, depeq ¼
1 S dEp ð1f Þseq ij ij
ð29Þ
namely, deps ¼
1 S dEp ð1f Þseq ij ij
ð30Þ
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
3.2. Pore evolution model
where
During the deformation of soft rock, pore evolution is determined by the volume change of voids in the material. It is assumed that the pore evolution is only determined by plastic pore deformation mechanism, not taking account of the pore nucleation. Therefore,
s0eq ¼
V ¼ Vvoid þ Vmatrix
ð31Þ
where Vvoid is total pore volume, Vmatrix is the volume of matrix, and V is the macroscopic volume of soft rock. Using f= Vvoid/V, the following equation can be obtained: V dVvoid Vvoid dV dVdVmatrix Vvoid dV df ¼ d void ¼ ¼ V V V V V V V dV dVmatrix Vmatrix Vvoid dV ð32Þ ¼ Vmatrix V V V V
1=2 Mpr 2kZ0 sm 2kZ0 2 ð1f Þpr
Z ¼ Z0 þ ðZm Z0 Þ
eps Rþ eps
Z0 ¼ Z0 þ ðZm Z0 Þ de p2 ¼
e p2 R þ e p2
1 s dep2 ð1f Þs0eq ij ij
389
ð43Þ
ð44Þ
ð45Þ
ð46Þ
4. Elasto-plastic damage constitutive model of soft rocks 4.1. Constitutive theory of irreversible thermodynamics
Using 1 f= Vmatrix/V, Eq. (32) can be written as dV dV dV dV dVmatrix df ¼ ð1f Þ matrix f ¼ ð1f Þ Vmatrix Vmatrix V V V
As a thermodynamic potential function, the Helmholtz energy ð33Þ dEpv ,
where dV/V is the macroscopic volumetric strain increment and dVmatrix/Vmatrix denotes volumetric strain increment depv of the matrix. Hence, Eq. (33) can be expressed as df ¼ ð1f ÞðdEpv depv Þ
ð34Þ
As the yield criteria for the matrix material can be transformed to one of the Drucker–Prager yield criteria under some conditions: FDP ¼ seq þ
6sin j ðsm C tan jÞ ¼ 0 3sin j
ð35Þ
According to the rule of orthogonal flow, the relative plastic strain increment can be expressed as depv ¼ dl
@FDP ; @sm
deps ¼ dl
@FDP @seq
Hence, we obtain @FDP 6sin j depv ¼ =@seq deps ¼ dep @sm 3sin j s
ð37Þ
Inserting Eq. (38) into Eq. (34), the pore evolution equation is given by p 6sin j Sij dEij df ¼ ð1f ÞdEpv ð39Þ 3sin j seq Finally, in order to keep the variables of the above models consistent with each other, s is chosen to denote the macroscopic stress, ep2 to denote the macroscopic plastic strain caused by pore deformation, e p2 to denote the plastic shearing strain of matrix, and s0eq to denote the equivalent stress of matrix. So, the above models can be rewritten as
p2 6sin j sij deij df ¼ ð1f Þdep2 v 0 3sin j seq
c ¼ cðee ; T; vk Þ
ð47Þ
From which it follows that
c_ ¼
@c e @c _ @c :Tþ : e_ þ : v_ k @ee @T @vk
ð48Þ
In the case of infinitesimal deformation, the total strain can be expressed as the sum of the elastic and the plastic strain ð49Þ
By differentiating Eqs. (48) and (49) with respect to time and substituting the resulting expression into the Clausius–Duhem inequality, the following elastic constitutive equation as a result of the non-negative requirement of entropy production can be obtained: @c @c ð50Þ ; S¼ @ee @T The thermodynamic conjugate forces Ak to vk can be defined as follows:
s¼r
Ak ¼ r
@c ðk ¼ 1; 2; . . . ; nÞ @vk
ð51Þ
Furthermore, the mechanical flux vector J and their thermodynamic conjugate force vector X can be defined as follows: n T g oT ð52Þ J ¼ r e_ p ; v_ k ; q ; X ¼ s; Ak ; T The rate of entropy production can be expressed as the scalar product of X and J as follows [38]:
r F ¼ XJ Z0
ð53Þ
ð40Þ
If a component Jm of J is assumed to be a functional of the corresponding thermodynamic force only, the existence of a dissipation potential may be proved
ð41Þ
and Jm is obtained as follows:
F ¼ FðXm ; vk ; TÞ
!
s2eq 3sm F2 ¼ 0 2 þ 2f cosh 1f 2 ¼ 0 2s0eq s eq
internal variables vk(k= 1,2, y, n), where ee and T are conventional internal variables, and vk denotes the changes in density and configuration of dislocations, and the development of microscopic cavities and other imperfections. By using these internal variables, the Helmholtz free energy of the materials can be written as [39]
e ¼ ee þ ep ð36Þ
Substituting Eq. (30) into Eq. (37), the relationship between matrix’s plastic deformation and macroscopic plastic deformation of soft rock is established as Sij dEpij 6sin j depv ¼ ð38Þ 3sin j ð1f Þseq
F1 ¼ s2eq þ 12kZM 2 pr ðsm pr Þ ¼ 0
c is related to elastic strain ee, absolute temperature T and other
Jm ¼ ð42Þ
@F @Xm
ð54Þ
ð55Þ
where the quantities before the semicolon (;) denote variables while those after (;) are parameters.
ARTICLE IN PRESS 390
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
The specific forms of Eq. (55) can be written as @F @F @F e_ p ¼ ; v_ k ¼ ; q¼ @s @Ak @ðg=TÞ
ð56Þ
These equations imply that the evolution equations of internal state variables can be derived from the dissipation potential F. It can be seen that the first item of Eq. (56) represents the inelastic constitutive equations, which is formally similar to the flow rule of plastic mechanics. The above thermodynamic theory provides a systematic framework for establishing constitutive and damage evolution equations. 4.2. Damage variable and damage evolution model When saturated soft rock materials are subjected to external loading, the elasto-plastic deformation is often accompanied by internal damage due to the development of various microscopic defects such as micropores and microcracks. Therefore, it is necessary to define a unified damage variable to describe this performance deterioration. According to Lemaitre [36], a scalar damage variable D is chosen, just considering the isotropic damage condition, and the effective stress tensor can be given by
sij
s~ ij ¼
ð57Þ 1D The damage influences of elastic compliance and plastic yield condition are taken into account. Hence the effective elastic compliance matrix is expressed as 1
1 C~ ¼
ð1DÞ2
C 1 2
1 6 6 v0 6 6 v0 1 6 ¼ 6 2 0 ð1DÞ E0 6 6 6 0 4 0
3
v0
v0
0
0
0
1
v0
0
0
0
v0
1
0
0
0
0
0
2ð1 þ v0 Þ
0
0
0 0
0 0
0 0
2ð1þ v0 Þ 0
0 2ð1 þ v0 Þ
7 7 7 7 7 7 7 7 7 5
ð58Þ By introducing the damage variable into the plastic friction yield function, the effective yield criterion can be obtained: F~ 1 ¼ ðs~ eq Þ2 þ 12kZM 2 pr ðs~ m pr Þ ¼ 0
ð59Þ
where s~ m ¼ sm =ð1DÞ, s~ eq ¼ seq =ð1DÞ. For the plastic pore deformation yield condition, it is only considered the damage effect of matrix induced by the initiation and development of microcracks. So the effective yield criterion can be written as " # s2eq 3ðsm Þ 1f 2 ¼ 0 ð60aÞ F~ 2 ¼ 0 2 þ 2f cosh 2s0eq s eq p2 6sin j sij deij df ¼ ð1f Þdep2 v 3sin j s0eq 1=2 Mpr 2kZ0 s~ m 2kZ0 2 ð1f Þpr
e p2 Z0 ¼ Z0 þ ðZm Z0 Þ R þ e p2 de
p2
1 ¼ s dep2 ð1f Þs0eq ij ij
s~ m ¼
sm 1D
Vk ¼ fvp ; D; bg
ð60bÞ
ð60cÞ
Ak ¼ fRp ; Y; Bg
ð62Þ
Hence, the Helmholtz free energy function per unit volume of material may be written as
rCðee ; D; vp ; bÞ ¼ rCe ðeeij ; DÞ þ rCp ðvp Þ þ rCd ðbÞ e
ð63Þ
e ij ; DÞ
denotes the elastic strain energy incorporating where rC ðe the effects of the damage variable D on the elastic response of the material, whereas rCp ðvp Þ and rCd(b) are the free energy induced by isotropic hardening and damage development, respectively. The term rCd(b) can be expressed as a simple form with internal variable b [35]: d
pc ðbÞ ¼ 12Cd b
2
ð64Þ
Combining Eqs. (51), (63), and (64), the thermodynamic conjugate force corresponding to b is @c ¼ Cd b @b
ð65Þ
where Cd is a material constant. The other thermodynamic conjugate force Y corresponding to damage variable D can be written as a form of strain energy release rate [50,51], so under the isotropic damage condition, Y may be expressed as Y¼
ð60dÞ
ð61Þ
The current damage state of the material is given by the damage variable D and the internal variable b has been introduced to prescribe the effects of current damage state on the further development of damage, which is similar to the hardening variable of plastic deformation. The thermodynamic conjugate forces corresponding to the internal variables are expressed as
B¼r
where
s0eq ¼
The method adopted for describing damage evolution regulation of rocks mainly concentrates in two aspects: (1) From the angle of phenomenal mechanics and based on statistical distribution theory, the damage variable is assumed to obey the regulation of Weibull distribution or average distribution, and the damage evolution function can be conducted [33,34]. (2) According to the continuum mechanics theory [38], some internal variables are chosen to deal with the damage problem. The former can describe the damage problem in some extent, but the model parameters are usually determined by fitting the experimental data, and their physical meaning are not clear. The latter is based on the irreversible thermodynamic theory, and the energy dissipation generated by plastic deformation or crack development during the damage evolution process is taken into account. The damage is seen as an internal variable in the thermodynamic system and a dissipation potential function is introduced to describe the damage surface. With this method, a systematic damage constitutive framework can be established. Therefore, the latter method is adopted here to analyze the damage evolution law of soft rock. Considering an isothermal condition, it will now employ a scalar isotropic hardening variable vp for plastic deformation, a scalar damage variable D to describe the anisotropic damage states, and another scalar variable b prescribing the current state of damage:
ðs~ Þ2 E0 ð1DÞ
ð66Þ
where ð60eÞ
ð60fÞ
s~ ¼ s~ eq 23ð1þ v0 Þ þ3ð12v0 Þðs~ m =s~ eq Þ
1=2
ð67Þ
The damage evolution equation can be obtained with a damage potential function. These experimental studies have revealed that a limited area, where the damage does not develop
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
and the material presents the elastic behavior, exists in the stress(strain) space [52,53]. Similar to plastic yield surface, some researchers proposed a damage surface Fd = 0 in the thermodynamic conjugate force space and established the relative damage criteria and dissipation potential function. However, there are no systematic results available for the determination of this damage surface so far. In this paper, based on the assumption of the existence of a damage surface and initial damage threshold, one of the simplest convex surfaces for the damage dissipation potential may be given as follows: Fd ðY; B; DÞ ¼ Yd BY0 ¼ 0
ð68aÞ
pffiffiffiffiffiffi Yd ¼ YY
ð68bÞ
where Y0 is the value of Yd corresponding to the initial damage threshold. Substituting Eq. (68) into (56), the incremental evolution equations of D and b can be obtained: dD ¼ dld
@Fd @Y
ð69Þ
db ¼ dld
@Fd @ðBÞ
ð70Þ
where ld Z0 denotes the damage consistence parameter. The loading and unloading conditions may be expressed as follows: (1) If Fd o0, the damage does not begin. Namely the multiplier satisfies dld = 0. (2) If Fd = 0, the new damage occurs. Namely the multiplier dld 40 and it can be obtained by using the damage consistence condition : @Fd dY @Y ð71Þ dld ¼ @F dB d @B db or dld ¼
@Fd =@Y dY Cd ð@Fd =@BÞ
ð72Þ
Based on the irreversible thermodynamic theory, one assumption has been made that the plastic deformation and damage evolution of soft rock are mutually independent considering isothermal condition. So, the Helmholtz free energy of the materials can be written as
rC ¼ rCe ðeeij ; DÞ þ rCp ðvp Þ þ rCd ðbÞ ¼ 12eeij C~ ijkl eekl þ rCp ðvp Þ þ rCd ðbÞ
ð73Þ
where C~ ijkl denotes effective elastic tensor and g(vp) is the potential function related to plastic dissipation. Substituting Eq. (73) into (50), the following equation can be obtained:
sij ¼ r
@C ¼ C~ ijkl eekl ¼ C~ ijkl ðekl epkl Þ @eeij
ð74Þ
The incremental form is dsij ¼ dC~ ijkl ðekl epkl Þ þ C~ ijkl ðdekl depkl Þ or 1 1 deij ¼ C~ ijkl dskl C~ ijkl dC~ ijkl ðeij epij Þ þ depij
1 1 ¼ C~ ijkl dskl þ dC~ ijkl skl þ depij 1
1
¼ C~ ijkl dskl þ
@C~ ijkl dD skl þ depij @D
1
ð75Þ
ð76Þ 1
where C~ ijkl dskl denotes the elastic strain increment, ð@C~ ijkl =@DÞdD skl denotes the damage strain increment and depij denotes the plastic strain increment. That means the total elasto-plastic strain increment of soft rock consists of three parts: elastic strain increment, damage strain increment, and plastic strain increment. Based on the study in the previous section, the plastic part is divided into a plastic friction component dep1 ij and a plastic shear : component dep2 ij p2 depij ¼ dep1 ij þ deij
ð77Þ
A potential function with a similar form of other mudstones [13] is chosen to describe the friction deformation of soft rock for the reason of lacking experimental data about plastic yield characteristics: sm pr ð78Þ Q1 ¼ seq þ aðsm pr Þln I0 where the variable I0 o0 denotes the intersection point between the potential surface and the axis sm, and the parameter a is the slope of the dividing line between the compression region and the dilation region. By using the potential function, the stress space is divided into two zones, corresponding to plastic compressibility and dilation, respectively. The boundary between these two zones is defined by the condition qQ1/qsm =0. In this model, the dividing line function can be expressed as [24] fb ¼ seq aðsm pr Þ ¼ 0
ð79Þ
Combining these two boundary functions, the expression of I0 can be obtained. Hence, considering the damage effect and nonassociated flow rule, the effective plastic potential function and dividing line can be expressed as follows: s~ m pr ð80Þ Q~ 1 ¼ s~ eq þ aðs~ m pr Þln I0 fb ¼ s~ eq aðs~ m pr Þ ¼ 0
4.3. Elasto-plastic damage constitutive law of soft rocks
391
ð81Þ
where s~ m ¼ sm =ð1DÞ, s~ eq ¼ seq =ð1DÞ. For plastic pore deformation, the effective plastic potential function is given formally similar to the yield criterion: " # s2eq 3sm ¼0 ð82Þ Q~ 2 ¼ 0 2 þ 2f cosh 2s0eq s eq where the parameters are the same as above. According to the non-associated rule of orthogonal flow, the following two equations can be obtained: dep1 ij ¼
dl1 @Q~ 1 @sij
ð83Þ
dep2 ¼ ij
dl2 @Q~ 2 @sij
ð84Þ
The plastic multiplier dl can be determined by the consistency principle dF~ 1 ¼
@F~ 1 @F~ 1 @F~ 1 @F~ 1 @F~ 1 dD þ dD dsij þ dZ ¼ dsij þ @sij @D @Z @sij @D @F~ 1 @Z p @F~ 1 @F~ 1 @F~ 1 @Z @Q~ 1 þ dD þ de ¼ dsij þ dl1 ¼0 @Z @eps s @sij @D @Z @eps @seq ð85Þ
@F~ 2 @F~ 2 @F~ 2 @F~ 2 @F~ 2 dD þ df ¼ dF~ 2 ¼ dsij þ dZ0 þ ds @sij @D @Z0 @f @sij ij
ARTICLE IN PRESS 392
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
þ
¼
@F~ 2 @D
dD þ
" p2 # @F~ 2 @Z0 @F~ 2 6sin j sij deij p2 p2 ð1f Þd d e þ e v @Z0 @e p2 @f 3sin j s0eq
@F~ 2 @F~ 2 @F~ 2 @Z0 1 @Q~ dD þ dsij þ sij dl2 2 0 0 p2 @sij @D @Z @e ð1f Þseq @sij " # @F~ 2 @Q~ 2 6sin j sij @Q~ 2 ¼0 ð1f Þdl2 þ dl2 @f @seq 3sin j s0eq @sij
ð86Þ
So,
dl1 ¼
@F~ 1 ds þ @F~ 1 dD ij @sij @D
dl2 ¼
ð87Þ
H1 @F~ 2 ds þ @F~ 2 dD ij @sij @D
ð88Þ
H2
where H1 ¼
@F~ 1 @Z @Q~ 1 @Z @eps @seq
ð89aÞ
@F~ 2 @Z0 1 @Q~ @F~ sij 2 2 0 0 p2 @Z @e ð1f Þseq @sij @f " # ~ ~ s @Q 2 6sin j ij @Q 2 ð1f Þ @seq 3sin j s0eq @sij
H2 ¼
ð89bÞ
Substituting dl1 and dl2 into Eqs. (83) and (84), two types of plastic deformation increment are obtained. During a general loading history, the plastic deformation of soft rock is determined by two plastic deformation mechanisms together (see Fig. 3). The present stress satisfies the compression state. In Fig. 3, the plane seq sm is divided into three areas by failure surface Fmax, friction yield surface F1 and pore deformation yield surface F2, including the elastic area A0, the plastic area A1, which is only relative to pore deformation mechanism and the plastic area A2, which is relative to two deformation mechanism simultaneously. At the initial stage of loading, the stress point falls into A0. With the increase of seq, the stress point reaches A1, and it represents the plastic pore deformation. And with continuous pore collapse, the contact area between grains will be increased. It may cause the enhancement of inner friction. So the friction mechanism is activated and plays a role in the pore deformation mechanism simultaneously, and the stress point reaches the A2 area; then more ligament fractures of skeleton occur with the further increase of seq and the deformation state of soft rock is basically determined by the transfixion of a fracture plane and development of a shear zone. Now the plastic friction mechanism plays a more important role than pore deformation mechanism, and the failure surface is reached finally.
eq Fmax: Failure Surface
O
m
Fig. 3. The plastic deformation determined by two plastic deformation mechanisms.
Three distinct constitutive domains can be identified as follows: (1) If F~ 1 o0; F~ 2 o 0, the applied stress state is fully inside the elastic domain or leads to an elastic unloading. No plastic flow occurs and therefore dl1 =dl2 = 0, namely depij ¼ 0. (2) If F~ 2 ¼ 0; dF~ 2 ¼ 0, but F~ 1 o 0 or F~ 1 ¼ 0; dF~ 1 o0, then the plastic pore deformation mechanism is activated while the plastic friction yield surface is not reached. The plastic multiplier dl2 40 is determined by Eq. (88) and dl =0. Namely, the plastic strain increment can be expressed as . depij ¼ dep2 ij ~ (3) If F 1 ¼ 0; dF~ 1 ¼ 0 and F~ 2 ¼ 0; dF~ 2 ¼ 0, the two plastic deformation mechanisms are both activated. The plastic multipliers dl1 40 and dl1 40 can be determined by Eqs. (87) and (88) together. The plastic strain increment can be expressed p2 as depij ¼ dep1 ij þdeij .
5. Verification of the proposed model The proposed model involves several parameters, which can be determined step by step from direct interpretation of triaxial compression test results of soft rocks and other similar mudstones. The two elastic parameters, E0 and v0, are obtained from the initial linear part of the stress–strain curve obtained from a conventional triaxial compression test. The values of E0 may depend on the confining pressure. This phenomenon is considered in this study. The cohesion C and inner friction angle j are obtained by drawing the locus of peak stresses in the p–q plane. The hardening parameters, Z0,Zm, R are determined by the study of other similar mudstones [13], and the parameter k is influenced intensively by the initial confining pressure. After fitting the test data, k is found to satisfy the logistic growth function: k =23.03 21.53/ (1+(sm0 /0.76)4.84), where sm0 is the initial confining pressure. The value of parameter a can be obtained by identifying the stress point corresponding to depv ¼ 0 in the stress–strain curve. By defining the locus of stresses at the onset of elastic strain as the initial damage stress point, and substituting this value into Eqs. (66) and (67), the initial damage threshold can be obtained. The value of material constant Cd can be determined by the experimental data of other elasto-plastic materials [35]. To verify the capacity of the proposed model for describing the mechanical behavior of soft rock, numerical simulations of several groups of triaxial compression tests by using the developed model are presented. These tests are performed on saturated silty mudstone, a kind of soft rock from the southern coastal region of China. The following values of parameters are used in Table 1. Note that f0 is the porosity of soft rock bearing no loading, ignoring the porosity variation during the process of applying the confining pressure in triaxial tests. Figs. 4–7 present the simulations of triaxial compression tests at four levels of confining pressures: 0, 0.5, 1.0, and 3.0 MPa. For convenience, compression is assumed to be positive. From Figs. 4–7, it can be seen that good agreement between the numerical simulations and experimental data is obtained. The mechanical behaviors of soft rock, such as plastic deformation, material damage, and transition from plastic compressibility to dilation, are correctly predicted by the proposed model. It can also be seen that the mechanical behavior is sensitive to confining pressure. At low confining pressures, say 0, 0.5, and 1 MPa, a brittle failure is observed, with a small softening regime and the material failure occurs with a small value of the axial
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
393
Table 1 Values of parameters in the model. E0 (MPa) Confining pressure 0 MPa
0.5 MPa
1 MPa
3 MPa
610
1280
1740
1765
v0
C (MPa)
j (deg)
Z0
Zm
R
a
f0 (%)
Cd
0.35
1.61
30.16
0.017
1
0.00008
1.41
15.1
5
Fig. 4. Simulation of a triaxial compression test on saturated soft rock with 0 MPa confining pressure.
Fig. 5. Simulation of a triaxial compression test on saturated soft rock with 0.5 MPa confining pressure.
strain. This means with the increase of seq, the friction mechanism is activated instead of the pore deformation mechanism in a short time. However, at a slightly higher confining pressure, say 3 MPa, the axial strain corresponding to the failure state is much higher, which means that the plastic friction yield surface is reached at a higher axial stress and the pore deformation mechanism plays a more important role in the whole process.
Fig. 6. Simulation of a triaxial compression test on saturated soft rock with 1 MPa confining pressure.
Fig. 7. Simulation of a triaxial compression test on saturated soft rock with 3 MPa confining pressure.
A falling phenomenon is shown at the last stage of the stress– strain curve. This means the transfixion of the fracture plane and development of the shear zone in soft rock finally causes the failure of material, and this strain softening behavior generally implies the formation of strain localization mode. Suitable regularization methods like non-local formulation and secondgradient modelling should be employed. However, this feature is beyond the scope of study of the present work and an asymptotic
ARTICLE IN PRESS 394
C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
at various confining pressure conditions. It has been shown that the proposed model correctly describes the main features of the mechanical behavior of saturated soft rock. In addition, the related results of other mudstone have been used in the determination of model parameters. The comparisons between the simulations and experimental data presented in this paper are only represent the first phase of validation of the model. Further advanced tests are still needed to adjust these parameters and check the general validity of the proposed model.
Acknowledgments
Fig. 8. Damage evolution curves at different confining pressures (0, 1, and 3 MPa).
perfect plastic regime is obtained to reflect the post-peak failure behavior of soft rock. The results of numerical simulations are appreciably greater than experimental data, mainly occurring after the yield locus. It is because the pore water pressure is not taken into account in the present model, inducing the effective stress of soft rock for numerical calculations greater than the actual situation. Nevertheless, the differentiation is still in the range of allowable error. In Fig. 8, the damage evolution curves at different confining pressures (0, 1, and 3 Mpa) are presented. The amount of elastic damage only happened at a small percentage of total structural damage of soft rock, and most parts of the damage took place during the plastic deformation phase, which mainly occurred after the activation of plastic friction mechanism. In this phase, the development of various microscopic defects such as micropores and microcracks can all be considered to be the performance of damage evolution. Moreover, by comparing the curves at different levels of confining pressures, it can also be found that the amount of damage at higher confining pressure is more than the one at lower confining pressure. That means the increase in confining pressure can strengthen the structure and reduce the amount of total damage of frozen soil.
6. Conclusion In this paper, an elasto-plastic damage constitutive model based on irreversible thermodynamics theory and damage mechanics theory is developed for saturated soft rock. Two plastic deformation mechanisms are taken into account: plastic friction and plastic pore deformation. By incorporating the volumetric deformation effect and considering the motion of yield function in stress space governed by its center position and rotation hardening rule, the plastic friction yield criterion is proposed. Meanwhile the plastic pore yield criterion is established as an adaptation of Gurson’s criterion using the proposed friction yield model to determine the matrix deformation of porous materials. Moreover, a damage variable describing the development of various microscopic defects is defined. Based on this variable, an elasto-plastic damage constitutive model for saturated soft rock is proposed for taking into account the two plastic mechanisms. The non-associated flow rule is adopted in the model, and the plastic strain increment-type expression is obtained. Numerical simulations have been compared with experimental data from triaxial compression tests that were performed on the saturated samples
The research is partly supported by the Natural Science Foundation of China (NSFC) through Grant no. 40672194, Guangdong Natural Science Foundation through Grant no. 06104932 and Doctorial Discipline Foundation through Grant no. 20060558060. The authors gratefully acknowledge the contribution of Engineer You-Rou Wang in completing the experiment. References [1] Mondo NH, Bjorlykke K, Jahren J. Experimental mechanical compaction of clay mineral aggregates—changes in physical properties of mudstones during burial. Mar. Petrol Geol. 2007;24:289–311. [2] Koncagul EC. Comparison of uniaxial compressive strength test and slake durability index values of shale samples from the breathtt formation, Kentucky. PhD thesis, University of Missouri, Rolla, 1998. [3] Botts ME. The effects of slaking on the engineering behavior of clay shale. PhD thesis, University of Colorado, Boulder, 1986. [4] Yamaguchi H, Yoshida K, Kuroshima I, Fukuda M. Slaking and shear properties of mudstone. In: Rock mechanics and power plants. Rotterdam: Balkema; 1988. p. 133–44. [5] Yoshinaka R, Tran TV, Osada M. Non-linear, stress and strain-dependent behavior of soft rocks under cyclic triaxial conditions. Int J Rock Mech Min Sci 1998;35:941–55. [6] Ye GL, Naito K, Sawada K. Experimental study on soft sedimentary rock under plane strain compression and creep tests. In: Proceedings of the third Asian rock mech symp, 2004. p. 865–70. [7] Singh TN, Verma AK. Prediction of creep characteristic of rock under varying environment. Environ Geol 2005:559–68. [8] Eseme E, Littke R, Krooss BM. Factors controlling the thermo-mechanical deformation of oil shales: implications for compaction of mudstones and exploitation. Mar Petrol Geol 2006;23:715–34. [9] Corkum AG, Martin CD. The mechanical behaviour of weak mudstone (Opalinus Clay) at low stresses. Int J Rock Mech Min Sci 2007;44:196–209. [10] Szymakowski J. Shear testing of soft rock masses. Geotech Test J 2007; 30:133–40. [11] Billiotte J, Yang D, Su K. Experimental study on gas permeability of mudstones. Phys Chem Earth 2008;33:231–6. [12] Kodaka T, Oka F, Kitahara H, Ohta H, Otani J. Observation of compaction bands under triaxial conditions for diatomaceous mudstone. In: Hyodo M, Murata H, Nakata Y, editors. Proc int symp geomech geotech particulate media, Ube Yamaguchi, Japan, 12–14 September 2006. Rotterdam: Balkema; 2006. p. 69–75. [13] Shao JF, Zhu QZ, Su K. Modeling of creep in rock materials in terms of material degradation. Comp Geotech 2003:549–55. [14] Li Hangzhou, Liao Hongjian, Sheng Qian. Study on statistical damage constitutive model of soft rock based on unified strength theory. Chin J Rock Mech Eng 2006;25:1331–6. [15] Dashnor H, Albert G, Francoise H, Christophe A. Saturated and unsaturated behaviour modelling of Meuse–Haute/Marne argillite. Int J Plast 2007;23: 733–66. [16] Adachi T, Oka F. An elasto-plastic constitutive model for soft rock with strain softening. Int J Numer Anal Meth Geonmech 1995;19:233–47. [17] Jun Yinglian, Tai Quanzhou, Yuan Hua, Guo Liangdai. Reliability analysis of railway tunnel construction process within soft and weak rock mass using the unified elasto-plastic strength criterion. Key Eng Mater 2007:2521–4. [18] Lade PV. Modelling the strengths of engineering materials in three dimensions. Mech Cohes-Frictional Mater 1997;2:339–56. [19] Desai CS. Mechanics of materials and interfaces: the disturbed state concept. Boca Raton: CRC Press; 2001. [20] Peric D, Ayari MA. On the analytical solutions for the three-invariant Camclay model. Int J Plast 2002;18:1061–82. [21] Schroeder C. Du coccolithe au re´servoir pe´trolier: approche phe´nome´nologique du comportement me´canique de la craie en vue de sa mode´lisation a´ diffe´rentes e´chelles. PhD thesis, Univ Lie´ge, France, 2003.
ARTICLE IN PRESS C.Y. Zhou, F.X. Zhu / International Journal of Rock Mechanics & Mining Sciences 47 (2010) 385–395
[22] Zhang H, Zhou L, Huang H. Strain localization analysis of anisotropic saturated and partially saturated porous media. Rock Soil Mech 2004;25: 675–80. [23] Dai Yonghao, Chen Weizhong, Wu Guojun, Zhou Xide. Study on elastoplastic damage model of unsaturated rock mass and its application. Chin J Rock Mech Eng 2008;27:728–35. [24] Xie SY, Shao JF. Elastoplastic deformation of a porous rock and water interaction. Int J Plast 2006:2195–225. [25] Jeong HY. A new yield function and a hydrostatic stress-controlled void nucleation model for porous solids with pressure-senstitive matrices. Int J Solids Struct 2002;39:1385–403. [26] Gurson AL. Porous rigid-plastic materials containing rigid inclusions-yield function: plastic potential and void nucleation. In: Taplin DMR, editor. Proc int conf fract. Oxford: Pergamon; 1977. p. 357–64. [27] Gurson AL. Continuum theory of ductile rupture by void nucleation and growth: part I—yield criteria and flow rules for porous ductile media. J Eng Mater Technol 1977;99:2–15. [28] Gurson AL. Plastic flow and fracture behavior of ductile materials incorporating void nucleation, growth and interaction. PhD thesis, Brown Univ, Providence, RI, 1975. [29] Tvergaard V. Material failure by void growth to coalescence. Adv Appl Mech 1990;27:83–147. [30] Perrin G, Leblond JB. Accelerated void growth in porous ductile solids containing two populations of cavities. Int J Plast 2000;16:91–120. [31] Leblond JB, Perrin G. Introduction a´ la mcanique de la rupture ductile des me´taux. Lecture notes. Paris: Ecole Polytechnique; 1996 [in French]. [32] Fu Zhiliang. Nonlinear creep damage model for soft rock under uniaxial compression. J Coal Sci Eng 2006;12:37–9. [33] Cao Wengui, Zhao Minghua, Liu Chengxue. Study on the model and its modifying method for rock softening and damage based on Weibull random distribution. Chin J Rock Mech Eng 2004;23:3226–31. [34] Xu Weiya, Wei Lide. Study on statistical damage constitutive model for rock. Chin J Rock Mech Eng 2002;21:787–91. [35] Hayakawa K, Murakami S, Liu Y. An irreversible thermodynamics theory for elastic–plastic-damage materials. Eur J Mech A Solids 1998;17:13–32. [36] Lemaitre J. A course of damage mechanics, 2nd ed. Berlin: Springer; 1992. [37] Murakami S. Anisotropic damage in metals. In: Boehler JP, editor. Failure criteria of structured media. Rotterdam: Balkema; 1993. p. 99–119.
395
[38] Krajcinovic D, Lemaitre J. Continumn damage mechanics—theory applications. Vienna: Springer; 1986. [39] Lemaitre J, Chaboche JL, Shrivastava B. Mechanics of solid materials. Cambridge: Cambridge University Press; 1990. [40] Collins IF, Kelly PA. A thermomechanical analysis of a family of soil models. Geotechnique 2002;52:507–18. [41] Collins IF, Muhunthan B. On the relationship between stress—dilatancy, anisotropy, and plastic dissipation for granular materials. Geotechnique 2003;53:611–8. [42] Collins IF, Hilder T. A theoretical framework for constructing elastic/plastic constitutive models of triaxial tests. Int J Numer Anal Meth Geomech 2002; 26:1313–47. [43] Shen Zhujiang. Summary on the failure criteria and yield functions. Chin J Geotech Eng 1995;17:1–8. [44] Collins IF, Houlsby GT. Application of thermomechanical principles to the modelling of geotechnical materials. Proc R Soc London 1997: 1975–2001. [45] Chiarelli AS. Experimental investigation and constitutive modeling of coupled elastoplastic damage in hard claystones. PhD thesis, Univ Lille, France, 2000. [46] Jeong HY, Pan J. A macroscopic constitutive law for porous solids with pressure-sensitive matrices and its implications to plastic flow localization. Int J Solids Struct 1995;32:3669–91. [47] Zhao Jisheng, Tao Xiaxin. A constitutive model of the porous media based on micro-mechanics damage. J Xian Min Inst 2002;17:73–8. [48] Yang qiang, Chen xin, Zhou weiyuan. Elasto-plastic damage model for geomaterials and strain localization analyses. Chin J Rock Mech Eng 2004; 23:3577–83. [49] Timoshenko SP, Goodier JN. Theory of elasticity. New York: McGraw-Hill; 1951. [50] George ZV, Deliktas B. A coupling anisotropic damage model for the inelastic response of composite materials. Comput Meth Appl Mech Eng 2000;183: 159–99. [51] Abu Al-Rub RK, Voyiadjis GZ. On the coupling of anisotropic damage and plasticity models for ductile materials. Int J Solids Struct 2003;40:2611–43. [52] Holcomb DJ, Costin LS. Detecting damage surface in brittle materials using acoustic emission. J Appl Mech 1986;53:536–44. [53] Holcomb DJ. Discrete memory in rock: a review. J Rheol 1994;28:725–58.