Analysis of a microcrack model and constitutive equations for time-dependent dilatancy of rocks Free

Summary

1 Introduction

Dilatancy of rocks plays an important role for tectonic movements inside the Earth as well as in many aspects of engineering (Brace 1966; Brace 1978; Tan Tjong Kie & Kang Wen-fa 1983). Many laboratory studies have been stimulated by the hypothesis that dilatancy is believed to play an important role in the mechanism of shallow earthquakes (Nur 1972; Scholz 1973). Basic studies on dilatancy are very helpful in gaining a better understanding of the complexity of fracturing in rock mechanics. The formation of cracks and loosening substances around a tunnel, which may lead to gradual instability, is also attributed to dilatancy.

To obtain a better understanding of the tectonics of the crust and the genesis of earthquakes, a profound knowledge concerning the brittle behaviour of rocks at various stresses and temperatures is desirable (Tan Tjong Kie 1989; Ngwenya 2001). A basic property of the rheology of rock is creep deformation. Under normal stresses, the typical creep curve may be subdivided into three parts: (1) transient creep in which the strain rate decreases with time; (2) steady creep in which the strain rate is constant; and (3) tertiary creep in which the strain rate increases progressively until failure. However, the microscopic mechanism of creep in rocks has not been studied comprehensively.

At low temperatures, it is observed that the volume strain of rocks under pressure increases with time due to the extension of microcracks in rocks (Kranz 1979). During experiments, it is also found that a crack often shows a significant extension rate at a value of stress intensity factor lower than the critical stress intensity factor K_c. Such a phenomenon may be attributed to several mechanisms including a dependence upon the chemical environment, among which stress corrosion is mentioned most frequently (Martin 1972; Anderson & Grew 1976; Atkinson 1984; Atkinson & Meredith 1987). Stress corrosion involves the weakening of strained bonds at crack tips with the chemical action of, for example, water, thereby facilitating crack growth. Other contributing mechanisms include dissolution, diffusion, ion exchange and microplasticity.

The concept of different scale media may be used to explain that a rock can exhibit different mechanical properties. For example, brittle rocks exhibit a creep-dilating behaviour under triaxial stresses, but when rock material around microcracks is viscoelastic there is no swelling. Such dilatancy is the result of an extension of microcracks, namely when dissipation of energy leads to subcritical growth of cracks over time. Therefore, the rock specimen as a whole has the properties of both viscoplasticity and time-dependent dilatancy, while the skeleton of the rock is viscoelastic.

In this paper, rheological fracture mechanics is used to analyse the process of crack growth for rocks during creep to explain the mechanism of creep dilatancy for rocks. A differential equation that describes the relationship between the deformation velocity of the rock and time is obtained. By solving this equation, a relationship between volumetric swelling and time is derived. It agrees well with the character of typical experimental creep curves.

Moreover, equations for damage constitution and damage evolution are presented. In these equations, the dilatancy volume strain is taken as the damage variable. The damage evolution equation can be determined from the creep curves obtained via an analysis of the course of microcrack extension. By this procedure we intend to find a way to establish and describe the macroscopic mechanical behaviour of rocks through an analysis of their microscopic structure.

2 Model of Microcracks

Most rocks contain inhomogeneities such as small cavities or cracks, particles with poor cohesion, as well as phases with different modulus or strengths from those of the matrix. Generally, the distribution of microcracks in rocks is very complicated. Therefore, it is very difficult to analyse the growth of various cracks in a complete way. However, based on experiments (Topponnier & Brace 1976; Waza 1980; Wong 2001), the sharp inclined crack shown in Fig. 1 is a model of a typical microcrack in rocks; two wings of the crack lie parallel to the maximum principal compressive stress. As for the mechanism of inelastic dilatancy in rocks, most studies focus on the sliding crack model first proposed by Brace (1966). The theoretical analysis by Nemat-Nasser & Horii (1982) and the numerical simulation by Tang & Kou (1998) show that under triaxal compressive stresses, the inclined crack in solids will extend in the direction at an angle of 0.39π from the crack inclination, and continue to extend parallel to the maximum compressive stress direction σ₁.

Figure 1

A wing crack.

Open in new tabDownload slide

