# A Creep Damage Constitutive Model for a Rock Mass with Nonpersistent Joints under Uniaxial Compression

Authors: Shuran Lv; Wanqing Wang; Hongyan Liu
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

As part of the rock mass, both the mesoscopic and macroscopic flaws will affect the creep mechanical behavior of the rock mass with nonpersistent joints. This study focuses on this kind of rock mass and establishes a creep damage model to account for the effect of the joint on its creep mechanical behavior. First, on basis of analyzing the rock element creep mechanism and the typical creep deformation curve, a new creep damage constitutive model for the rock element is set up by introducing the damage theory and Kachanov damage evolution law into the classic creep constitutive model such as J body model. Second, the determination method of the proposed model parameters is studied in detail. Third, the calculation method of the macroscopic damage caused by the joint proposed by others is introduced which can consider the joint geometry, strength, and deformation parameters at the same time. Finally, the creep damage model for the rock mass with nonpersistent joints under uniaxial compression is proposed. The calculation examples indicate that it can present the effect of the joint on the rock mass creep mechanical behavior.


---

## Body

## 1. Introduction

The time-dependent deformation behavior and failure process of the rock mass has a very significant effect on the stability and failure of the rock slopes and tunnels [1, 2]. This deformation might begin few minutes after excavation and continue for up to several years. Thus, the identification and modeling of the rock mass time-dependent mechanical behavior is very important in rock engineering. In general, the rock mass time-dependent mechanical behavior includes creep, stress relaxation, and elastic aftereffect, among which the rock mass creep mechanical behavior is one of the most important and challenging issues in many practical rock projects [3, 4]. Meanwhile, the rock mass is always intersected by the joints (here we call them the macroscopic flaws) into the rock blocks (here we call them the rock element) and at the same time the rock element also contains many microcracks (here we call them the microscopic flaws). As is known to all, these two different scale flaws both affect the jointed rock mass strength and deformability. Now because of its complexity, the existing studies on the creep mechanical behavior mainly focused on the rock element involving theoretical, experimental, and numerical methods [5–7]. Especially, the study on the creep constitutive model for the rock element is one of the most fundamental and important issues in the rock mass rheological mechanics which attracts much attention [8–10]. In general, these models can be divided into the following two categories according to their establishment method. One is the empirical model which is set up by fitting with the experimental creep curve of the rock element [11], in which the proposed models were in good agreement with the practical engineering. However, they varied greatly with different engineering backgrounds and were not universal. The other is the theoretical model which is mainly established by combining the following three ideal elements namely the Hooke spring for an elastic solid, the Newton dashpot for a viscous fluid, and the St. Venant frictional element for a perfect plastic body in series or parallels [12]. From this viewpoint, many creep models for the rock element are set up, for example, Maxwell model, Kelvin model, Burgers model, and Bingham model. Although they can perfectly simulate the rock element deformation behavior from the instantaneous elastic deformation to the steady-state creep deformation, these models cannot reflect the accelerating creep deformation. In many engineering cases, the accelerating creep is often observed before the rock mass rheological failure occurs [3, 13]. Therefore, many works have been carried out in order to perfectly reflect the accelerating creep behavior of the rock element. In sum, these works can also be classified into two categories namely the nonlinear rheological model and the damage rheological model. Firstly, researchers have set up many nonlinear rheological models to describe the accelerating creep stage. Yang et al. [14] proposed a new nonlinear visco-elasto-plastic creep model with creep threshold and long-term strength by connecting an instantaneous elastic Hooke body, a visco-elasto-plastic Schiffman body, and a nonlinear visco-plastic body in series, which could describe the typical creep behavior including the primary creep stage, the secondary creep stage, and the tertiary creep stage. Jiang et al. [15] modified the quasi-static Nishihara model for transient and steady-state creep by adding a strain-triggered inertial element in series, which can simulate a phase of quadratic accelerating creep. Zhao et al. [16] proposed a nonlinear elasto-visco-plastic rheological model by connecting a Hooke body, a parallel combination of Hooke and plastic slide bodies, a Kelvin body, and a generalized Bingham body, which can perfectly describe the loading and unloading creep behavior and the full stages of creep, especially the tertiary creep stage. Besides adopting the nonlinear components to simulate the accelerating creep process, many researchers assumed it to be a rapid coalescence and growth of the microcracks, which can be called the damage evolution. Therefore, the damage mechanics was adopted to study the creep property of the rock element especially the accelerating creep process. For instance, Golshani et al. [17] took into account the time-dependent behavior of brittle material to extend the micromechanics-based damage model proposed by Golshani et al. [18] and established a creep damage model in capable of reproducing the three stages of creep behavior. By introducing a damage factor into the Burgers model, Yang et al. [19] proposed a new nonlinear creep damage model to describe the entire creep process of coal specimens under three-dimensional stress, which agrees very well with the experimental results. Based on the nonlinear damage creep characteristics of the rock element and damage variable, Cao et al. [20] proposed a new nonlinear damage creep constitutive model of high-stress soft rock element, which can reasonably describe the soft rock creep deformation.On basis of the existing research studies on the creep mechanical behavior of the rock element, some investigations on that of the rock mass is also in progress. Wong et al. [21] experimentally investigated the creep damage and failure mechanism around an opening in rock-like material containing nonpersistent joints and found that λ=σ3/σ1 had much effect on the rock mass creep failure mode, time, and stability. Through the laboratory experiment on the sandstone samples with nonpersistent joints subjected to shear creep, Wu et al. [22] developed a new creep model considering the effect of joint length on the rock time-dependent behavior. Li et al. [6] developed a damage rheology model for jointed rock masses in which the damage caused by the joints to the rock mass is described with a second-order tensor.Overall, more and more creep models for the rock element and jointed rock mass are proposed which can perfectly simulate the whole creep deformation process especially the accelerating creep one. However, some improvements should be made to the previous models. First, in order to simulate the whole creep deformation process especially the accelerating creep stage, one approach is to adopt more and more basic components in series or parallels, which will lead to the difficulty in solving the constitutive equation. The other approach is to introduce the damage variable into the creep equation, which has no corresponding physical background and theoretical basis. Second, the determination of the model parameters in the previous models is too difficult. Some of them are determined only by experience, fitting with the experimental curves or other complicated theory, which greatly limits their practical applications [23]. Third, the study on the creep damage model mainly focuses on the rock element and that for the jointed rock mass is very rare. Hence, Kachanov creep damage theory is firstly introduced into the classic creep model such as J body model to set up the creep damage model for the rock element which can simulate the accelerating creep stage. Second, the determination of the model parameters in the proposed models is studied, and then the validity of the proposed model for the rock element is verified with the experimental data. Finally, based on the coupled damage theory, the creep damage model for the rock mass with nonpersistent joints is proposed.


## 2. Rock Element Creep Mechanism and Creep Deformation Curve

In order to perfectly study the rock element creep mechanism and its creep deformation characteristic, lots of rock element creep deformation experiments are made by the researchers [14, 24], and the typical creep deformation curve is shown in Figure 1. It can be seen that the rock element creep deformation curves can be classified into the following two types. One is the steady creep deformation for instance curve I in which the creep rate gradually decreases to 0, and accordingly the creep deformation gradually tends to be constant when the applied load is less than the rock element long-term strength. The other is the unsteady creep deformation for instance curve II in which the creep rate gradually increases, and accordingly the creep deformation gradually increases and eventually the rock element failure occurs when the applied load exceeds the rock long-term strength. Therefore, the creep strain is usually described as a unitary process represented by the sequence of specific time-dependent deformations: the instantaneous elastic strain, the transient creep strain, the steady creep strain, and the accelerating creep strain [12], shown as curve II in Figure 1. Now the previous studies mainly focus on the transient and steady creep stages; however, in fact the accelerating creep stage is more significant and assumed to be much related to the rock stability.Figure 1
The typical creep curve of the rock element. It has four typical stages: (1) The instantaneous elastic strain (OA): when the stress is applied on the rock element, an instantaneous elastic strainε0 is generated immediately, which is often modeled by Hooke’s law. (2) The transient creep strain (AB): the strain increases with time, but the strain rate gradually decreases to a constant. (3) The steady creep strain (BC): the strain rate nearly keeps constant, and the strain increases linearly with time. (4) The accelerating creep strain (CD): the strain rate increases rapidly with time and finally followed by the failure of the rock element.Therefore, through the analysis of Figure1 and the previous research studies, the following conclusions can be obtained.(i)
Many creep models, for example, the Maxwell model, the Kelvin model, the J body model, and the generalized Kelvin model can only simulate the first three stages of the creep curve and cannot simulate the accelerating creep stage. It is assumed that the components adopted in the above models are linear, namely, the elastic modulus and viscous coefficient keep constant during the whole creep deformation process. However, during the accelerating creep stage, the nonlinearly permanent plastic deformation will occur which indicates that the irreversible damage happens. Therefore, the component which can reflect the rock element creep damage should be adopted to describe the accelerating creep stage.(ii)
The steady creep stage of the rock element in practice is approximately linear, and therefore its corresponding deformation rate is assumed to be the average of the curve’s slope.(iii)
During the accelerating creep stage, the rock element damage will become worse and worse and eventually fail. The timetF when the rock element creep failure occurs is called the rock element creep lifetime where the slope of the creep curve is infinite and the rock element damage is 1.


## 3. A Creep Damage Constitutive Model for a Rock Element Based on J Body Model

### 3.1. Establishment of the Rock Element Creep Damage Constitutive Model

As stated above, researchers have set up many classic rock element creep constitutive models such as Maxwell model, Kelvin model, and J body model. Here, J body model is taken for an example. It is composed of an elastic body and a Maxwell body in series, shown in Figure2(a). Assume the elastic moduli of these two elastic components and the viscous coefficient of the viscous component are E1, E2, and η, respectively. According to the combination relationship of these three components shown in Figure 2(a), the following stress-strain relationship can be obtained:(1)σ=σ1+σ2,ε=ε1=ε2,ε1=σ1E1,ε˙2=σ˙2E2+σ2η,where σ1andε1 are the stress and strain of the elastic component, respectively; σ2andε2 are the stress and strain of the Maxwell body, respectively; and σandε are the stress and strain of the J body, respectively.Figure 2
The representations of the rock element creep model with damage. (a) The original J body model. (b) The stress-strain relation of the J body model. (c) The improved J body model with damage.
(a)
(b)
(c)Therefore, the constitutive model for the J body can be obtained from equation (1):(2)E1+E2E2ε˙+E1ηε=ση+σ˙E2.Because the stress is constant, namely,σ=σ0=const on the condition of creep deformation, the creep equation for the J body can be derived from equation (2):(3)ε=σ0E1+E2−σ0E1exp−E1E2tηE1+E2+σ0E1.From its corresponding creep curve shown in Figure2(b), we find that the instantaneous elastic strain is σ0/E1+E2 at time t = 0 s, and then with time going on, the strain becomes larger and larger and finally tends to be a constant value σ0/E1. It can be seen that the J body model can reflect the four typical stages of the whole creep deformation except for the accelerating creep stage. This is because the irreversible plastic damage and failure will occur in the rock element at the accelerating creep stage. It cannot be simulated with the J body model because all the components in this model are linear. As stated above, the often adopted methods to set up the rock element creep model are classified into the following two types. One is to increase the number of the components, for example, the famous traditional Nishihara model composed of a Hooke body, Kelvin body, and Bingham body in series [25]. Xu et al. [26] proposed a nonlinear visco-elasto-plastic rheological model (Hohai Model) for a rock with seven components in series or parallels. However, with increase in the number of the components in the model, its corresponding constitutive equation also becomes more and more complicated and is not suitable for practical engineering application. The other is to introduce the damage to the component. It assumes that the component’s elastic modulus or viscous coefficient is not constant anymore, and it will decrease with the creep deformation [19, 20]. Therefore, in view of this idea, the damage is introduced into the J body model to reflect the plastic failure occurred in the accelerating creep stage. The elastic component E1 is revised to the damage one. So, the improved J body model is assumed to be composed of a Maxwell body and a damage body E1t in parallels, shown in Figure 2(c).Then, the next task is to obtain the damage evolution law of the damage componentE1t. From Lemaitre strain equivalence hypothesis [27], we obtain(4)E1t=E11−Dt,where Dt is the damage variable, a function of the time. The other parameters are stated as above.Kachanov [28] investigated the creep damage and assumed that the rock element creep damage rate could be expressed as(5)dDtdt=Bσ1−Dtn,where B and n are the material constants of the rock element.Whent = 0, the rock element is undamaged, Dt=0. By integrating equation (5), we obtain(6)Dt=1−1−B1+nσnt1/n+1.Assume the critical failure time of creep damage istF and the corresponding damage variable Dt of the rock element is 1. Then, we obtain(7)tF=1B1+nσn.From equations (6)∼(7), we obtain(8)Dt=1−1−ttF1/n+1.Then, the constitutive equation of the damage body can be obtained:(9)ε=σ0E11−ttF−1/n+1.Finally, we obtain the creep constitutive equation of the J body model with damage(10)ε=σ0E11−Dt+E2−σ0E11−Dt⋅exp−E11−DtE2tηE11−Dt+E2+σ0E11−Dt.