Based on experimental observations and a theoretical analysis, the author introduces an ideal model for microcracks, which is quite similar to the model proposed by Ashby & Sammis (1990). In this model, an array of cracks with the same shape and initial size are distributed evenly in rocks (Fig. 2). Under triaxial compressive stresses, the crack will extend parallel to the direction of the maximum compressive stress. In accordance with Ashby & Sammis (1990), the stress intensity factor at tips of a wing crack is analysed as follows.

Figure 2

Model of extending microcracks.

Open in new tabDownload slide

The remote field stress σ₁, σ₃ (positive for dilatation and negative for compression) will create a shear stress τ and a normal stress σ on the surface of initial cracks. The sliding of cracks (resisted by the frictional force, with friction coefficient μ), opens the mouth of each wing crack. The wedging can be deemed to be a result of the force F₃ that is parallel to X₃ and acts at the midpoint of the crack. Stresses τ and σ are given by

1

where ψ is the angle between the initial crack and the main stress σ₁.

F₃ is the component of the sliding force acting parallel to σ₃,

2

It creates a stress intensity factor tending to open the crack (Tada 1985).

Also there is an average internal stress, σⁱ₃, acting on the wing cracks that describes the interaction among cracks:

3

where S is the distance between the central points of two wing cracks and 2a is the length of the crack. If the number of cracks in unit area is N, the relation between S and N is

4

Therefore, the stress intensity factor at tips of wing cracks is given by

5

We then choose β so that the first term in eq. (5) becomes equal to that for the inclined crack when L= 0. The remote field compressive stress σ₃ creates an additional stress intensity factor tending to close the crack. The first term in the brackets of eq. (5) describes the interaction between the wedging and cracks, and the second term describes the closing effect due to the lateral compressive stress.

3 Analysis of Microcrack Extension for Creep Dilatancy

Griffith (1924) has indicated that if a crack grows to length L, the release rate of elastic energy is greater than or equal to the increasing rate of surface energy, so the crack will be in an unstable state and must extend. If we suppose that the medium around the crack is viscoelastic, due to the stress singularity at the tips of cracks, a higher dissipation of energy will occur, which will cause a subcritical extension of microcracks. In this paper the law of conservation is applied directly as the failure criterion. This method is a generalization of the classical Griffith approach for cracks in elastic materials.

Under the isothermal condition, the appropriate statement of the global conservation of energy for extension of cracks in viscoelastic media (Christensen 1982) is given by

6

Herein U is the elastic energy in the infinite plate containing one crack, D_p is the dissipation energy due to the rheology of the medium and S_e is the surface energy. Eq. (6) may also be written as

7

The expression for the energy change rate (Jaeger & Cook 1976) may be used for the model mentioned above. When neglecting K_Π,

8

where ν is the Poisson ratio and E represents the Young’s modulus.

To make the analysis easier, we suppose the extension of cracks will take place step by step. Accordingly, when the stresses act on rocks for the time period T(L= 0), the crack begins growing to a length L and then stops. At another time period T(L) later on, the crack continues to grow to a new length L and stops again. Furthermore, there exists a critical length L_c. When L is larger then L_c, T(L) is equal to zero, which means crack growth will accelerate until failure. T(L) is the time period from the stop of the last crack extension up to the start of the next crack extension.

Considering the actual process of crack extension described above, the stress around tips of cracks will be constant in duration of crack growth stopping temporarily, in accordance with the correspondence principle. The stress formula for a type I crack is

9

where r and θ are polar coordinates, and the origin is at the tip of the wing crack (which will move with the extension of the crack). σ_ij is the stress tensor and K_I is the stress intensity factor. f_ij(θ) is a function depending on loading type.

The viscoelastic constitutive equation for creep may be written as follows:

10

where D(t) and K(t) are the deviatoric strain creep compliance and the volume strain creep compliance, respectively. σ_ij and ɛ_ij are the stress tensor and strain tensor, respectively, in a coordinate system that moves with the crack extension. The origin of the coordinate system is again at the tips of cracks. While S_ij is the deviatoric stress tensor, p is the hydrostatic stress, e_ij is the deviatoric strain tensor and e_v is the volumetric strain.

The dissipative energy Q at some point of the medium may be calculated as

11

where time T counts from the halting of crack growth. It is assumed that the dissipation of energy is caused by the stress singularity at tips of the cracks.

Taking a circle around the tip of a crack, the total dissipative energy may be calculated as follows:

12

13

14

where negative signs indicate loss of energy.