### 3.2. Determination of the Parameters in the Proposed Model

The proposed creep damage constitutive model based on the J body model has five parameters namelyE1, E2, n, tF, and η, and how to determine them is also important. Now the often adopted methods are fitting with the experimental data, which are not only complicated but also cannot clearly explain their physical meaning. Therefore, it is very necessary to seek other methods to determine them. Now the determination of these parameters is discussed from the physical meaning of the proposed model.


#### 3.2.1. The Critical Creep Failure TimetF

As stated above, the creep damage of the rock element is 1 at the critical creep failure timetF, and the corresponding creep strain and strain rate both tend to be infinite. As shown in Figure 1, the critical creep fracture time tF is easily obtained from the rock element creep experimental curve.


#### 3.2.2. The Elastic Moduli ofE1 and E2

From equation (10), it can be seen that Dt=0 at the beginning time t = 0 s, and then under the stress σ0, the initial elastic strain ε0 of the proposed model is(11)ε0=σ0E1+E2.After the rock element goes into the steady creep stage, the creep curve exhibits a linear increase and no damage occurs in the rock element. Therefore, the elastic modulus of the elastic componentE1 is still E1. From the structure characteristic of the J body, it can be seen that the intersection point of the oppositely elongated line of the steady creep stage with the strain axial is the final strain ε1 of the elastic component E1 (shown in Figure 1), namely,(12)ε1=σ0E1.Combining equations (11) and (12), the elastic moduli E1 and E2 of these two elastic bodies can be solved.


#### 3.2.3. Determination ofn

From the rock element creep deformation process and its mechanism, it can be seen that the viscoelastic deformation is finished after the rock element goes into the accelerating creep stage, and only the permanent deformation occurs with time going on. Then, from equation (9), the permanent strain εr of the sample is(13)εr=σ0E11−ttF−1/n+1.Assumem points (ti, εi) (i = 1, 2, …, m) are selected from the accelerating creep stage, and at the same time σ0 is known and tF and E1 have been solved; therefore, the calculation method of ni is(14)ni=ln1−ti/tFlnσ0/E1εri−1.The average value ofni is assumed to be the value of n:(15)n=1m∑i=1mni.


#### 3.2.4. Determination ofη

The strain-time curve of the viscous component is linear which indicates that the viscosity is mainly reflected in the stable creep stage. Assumem points (ti, εi) (i = 1, 2, …, m) are selected from the stable creep stage, and there is(16)εi=σ0E11−Dti+E2−σ0E11−Dti⋅exp−E11−DtiE2tiηiE11−Dti+E2+σ0E11−Dti.So the corresponding viscosity coefficientηi (i = 1, 2, …, m) can be obtained. Finally, the average value of ηi is assumed to be the value of η(17)η=1m∑i=1mηi.


## 3.1. Establishment of the Rock Element Creep Damage Constitutive Model

As stated above, researchers have set up many classic rock element creep constitutive models such as Maxwell model, Kelvin model, and J body model. Here, J body model is taken for an example. It is composed of an elastic body and a Maxwell body in series, shown in Figure2(a). Assume the elastic moduli of these two elastic components and the viscous coefficient of the viscous component are E1, E2, and η, respectively. According to the combination relationship of these three components shown in Figure 2(a), the following stress-strain relationship can be obtained:(1)σ=σ1+σ2,ε=ε1=ε2,ε1=σ1E1,ε˙2=σ˙2E2+σ2η,where σ1andε1 are the stress and strain of the elastic component, respectively; σ2andε2 are the stress and strain of the Maxwell body, respectively; and σandε are the stress and strain of the J body, respectively.Figure 2
The representations of the rock element creep model with damage. (a) The original J body model. (b) The stress-strain relation of the J body model. (c) The improved J body model with damage.
(a)
(b)
(c)Therefore, the constitutive model for the J body can be obtained from equation (1):(2)E1+E2E2ε˙+E1ηε=ση+σ˙E2.Because the stress is constant, namely,σ=σ0=const on the condition of creep deformation, the creep equation for the J body can be derived from equation (2):(3)ε=σ0E1+E2−σ0E1exp−E1E2tηE1+E2+σ0E1.From its corresponding creep curve shown in Figure2(b), we find that the instantaneous elastic strain is σ0/E1+E2 at time t = 0 s, and then with time going on, the strain becomes larger and larger and finally tends to be a constant value σ0/E1. It can be seen that the J body model can reflect the four typical stages of the whole creep deformation except for the accelerating creep stage. This is because the irreversible plastic damage and failure will occur in the rock element at the accelerating creep stage. It cannot be simulated with the J body model because all the components in this model are linear. As stated above, the often adopted methods to set up the rock element creep model are classified into the following two types. One is to increase the number of the components, for example, the famous traditional Nishihara model composed of a Hooke body, Kelvin body, and Bingham body in series [25]. Xu et al. [26] proposed a nonlinear visco-elasto-plastic rheological model (Hohai Model) for a rock with seven components in series or parallels. However, with increase in the number of the components in the model, its corresponding constitutive equation also becomes more and more complicated and is not suitable for practical engineering application. The other is to introduce the damage to the component. It assumes that the component’s elastic modulus or viscous coefficient is not constant anymore, and it will decrease with the creep deformation [19, 20]. Therefore, in view of this idea, the damage is introduced into the J body model to reflect the plastic failure occurred in the accelerating creep stage. The elastic component E1 is revised to the damage one. So, the improved J body model is assumed to be composed of a Maxwell body and a damage body E1t in parallels, shown in Figure 2(c).Then, the next task is to obtain the damage evolution law of the damage componentE1t. From Lemaitre strain equivalence hypothesis [27], we obtain(4)E1t=E11−Dt,where Dt is the damage variable, a function of the time. The other parameters are stated as above.Kachanov [28] investigated the creep damage and assumed that the rock element creep damage rate could be expressed as(5)dDtdt=Bσ1−Dtn,where B and n are the material constants of the rock element.Whent = 0, the rock element is undamaged, Dt=0. By integrating equation (5), we obtain(6)Dt=1−1−B1+nσnt1/n+1.Assume the critical failure time of creep damage istF and the corresponding damage variable Dt of the rock element is 1. Then, we obtain(7)tF=1B1+nσn.From equations (6)∼(7), we obtain(8)Dt=1−1−ttF1/n+1.Then, the constitutive equation of the damage body can be obtained:(9)ε=σ0E11−ttF−1/n+1.Finally, we obtain the creep constitutive equation of the J body model with damage(10)ε=σ0E11−Dt+E2−σ0E11−Dt⋅exp−E11−DtE2tηE11−Dt+E2+σ0E11−Dt.


## 3.2. Determination of the Parameters in the Proposed Model

The proposed creep damage constitutive model based on the J body model has five parameters namelyE1, E2, n, tF, and η, and how to determine them is also important. Now the often adopted methods are fitting with the experimental data, which are not only complicated but also cannot clearly explain their physical meaning. Therefore, it is very necessary to seek other methods to determine them. Now the determination of these parameters is discussed from the physical meaning of the proposed model.


### 3.2.1. The Critical Creep Failure TimetF

As stated above, the creep damage of the rock element is 1 at the critical creep failure timetF, and the corresponding creep strain and strain rate both tend to be infinite. As shown in Figure 1, the critical creep fracture time tF is easily obtained from the rock element creep experimental curve.


### 3.2.2. The Elastic Moduli ofE1 and E2

From equation (10), it can be seen that Dt=0 at the beginning time t = 0 s, and then under the stress σ0, the initial elastic strain ε0 of the proposed model is(11)ε0=σ0E1+E2.After the rock element goes into the steady creep stage, the creep curve exhibits a linear increase and no damage occurs in the rock element. Therefore, the elastic modulus of the elastic componentE1 is still E1. From the structure characteristic of the J body, it can be seen that the intersection point of the oppositely elongated line of the steady creep stage with the strain axial is the final strain ε1 of the elastic component E1 (shown in Figure 1), namely,(12)ε1=σ0E1.Combining equations (11) and (12), the elastic moduli E1 and E2 of these two elastic bodies can be solved.


### 3.2.3. Determination ofn

From the rock element creep deformation process and its mechanism, it can be seen that the viscoelastic deformation is finished after the rock element goes into the accelerating creep stage, and only the permanent deformation occurs with time going on. Then, from equation (9), the permanent strain εr of the sample is(13)εr=σ0E11−ttF−1/n+1.Assumem points (ti, εi) (i = 1, 2, …, m) are selected from the accelerating creep stage, and at the same time σ0 is known and tF and E1 have been solved; therefore, the calculation method of ni is(14)ni=ln1−ti/tFlnσ0/E1εri−1.The average value ofni is assumed to be the value of n:(15)n=1m∑i=1mni.


### 3.2.4. Determination ofη

The strain-time curve of the viscous component is linear which indicates that the viscosity is mainly reflected in the stable creep stage. Assumem points (ti, εi) (i = 1, 2, …, m) are selected from the stable creep stage, and there is(16)εi=σ0E11−Dti+E2−σ0E11−Dti⋅exp−E11−DtiE2tiηiE11−Dti+E2+σ0E11−Dti.So the corresponding viscosity coefficientηi (i = 1, 2, …, m) can be obtained. Finally, the average value of ηi is assumed to be the value of η(17)η=1m∑i=1mηi.


## 3.2.1. The Critical Creep Failure TimetF

As stated above, the creep damage of the rock element is 1 at the critical creep failure timetF, and the corresponding creep strain and strain rate both tend to be infinite. As shown in Figure 1, the critical creep fracture time tF is easily obtained from the rock element creep experimental curve.


## 3.2.2. The Elastic Moduli ofE1 and E2

From equation (10), it can be seen that Dt=0 at the beginning time t = 0 s, and then under the stress σ0, the initial elastic strain ε0 of the proposed model is(11)ε0=σ0E1+E2.After the rock element goes into the steady creep stage, the creep curve exhibits a linear increase and no damage occurs in the rock element. Therefore, the elastic modulus of the elastic componentE1 is still E1. From the structure characteristic of the J body, it can be seen that the intersection point of the oppositely elongated line of the steady creep stage with the strain axial is the final strain ε1 of the elastic component E1 (shown in Figure 1), namely,(12)ε1=σ0E1.Combining equations (11) and (12), the elastic moduli E1 and E2 of these two elastic bodies can be solved.


## 3.2.3. Determination ofn

From the rock element creep deformation process and its mechanism, it can be seen that the viscoelastic deformation is finished after the rock element goes into the accelerating creep stage, and only the permanent deformation occurs with time going on. Then, from equation (9), the permanent strain εr of the sample is(13)εr=σ0E11−ttF−1/n+1.Assumem points (ti, εi) (i = 1, 2, …, m) are selected from the accelerating creep stage, and at the same time σ0 is known and tF and E1 have been solved; therefore, the calculation method of ni is(14)ni=ln1−ti/tFlnσ0/E1εri−1.The average value ofni is assumed to be the value of n:(15)n=1m∑i=1mni.


## 3.2.4. Determination ofη