The surface energy is

15

with η representing the surface energy in unit area.

Combining eqs (8), (13), (15) and (7), the following formula can be derived:

16

Eq. (16) must be satisfied by T(L).

If the viscoelastic creep compliance is taken as

17

the following equation can be found:

18

In accordance with the definition of T(L), the relation between T(L) and the extending velocity of cracks can be described as

19

Note that the time t herein is a scale for the total process of crack extension, and that it is different from T.

If the crack is assumed to be planar, the relation between the dilatancy volume strain and the extending length of cracks is

20

where α_c denotes the width of cracks and N is the number of microcracks in a planar unit area. ɛ_vp is the dilatancy volume strain.

4 Computing Analysis

Solving the differential eq. (19) using the Runge–Kutta algorithm, the relation between the length L and time t, and the relation between the creep volume strain and time may be obtained. Plausible values of parameters required for computational purposes are taken as a particular example. Fig. 3 displays a theoretical curve that is in line with a typical creep curve. The curve also includes three stages: primary, secondary and tertiary creep. Referring to eq. (5), it is obvious that in the primary creep stage, the stress intensity factor at the tips of cracks decreases with time due to the confining pressure; in the secondary creep stage, the stress intensity factor at the tips of cracks remains approximately constant with time due to the interaction of cracks and the confining pressure; while in the tertiary creep stage, the stress intensity factor at the tips of cracks increases quickly due to the interaction of cracks. These results show that the model established in this paper is an effective approach to describing the mechanism of creep dilatancy in rocks.

Figure 3

Creep volumetric strain curve computed by the model. The parameters for the computation are as follows: σ₁= 120 MPa, σ₃= 35 MPa, a= 0.2 cm, μ= 0.1, E= 50 000 MPa, ν= 0.25, β= 0.1, η= 1000 MPa, a₀= 0.05/h, β_s= 1 (MPa)⁻¹, N= 3.3, α_c= 0.05 cm, ψ= 45°.

Open in new tabDownload slide

The influences upon the crack extension velocity of the initial length of a crack, the inclination of a crack and the field stresses, are as follows.

4.1 Initial length of crack

Apparently, the crack extension velocity increases with the initial length of a crack. Fig. 4 shows a group of curves obtained with initial lengths of cracks ranging from 0.23 to 0.26 cm. Other parameters were unchanged (same as in Fig. 3). Fig. 4 reveals that the time of failure for rocks with shorter cracks is longer than that for rocks with longer cracks. The strain rate increases with increasing initial length of pre-existing cracks.

Figure 4

Creep volumetric strain curves computed for different lengths of crack.

Open in new tabDownload slide

4.2 Inclination of crack

Fig. 5 shows a group of curves obtained with crack inclinations changing from π/3 to π/5 but with other parameters as in Fig. 3 during creep. Obviously, the strain rate increases when the inclination of cracks changes from π/3 to π/4, while it decreases when the inclination of cracks changes from π/4 to π/5. So, it is interesting that the extending velocity of cracks does not depend monotonically on the inclination of cracks. This also means that the flat or steep inclination of cracks will exhibit a longer failure time, or say the cracks will extend with more difficulty.

Figure 5

Creep volumetric strain curves computed for different inclinations of crack.

Open in new tabDownload slide

4.3 Confining field stress

The influence of confining field stresses upon the extension velocity of cracks falls in line with the theory and results we have expected. With increasing confining pressure the crack will extend with more difficulty and will require a longer failure time. This is shown in Fig. 6.

Figure 6

Creep volumetric strain curves computed for different confining pressures.

Open in new tabDownload slide

5 Constitutive Equation and Damage Evolution Equation

Another major objective of this paper is to establish and describe macroscopic mechanical behaviour through an analysis of the microscopic structure. Rice’s internal variable theory including a micro-to-macro transition has been employed to formulate a micromechanical, 2-D damage model of brittle deformation in compression (Basista & Gross 1998). Yet, our analysis method will be different from that of Basista & Gross. We introduce the associated constitutive equation suggested by Tan Tjong Kie & Kang Wen-fa (1983), namely

21

22

with