The strain-time curve of the viscous component is linear which indicates that the viscosity is mainly reflected in the stable creep stage. Assumem points (ti, εi) (i = 1, 2, …, m) are selected from the stable creep stage, and there is(16)εi=σ0E11−Dti+E2−σ0E11−Dti⋅exp−E11−DtiE2tiηiE11−Dti+E2+σ0E11−Dti.So the corresponding viscosity coefficientηi (i = 1, 2, …, m) can be obtained. Finally, the average value of ηi is assumed to be the value of η(17)η=1m∑i=1mηi.


## 4. Verification of the Proposed Model for the Rock Element

In order to validate the proposed model, the rock element creep experiment made by Kang et al. [10] is adopted. The rock is selected from Pingdingshan coal mine located in Henan Province, China. All the tests are carried out on the cylindrical specimens at room temperature. The target size of cylindrical specimens is 5 cm in diameter and 10 cm high according to the international standard for rock mechanics tests, and its average mass density is 1.59 g/cm3. All creep tests are conducted on an MTS815.02 rock servo-controlled equipment, and the other details of the experiment can be found in Reference [10]. When the axial stress is little, only the instant elastic deformation and transient creep deformation occur. With increase in the axial stress to a certain value, the creep damage will happen, and accordingly the accelerating creep stage occurs. Therefore, the verification of the proposed model can be classified into two categories namely the steady creep deformation model and the unsteady creep deformation model.


### 4.1. The Steady Creep Deformation Model

The steady creep deformation curve can be obtained when the axial stress is little; therefore, the typical strain-time curve under the axial stressσ0 = 3 MPa obtained by Kang et al. [10] is shown in Figure 3. For this condition, tF = ∞, and then from equation (8), we can obtain Dt=0. Therefore, equation (10) can be simplified as the classic J body model. Then, according to the proposed method, we can obtain the parameters of the rock element steady creep model are as follows: σ0 = 3.0 MPa, ε0 = 0.0011, E1 = 2143 MPa, E2 = 584 MPa, and η = 180 MPa·h. Figure 3 compares the strain-time curve tested by Kang et al. [10] and the one predicted by the proposed model. It can be seen that they fit with each other very well. It indicates that the proposed model can perfectly simulate the steady creep deformation of the rock element.Figure 3
Comparison of the predicted results and experimental data under the axial stressσ0 = 3 MPa [10].


### 4.2. The Unsteady Creep Deformation Model

The unsteady creep deformation curve can be obtained when the axial stress is large. The strain-time curve under the axial stressσ0 = 12 MPa obtained by Kang et al. [10] is shown in Figure 4, which has the accelerating creep stage. Then, according to the proposed method, we obtain the parameters of the rock element unsteady creep model as follows: σ0 = 12.0 MPa, ε0 = 0.00183, E1 = 3333 MPa, E2 = 3217 MPa, η = 120 MPa·h, tF = 1.001 h, and n = 2.21. Figure 4 compares the strain-time curve tested by Kang et al. [10] and the one predicted by the proposed model. It can be seen that they fit with each other. It indicates that the proposed model can perfectly simulate the unsteady creep deformation of the rock element. Meanwhile, we also give the predicted result by the J body model in Figure 4, which indicates that the original J body model cannot simulate the unsteady creep deformation of the rock element.Figure 4
Comparison of the predicted results and experimental data under the axial stressσ0 = 12 MPa [10].


## 4.1. The Steady Creep Deformation Model

The steady creep deformation curve can be obtained when the axial stress is little; therefore, the typical strain-time curve under the axial stressσ0 = 3 MPa obtained by Kang et al. [10] is shown in Figure 3. For this condition, tF = ∞, and then from equation (8), we can obtain Dt=0. Therefore, equation (10) can be simplified as the classic J body model. Then, according to the proposed method, we can obtain the parameters of the rock element steady creep model are as follows: σ0 = 3.0 MPa, ε0 = 0.0011, E1 = 2143 MPa, E2 = 584 MPa, and η = 180 MPa·h. Figure 3 compares the strain-time curve tested by Kang et al. [10] and the one predicted by the proposed model. It can be seen that they fit with each other very well. It indicates that the proposed model can perfectly simulate the steady creep deformation of the rock element.Figure 3
Comparison of the predicted results and experimental data under the axial stressσ0 = 3 MPa [10].


## 4.2. The Unsteady Creep Deformation Model

The unsteady creep deformation curve can be obtained when the axial stress is large. The strain-time curve under the axial stressσ0 = 12 MPa obtained by Kang et al. [10] is shown in Figure 4, which has the accelerating creep stage. Then, according to the proposed method, we obtain the parameters of the rock element unsteady creep model as follows: σ0 = 12.0 MPa, ε0 = 0.00183, E1 = 3333 MPa, E2 = 3217 MPa, η = 120 MPa·h, tF = 1.001 h, and n = 2.21. Figure 4 compares the strain-time curve tested by Kang et al. [10] and the one predicted by the proposed model. It can be seen that they fit with each other. It indicates that the proposed model can perfectly simulate the unsteady creep deformation of the rock element. Meanwhile, we also give the predicted result by the J body model in Figure 4, which indicates that the original J body model cannot simulate the unsteady creep deformation of the rock element.Figure 4
Comparison of the predicted results and experimental data under the axial stressσ0 = 12 MPa [10].


## 5. Establishment of the Creep Damage Model for the Rock Mass with Nonpersistent Joints

### 5.1. The Damage Caused by One Nonpersistent Joint to the Rock Mass

Under uniaxial compression shown in Figure5, Li et al. [29] obtained the damage variable Dj along the loading direction caused by one nonpersistent joint to the rock mass according to fracture mechanics for a plane stress issue:(18)Dj=1−11+2/V1/σ2∫0AKΙ2+KΙΙ2dA,where V is the volume of the rock mass, KI and KII are the first and second stress intensity factors (SIFs) at the joint tip, respectively, A is the joint area, for one joint, A = Ba (unilateral joint) or 2Ba (central joint). B is the joint depth and 2a is the joint length, as shown in Figure 5.Figure 5
Model of a rock mass with one nonpersistent joint: joint length 2a, joint depth b, joint dip angle α, and wing joint initiation angle θ.Then, based on this method, Liu et al. [30] have deduced the calculation formula of the damage variable caused by one nonpersistent joint, which can consider the effect of three kinds of joint parameters, namely, the joint geometry parameter (such as the joint length and dip angle), strength parameter (such as the joint friction angle ϕ), and deformation parameter (such as the joint normal stiffness kn and shear stiffness ks) on the macroscopic damage variable. And we also adopt this method and do not state it again.


### 5.2. A Creep Damage Constitutive Model for the Jointed Rock Mass by Coupling the Mesoscopic and Macroscopic Flaws

As stated above, the rock mass contains both two different scales of flaws which will affect the rock mass mechanical behavior at the same time. So Liu et al. [31] proposed the coupled damage caused by the macroscopic and mesoscopic flaws. If assume the damage variable along the loading direction induced by the macroscopic and mesoscopic flaws is Dj and D, respectively, their coupled damage variable Dc is(19)Dc=1−1−Dj1−D1−DjD.It is known from the damage theory that the different scales of flaws in a rock mass are corresponding to the different damage variables. So, substituting the damage variableD in equation (10) with the coupled damage variable Dc in equation (19), the uniaxial creep damage constitutive model for a rock mass with nonpersistent joints can be obtained:(20)εj=σ0E1χ+E2−σ0E1χ⋅exp−E1χE2tηiE1χ+E2+σ0E1χ,where χ=1−Dj1−Dt/1−DjDt, and Dt is solved with equation (8); εj is the rock mass strain.


## 5.1. The Damage Caused by One Nonpersistent Joint to the Rock Mass

Under uniaxial compression shown in Figure5, Li et al. [29] obtained the damage variable Dj along the loading direction caused by one nonpersistent joint to the rock mass according to fracture mechanics for a plane stress issue:(18)Dj=1−11+2/V1/σ2∫0AKΙ2+KΙΙ2dA,where V is the volume of the rock mass, KI and KII are the first and second stress intensity factors (SIFs) at the joint tip, respectively, A is the joint area, for one joint, A = Ba (unilateral joint) or 2Ba (central joint). B is the joint depth and 2a is the joint length, as shown in Figure 5.Figure 5
Model of a rock mass with one nonpersistent joint: joint length 2a, joint depth b, joint dip angle α, and wing joint initiation angle θ.Then, based on this method, Liu et al. [30] have deduced the calculation formula of the damage variable caused by one nonpersistent joint, which can consider the effect of three kinds of joint parameters, namely, the joint geometry parameter (such as the joint length and dip angle), strength parameter (such as the joint friction angle ϕ), and deformation parameter (such as the joint normal stiffness kn and shear stiffness ks) on the macroscopic damage variable. And we also adopt this method and do not state it again.


## 5.2. A Creep Damage Constitutive Model for the Jointed Rock Mass by Coupling the Mesoscopic and Macroscopic Flaws

As stated above, the rock mass contains both two different scales of flaws which will affect the rock mass mechanical behavior at the same time. So Liu et al. [31] proposed the coupled damage caused by the macroscopic and mesoscopic flaws. If assume the damage variable along the loading direction induced by the macroscopic and mesoscopic flaws is Dj and D, respectively, their coupled damage variable Dc is(19)Dc=1−1−Dj1−D1−DjD.It is known from the damage theory that the different scales of flaws in a rock mass are corresponding to the different damage variables. So, substituting the damage variableD in equation (10) with the coupled damage variable Dc in equation (19), the uniaxial creep damage constitutive model for a rock mass with nonpersistent joints can be obtained:(20)εj=σ0E1χ+E2−σ0E1χ⋅exp−E1χE2tηiE1χ+E2+σ0E1χ,where χ=1−Dj1−Dt/1−DjDt, and Dt is solved with equation (8); εj is the rock mass strain.


## 6. Analysis of the Calculation Examples

### 6.1. The Creep Mechanical Behavior of a Rock Mass with One Nonpersistent Joint

In order to illustrate the effect of the nonpersistent joint on the rock mass creep mechanical behavior, the rock element creep experiment made by Kang et al. [10] is also adopted here. According to the experiment made by Kang et al. [10], the elastic modulus of the rock element under the axial stress of 3 MPa and 12 MPa are 2727 MPa and 6550 MPa, respectively. Meanwhile, assume the Poisson’s ratio of the rock element is 0.3 which is assumed to be not affected by the axial stress. Assume there is one nonpersistent joint in the sample (Figure 6), the damage variables along the loading direction caused by the nonpersistent joint to the rock mass under the axial stress 3 MPa and 12 MPa are 0.0054 and 0.013, respectively, which can be calculated with equation (18).Figure 6
The calculation model. The joint parameters are as follows: 2a = 2 cm, B = 1 cm, α = 45°, kn = 2000 MPa/cm, ks = 1000 MPa/cm, and ϕ = 10°, respectively.Then, the predictions of the creep strain-time curves for the rock mass sample under the axial stressσ0 = 3 MPa and 12 MPa are shown in Figures 6 and 7. The following results can be obtained. For the axial stress σ0 = 3 MPa, only the steady creep deformation occurs in the rock mass. The only difference is that the strain of the rock mass with one nonpersisent joint is gently larger than that of the rock element. This is because the existence of the joint causes the damage to the rock mass, which reduces its strength and increases its deformability. However, when the axial stress increases to 12 MPa, the unsteady creep deformation will occur in both the rock element and the rock mass. Meanwhile, the strain-time curve of the rock mass becomes much complicated comparing with that of the rock element. At the stage of transient creep, the strain of the rock mass is a little larger than that of the rock element because of the joint. However, with time going on, the strain of the rock mass gradually becomes less than that of the rock element at the steady creep stage, and finally the strain of the rock mass gradually becomes much larger than that of the rock element at the accelerating creep stage. It indicates that there exists the complicated interaction between the microcracks and the joint.Figure 7
Comparison of the creep strain-time curves for the rock element and the rock mass with one nonpersistent joint at different axial stressσ0.
(a)
(b)


### 6.2. The Effect of the Joint Length on the Rock Mass Creep Mechanical Behavior