Herein σ_ij is the stress tensor, S_ij is the deviatoric stress tensor, p is the hydrostatic stress, e_ij is the deviatoric strain tensor, e_v is the volumetric strain, δ_ij is the Kronecker tensor, E, G and ν are the elastic parameters: Young’s modulus, shear modulus and Poisson’s ratio; D, n and f are material constants, where f is the threshold value for dilatancy in shear, and F(σ_oct/f) is a function describing the dilatancy strain.

Eqs (21) and (22) include two parts: elasticity and plasticity. The first two terms in eq. (21) and the first term in eq. (22) represent the elastic strain characterized by the elastic parameters E, G and ν; the last terms in eqs (21) and (22) show the dilatant part. The parameter D characterizes the rock behaviour in the dilatancy regime and is called the dilatancy parameter. It can be seen from the above expressions that the dilatant strain is limited by the factor (σ_oct/f) and by the function F(σ_oct/f).

If the dilatant volume strain is taken as the damage variable ω, it indicates the internal irrecoverable state for rocks. It is easily found that the second term in eq. (22) may be replaced by a damage variable ω. Compared with the third term in eq. (21), the following damage constitutive equation may be obtained from eqs (21) and (22):

23

where C is the reciprocal of D. Apparently, the elasticity property of rock here is supposed without any effect of dilatancy, and the plasticity property of rocks is relative to the dilatant strain or the damage variable. The parameter C also characterizes the plastic deformation of the rock shearing mechanism.

In addition to the above constitutive equations, the evolution equation for damage is required to be set up. The following is the suggested evolution equation:

24

with H_i being a positive material constant. The first term in the function H(ω) represents an exponentially descending curve corresponding to the strain softening portion of the stress–strain curve; the second term indicates the exponentially ascending curve corresponding to the strain hardening portion prior to the peak strength; the third term reveals the yield strength in the residual state. When H_i is selected and satisfies the conditions H₄ > H₂ and H₃H₄ > H₁H₂, the second term will lose the strain hardening features more rapidly than the first item. H(ω) can be obtained by fitting the creep curve of volumetric dilatancy strain with a non-linear least-squares method. Eq. (23) is meaningful when σ_oct > σ_y. When ω= 0, the formula σ_y=H₁−H₃+H₅, is the threshold value at the starting point of dilatancy. H₁, H₂, H₃, H₄ and H₅ have different values under different confining pressures.

Fig. 7 shows the strain hardening–strain softening function H(ω), which is obtained from the creep curve shown in Fig. 3. Solving eqs (13) and (24) under the condition of constant strain rate 2.0 × 10⁻⁴ s⁻¹, the stress–strain relation for rocks can be derived numerically. The computational steps are shown below:

25

26

27

28

where k is the number of iterations. As the parameter C in eq. (26) represents the shearing behaviour related to the dilatancy of rocks, it cannot be derived from current theories. However, some experiments indicate that C is close to 0.2 (Tan Tjong Kie 1989). Therefore, the value C= 0.2 is applied in the calculations. The computed stress–strain relationship is shown in Fig. 8. The characteristics of the curves during the deformation process are well indicated.

Figure 7

Hardening–softening H(ω) for the same parameters as Fig. 3.

Open in new tabDownload slide

Figure 8

The predicated behaviour of stress–strain under a constant strain rate of 2 × 10⁻⁴ s⁻¹.

Open in new tabDownload slide

6 Conclusion

An approximate physical model for time-dependent microcrack evolution in rocks under triaxial compressive stresses has been developed and analysed. The model is based on the subcritical growth of an array of small and inclined wing cracks as well as their interactions. The mechanism of creep dilatancy in rocks is well described through this model. The general agreement of the computational results with the typical creep curve for rocks indicates that this model is plausible. The analysis further indicates that the concentration of stresses at tips of microcracks, interaction of microcracks and energy dissipation are three important factors for crack extension in creeping rocks. The iterative analysis method for crack extension constitutes a distinguishing feature of this paper.

However, our final objective is to derive a macroscopic model for the mechanical behaviour of rocks through the analysis of microscopic structures. The work in this paper is only an initial attempt in this respect. It is concluded that the description of the mechanical property of rocks using the damaging mechanical method is a possible way of establishing a bridge between the macroscopic behaviour and the microscopic structure.

Acknowledgments

We thank Professors Bai Wuming and Yang Wenxiang for their careful reviews and constructive suggestions with regard to this paper. Thanks should also be given to the National Foundation Committee for Natural Science of China, which funded this study. The project code is 10032040.

References