Here the calculation model in Figure6 is also adopted to discuss the effect of the joint length on the sample creep mechanical behavior. Assume the axial stress is also 12 MPa. It can be seen from Figure 8 that when the joint length increases from 2 cm to 3 cm and 4 cm, respectively, the creep strain will increase. Here, the creep strain at time t = 0.81 h is taken for an example, and the creep strain will increase from 0.00536 to 0.00565 and 0.00597, respectively. It indicates that increase in the joint length will lead to the increase in the rock mass damage, which accordingly increases the sample creep deformability.Figure 8
Effect of the joint length on the creep stress-strain curves of the rock mass with one nonpersistent joint at axial stressσ0 = 12 MPa.


### 6.3. The Effect of the Joint Friction Angle on the Rock Mass Creep Mechanical Behavior

Here, the calculation model in Figure6 is also adopted to discuss the effect of the joint friction angle on the sample creep mechanical behavior. Assume the axial stress is also 12 MPa. It can be seen from Figure 9 that when the joint friction angle increases from 10° to 20° and 30°, respectively, the creep strain will decrease. Here, the creep strain at time t = 0.81 h is taken for an example, and the creep strain will increase from 0.00536 to 0.00525 and 0.00515, respectively. It indicates that increase in the joint friction angle will lead to the decrease in the rock mass damage, which accordingly reduces the sample creep deformability.Figure 9
Effect of the joint friction angle on the creep stress-strain curves of the rock mass with one nonpersistent joint at axial stressσ0 = 12 MPa.


### 6.4. The Effect of the Joint Normal and Shear Stiffness on the Rock Mass Creep Mechanical Behavior

Here, the calculation model in Figure6 is also adopted to discuss the effect of the joint deformation parameter, namely, the joint normal and shear stiffness on the sample creep mechanical behavior. Assume the axial stress is 12 MPa. The following conclusions can be obtained from Figure 10. Firstly, when the joint normal stiffness kn increases from 200 MPa/cm to 2000 MPa/cm and 20000 MPa/cm, respectively, the creep strain will increase. Here, the creep strain at time t = 0.81 h is taken for an example, and the creep strain will increase from 0.00531 to 0.00536 and 0.00545, respectively. It indicates that increase in the joint normal stiffness will lead to the increase in the sample creep strain. This is because the resolved normal stress on the joint face will decrease with increase in the joint normal stiffness, and therefore the friction to resist the slide of the rock block along the joint face will decrease which will easily cause the rock mass to shear failure along the joint face. Therefore, it leads to the increase in damage and the rock mass creep strain. However, when the joint shear stiffness increases from 100 MPa/cm to 1000 MPa/cm and 10000 MPa/cm, the sample creep strain at time t = 0.81 h decreases from 0.00553 to 0.00536 and 0.00505, respectively. It indicates that increase in the joint shear stiffness will lead to the decrease in the sample creep strain. In sum, the joint normal stiffness and shear stiffness has some effect on the sample creep strain. But this result does not agree with that obtained by Prudencio and Van Sint Jan [32] who found that the jointed rock mass flexibility would decrease with increase in the joint normal stiffness and shear stiffness. However, the results in Figures 10(a) and 3(b) show that the jointed rock mass flexibility increases and decreases with increase in the joint normal stiffness and shear stiffness, respectively. Through the thorough research, we assume that the arbitrary selection of the joint normal stiffness and shear stiffness is the primary reason leading to this result. Bandis et al. [33] found that joint shear stiffness was much lower than its normal stiffness. However, through the study on the variation of the Kn/Ks ratio with the normal stress, they found that the Kn/Ks ratio keeps constant for a given normal stress. Therefore, the effect of the joint normal stiffness and shear stiffness on the rock mass creep deformation behavior should be studied by varying them in proportion. Here, the initial joint normal stiffness and shear stiffness is assumed to be 200 MPa/cm and 100 MPa/cm, respectively, as before, and then we increase them by 10 and 100 times, respectively. It can be seen from Figure 10(c) that the rock mass creep strain at time t = 0.81 h decreases from 0.00548 to 0.00536 and 0.00507, respectively. It shows that the rock mass creep strain decreases with increase in the joint normal and shear stiffness and fits with the research conclusion obtained by Prudencio and Van Sint Jan [32]. It indicates that the effect of the joint shear stiffness on the rock mass creep strain is much more than that of the joint normal stiffness, and when they increase in proportion, the rock mass creep strain decreases.Figure 10
Effect of the joint normal and shear stiffness on the creep strain-time curves of the rock mass with one nonpersistent joint at axial stressσ0 = 12 MPa. (a) The joint normal stiffness is different; (b) the joint shear stiffness is different; (c) the joint normal and shear stiffness varies at the same Kn/Ks ratio.
(a)
(b)
(c)


## 6.1. The Creep Mechanical Behavior of a Rock Mass with One Nonpersistent Joint

In order to illustrate the effect of the nonpersistent joint on the rock mass creep mechanical behavior, the rock element creep experiment made by Kang et al. [10] is also adopted here. According to the experiment made by Kang et al. [10], the elastic modulus of the rock element under the axial stress of 3 MPa and 12 MPa are 2727 MPa and 6550 MPa, respectively. Meanwhile, assume the Poisson’s ratio of the rock element is 0.3 which is assumed to be not affected by the axial stress. Assume there is one nonpersistent joint in the sample (Figure 6), the damage variables along the loading direction caused by the nonpersistent joint to the rock mass under the axial stress 3 MPa and 12 MPa are 0.0054 and 0.013, respectively, which can be calculated with equation (18).Figure 6
The calculation model. The joint parameters are as follows: 2a = 2 cm, B = 1 cm, α = 45°, kn = 2000 MPa/cm, ks = 1000 MPa/cm, and ϕ = 10°, respectively.Then, the predictions of the creep strain-time curves for the rock mass sample under the axial stressσ0 = 3 MPa and 12 MPa are shown in Figures 6 and 7. The following results can be obtained. For the axial stress σ0 = 3 MPa, only the steady creep deformation occurs in the rock mass. The only difference is that the strain of the rock mass with one nonpersisent joint is gently larger than that of the rock element. This is because the existence of the joint causes the damage to the rock mass, which reduces its strength and increases its deformability. However, when the axial stress increases to 12 MPa, the unsteady creep deformation will occur in both the rock element and the rock mass. Meanwhile, the strain-time curve of the rock mass becomes much complicated comparing with that of the rock element. At the stage of transient creep, the strain of the rock mass is a little larger than that of the rock element because of the joint. However, with time going on, the strain of the rock mass gradually becomes less than that of the rock element at the steady creep stage, and finally the strain of the rock mass gradually becomes much larger than that of the rock element at the accelerating creep stage. It indicates that there exists the complicated interaction between the microcracks and the joint.Figure 7
Comparison of the creep strain-time curves for the rock element and the rock mass with one nonpersistent joint at different axial stressσ0.
(a)
(b)


## 6.2. The Effect of the Joint Length on the Rock Mass Creep Mechanical Behavior

Here the calculation model in Figure6 is also adopted to discuss the effect of the joint length on the sample creep mechanical behavior. Assume the axial stress is also 12 MPa. It can be seen from Figure 8 that when the joint length increases from 2 cm to 3 cm and 4 cm, respectively, the creep strain will increase. Here, the creep strain at time t = 0.81 h is taken for an example, and the creep strain will increase from 0.00536 to 0.00565 and 0.00597, respectively. It indicates that increase in the joint length will lead to the increase in the rock mass damage, which accordingly increases the sample creep deformability.Figure 8
Effect of the joint length on the creep stress-strain curves of the rock mass with one nonpersistent joint at axial stressσ0 = 12 MPa.


## 6.3. The Effect of the Joint Friction Angle on the Rock Mass Creep Mechanical Behavior

Here, the calculation model in Figure6 is also adopted to discuss the effect of the joint friction angle on the sample creep mechanical behavior. Assume the axial stress is also 12 MPa. It can be seen from Figure 9 that when the joint friction angle increases from 10° to 20° and 30°, respectively, the creep strain will decrease. Here, the creep strain at time t = 0.81 h is taken for an example, and the creep strain will increase from 0.00536 to 0.00525 and 0.00515, respectively. It indicates that increase in the joint friction angle will lead to the decrease in the rock mass damage, which accordingly reduces the sample creep deformability.Figure 9
Effect of the joint friction angle on the creep stress-strain curves of the rock mass with one nonpersistent joint at axial stressσ0 = 12 MPa.


## 6.4. The Effect of the Joint Normal and Shear Stiffness on the Rock Mass Creep Mechanical Behavior

Here, the calculation model in Figure6 is also adopted to discuss the effect of the joint deformation parameter, namely, the joint normal and shear stiffness on the sample creep mechanical behavior. Assume the axial stress is 12 MPa. The following conclusions can be obtained from Figure 10. Firstly, when the joint normal stiffness kn increases from 200 MPa/cm to 2000 MPa/cm and 20000 MPa/cm, respectively, the creep strain will increase. Here, the creep strain at time t = 0.81 h is taken for an example, and the creep strain will increase from 0.00531 to 0.00536 and 0.00545, respectively. It indicates that increase in the joint normal stiffness will lead to the increase in the sample creep strain. This is because the resolved normal stress on the joint face will decrease with increase in the joint normal stiffness, and therefore the friction to resist the slide of the rock block along the joint face will decrease which will easily cause the rock mass to shear failure along the joint face. Therefore, it leads to the increase in damage and the rock mass creep strain. However, when the joint shear stiffness increases from 100 MPa/cm to 1000 MPa/cm and 10000 MPa/cm, the sample creep strain at time t = 0.81 h decreases from 0.00553 to 0.00536 and 0.00505, respectively. It indicates that increase in the joint shear stiffness will lead to the decrease in the sample creep strain. In sum, the joint normal stiffness and shear stiffness has some effect on the sample creep strain. But this result does not agree with that obtained by Prudencio and Van Sint Jan [32] who found that the jointed rock mass flexibility would decrease with increase in the joint normal stiffness and shear stiffness. However, the results in Figures 10(a) and 3(b) show that the jointed rock mass flexibility increases and decreases with increase in the joint normal stiffness and shear stiffness, respectively. Through the thorough research, we assume that the arbitrary selection of the joint normal stiffness and shear stiffness is the primary reason leading to this result. Bandis et al. [33] found that joint shear stiffness was much lower than its normal stiffness. However, through the study on the variation of the Kn/Ks ratio with the normal stress, they found that the Kn/Ks ratio keeps constant for a given normal stress. Therefore, the effect of the joint normal stiffness and shear stiffness on the rock mass creep deformation behavior should be studied by varying them in proportion. Here, the initial joint normal stiffness and shear stiffness is assumed to be 200 MPa/cm and 100 MPa/cm, respectively, as before, and then we increase them by 10 and 100 times, respectively. It can be seen from Figure 10(c) that the rock mass creep strain at time t = 0.81 h decreases from 0.00548 to 0.00536 and 0.00507, respectively. It shows that the rock mass creep strain decreases with increase in the joint normal and shear stiffness and fits with the research conclusion obtained by Prudencio and Van Sint Jan [32]. It indicates that the effect of the joint shear stiffness on the rock mass creep strain is much more than that of the joint normal stiffness, and when they increase in proportion, the rock mass creep strain decreases.Figure 10
Effect of the joint normal and shear stiffness on the creep strain-time curves of the rock mass with one nonpersistent joint at axial stressσ0 = 12 MPa. (a) The joint normal stiffness is different; (b) the joint shear stiffness is different; (c) the joint normal and shear stiffness varies at the same Kn/Ks ratio.
(a)
(b)
(c)


## 7. Conclusions

By introducing the damage component into the J body model, this work proposes a new creep damage model for the rock element which can perfectly reflect the accelerating creep stage.Then, a new creep damage constitutive model for the rock mass with one nonpersistent joint is proposed, which can consider the effect of the macroscopic and mesoscopic flaws on the rock mass creep mechanical behavior.The calculation examples indicate that the proposed creep damage model for the jointed rock mass can reasonably reflect the effect of joint geometry, strength, and deformation parameters on the rock mass creep mechanical behavior at the same time. Overall, the proposed model provides a new way to investigate the creep mechanical behavior of the rock mass with nonpersistent joints.



---

*Source: 4361458-2019-11-11.xml*