# Rock mass plasticity

Boudinaged quartz vein(with strain fringe) showing sinistral shear sense, Starlight Pit, Fortnum Gold Mine, Western Australia

Plasticity theory for rocks is concerned with the response of rocks to loads beyond the elastic limit. Historically, conventional wisdom has it that rock is brittle and fails by fracture while plasticity is identified with ductile materials. In field scale rock masses, structural discontinuities exist in the rock indicating that failure has taken place. Since the rock has not fallen apart, contrary to expectation of brittle behavior, clearly elasticity theory is not the last work.[1]

Theoretically, the concept of rock plasticity is based on soil plasticity which is different from metal plasticity. In metal plasticity, for example in steel, the size of a dislocation is sub-grain size while for soil it is the relative movement of microscopic grains. The theory of soil plasticity was developed in the 1960s at Rice University to provide for inelastic effects not observed in metals. Typical behaviors observed in rocks include strain softening, perfect plasticity, and work hardening.

Application of continuum theory is possible in jointed rocks because of the continuity of tractions across joints even through displacements may be discontinuous. The difference between an aggregate with joints and a continuous solid is in the type of constitutive law and the values of constitutive parameters.

## Experimental evidence

Experiments are usually carried out with the intention of characterizing the mechanical behavior of rock in terms of rock strength. The strength is the limit to elastic behavior and delineates the regions where plasticity theory is applicable. Laboratory tests for characterizing rock plasticity fall into four overlapping categories: confining pressure tests, pore pressure or effective stress tests, temperature-dependent tests, and strain rate-dependent tests. Plastic behavior has been observed in rocks using all these techniques since the early 1900s.[2]

The Boudinage experiments [3] show that localized plasticity is observed in certain rock specimens that have failed in shear. Other examples of rock displaying plasticity can be seen in the work of Cheatham and Gnirk.[4] Test using compression and tension show necking of rock specimens while tests using wedge penetration show lip formation. The tests carried out by Robertson [5] show plasticity occurring at high confining pressures. Similar results are observable in the experimental work carried out by Handin and Hager,[6] Paterson,[7] and Mogi.[8] From these results it appears that the transition from elastic to plastic behavior may also indicate the transition from softening to hardening. More evidence is presented by Robinson [9] and Schwartz.[10] It is observed that the higher the confining pressure, the greater the ductility observed. However, the strain to rupture remains roughly the same at around 1.

The effect of temperature on rock plasticity has been explored by several teams of researchers.[11] It is observed that the peak stress decreases with temperature. Extension tests (with confining pressure greater than the compressive stress) show that the intermediate principal stress as well as the strain rate has an effect on the strength. The experiments on the effect of strain rate by Serdengecti and Boozer [12] show that increasing the strain rate makes rock stronger but also makes it appear more brittle. Thus dynamic loading may actually cause the strength of the rock to increase substantially. Increase in temperature appears to increase the rate effect in the plastic behavior of rocks.

After these early explorations in the plastic behavior of rocks, a significant amount of research has been carried out on the subject, primarily by the petroleum industry. From the accumulated evidence, it is clear that rock does exhibit remarkable plasticity under certain conditions and the application of a plasticity theory to rock is appropriate.

## Governing equations

The equations that govern the deformation of jointed rocks are the same as those used to describe the motion of a continuum:[13]

{ \begin{align} \dot{\rho} + \rho~\boldsymbol{\nabla} \cdot \mathbf{v} & = 0 & & \qquad\text{Balance of Mass} \\ \rho~\dot{\mathbf{v}} - \boldsymbol{\nabla} \cdot \boldsymbol{\sigma} - \rho~\mathbf{b} & = 0 & & \qquad\text{Balance of Linear Momentum} \\ \boldsymbol{\sigma} & = \boldsymbol{\sigma}^T & & \qquad\text{Balance of Angular Momentum} \\ \rho~\dot{e} - \boldsymbol{\sigma}:(\boldsymbol{\nabla}\mathbf{v}) + \boldsymbol{\nabla} \cdot \mathbf{q} - \rho~s & = 0 & & \qquad\text{Balance of Energy.} \end{align} }

where $\rho(\mathbf{x},t)$ is the mass density, $\dot{\rho}$ is the material time derivative of $\rho$, $\mathbf{v}(\mathbf{x},t) = \dot{\mathbf{u}}(\mathbf{x},t)$ is the particle velocity, $\mathbf{u}$ is the particle displacement, $\dot{\mathbf{v}}$ is the material time derivative of $\mathbf{v}$, $\boldsymbol{\sigma}(\mathbf{x},t)$ is the Cauchy stress tensor, $\mathbf{b}(\mathbf{x},t)$ is the body force density, $e(\mathbf{x},t)$ is the internal energy per unit mass, $\dot{e}$ is the material time derivative of $e$, $\mathbf{q}(\mathbf{x},t)$ is the heat flux vector, $s(\mathbf{x},t)$ is an energy source per unit mass, $\mathbf{x}$ is the location of the point in the deformed configuration, and t is the time.

In addition to the balance equations, initial conditions, boundary conditions, and constitutive models are needed for a problem to be well-posed. For bodies with internal discontinuities such as jointed rock, the balance of linear momentum is more conveniently expressed in the integral form, also called the principle of virtual work:

$\int_{\Omega} [\boldsymbol{\sigma}\cdot\nabla{\mathbf{w}} - \rho\,\mathbf{b}\cdot\mathbf{w} + \rho\,\dot{\mathbf{v}}\cdot\mathbf{w}]\,\text{dV} = \int_{\partial\Omega} \mathbf{t}\cdot\mathbf{w}\,\text{dS}$

where $\Omega$ represents the volume of the body and $\partial\Omega$ is its surface (including any internal discontinuities), $\mathbf{w}$ is an admissible variation that satisfies the displacement (or velocity) boundary conditions, the divergence theorem has been used to eliminate derivatives of the stress tensor, and $\mathbf{t}$ are surface tractions on the surfaces $\partial\Omega$. The jump conditions across stationary internal stress discontinuities require that the tractions across these surfaces be continuous, i.e.,

$\mathbf{n}\cdot\boldsymbol{\sigma}^{+} + \mathbf{n}\cdot\boldsymbol{\sigma}^{-1} = \mathbf{0} \qquad \text{or} \qquad \mathbf{n}\cdot[[\boldsymbol{\sigma}]] = \mathbf{0}$

where $\boldsymbol{\sigma}^{+},\boldsymbol{\sigma}^{-}$ are the stresses in the sub-bodies $\Omega^{+},\Omega^{-}$, and $\mathbf{n}$ is the normal to the surface of discontinuity.

### Constitutive relations

Stress-strain curve showing typical plastic behavior of rocks in uniaxial compression. The strain can be decomposed into a recoverable elastic strain ($\varepsilon_e$) and an inelastic strain ($\varepsilon_p$). The stress at initial yield is $\sigma_0$. For strain hardening rocks (as shown in the figure) the yield stress increases with increasing plastic deformation to a value of $\sigma_y$.

For small strains, the kinematic quantity that is used to describe rock mechanics is the small strain tensor $\boldsymbol{\varepsilon} = \tfrac{1}{2}\left[\nabla\mathbf{u} + (\nabla\mathbf{u})^T\right] \,.$ If temperature effects are ignored, four types of constitutive relations are typically used to describe small strain deformations of rocks. These relations encompass elastic, plastic, viscoelastic, and viscoplastic behavior and have the following forms:

1. Elastic material: $\,\,\boldsymbol{\sigma} = \mathsf{H}:\boldsymbol{\varepsilon}\,\,$ or $\,\,\sigma_{ij} = H_{ijkl}\,\varepsilon_{kl}\,\,$. For an isotropic, linear elastic, material this relation takes the form $\,\,\boldsymbol{\sigma} = 2\mu\,\boldsymbol{\varepsilon} + \lambda\,\text{tr}(\boldsymbol{\varepsilon})\,\boldsymbol{I}\,\,$ or $\,\,\sigma_{ij} = 2\mu\varepsilon_{ij} + \lambda\varepsilon_{kk}\delta_{ij}$. The quantities $\mu,\lambda$ are the Lamé parameters.
2. Viscous fluid: For isotropic materials, $\,\,\boldsymbol{\sigma} = -p\,\boldsymbol{I} + 2\mu\,\dot{\boldsymbol{\varepsilon}} + \lambda\,\text{tr}(\dot{\boldsymbol{\varepsilon}})\,\boldsymbol{I}\,\,$ or $\,\,\sigma_{ij} = -P\,\delta_{ij} + 2\mu\dot{\varepsilon}_{ij} + \lambda\dot{\varepsilon}_{kk}\delta_{ij}$ where $\mu$ is the shear viscosity and $\lambda$ is the bulk viscosity.
3. Nonlinear material: Isotropic nonlinear material relations take the form $\,\,\boldsymbol{\sigma} = 2\mu\,\boldsymbol{\varepsilon} + \lambda\,\text{tr}(\boldsymbol{\varepsilon})\,\boldsymbol{I} + \lambda'\,\boldsymbol{\varepsilon}\cdot\boldsymbol{\varepsilon}\,\,$ or $\,\,\sigma_{ij} = 2\mu\varepsilon_{ij} + \lambda\varepsilon_{kk}\delta_{ij} + \lambda'\,\varepsilon_{ik}\,\varepsilon_{kj}$. This type of relation is typically used to fit experimental data and may include inelastic behavior.
4. Quasi-linear materials: Constitutive relations for these materials are typically expressed in rate form, e.g., $\,\,\dot{\boldsymbol{\sigma}} = \mathsf{H}(\boldsymbol{\sigma}):\dot{\boldsymbol{\varepsilon}}\,\,$ or $\,\,\dot{\sigma}_{ij} = H_{ijkl}(\sigma_{mn})\,\dot{\varepsilon}_{kl}\,\,$.

A failure criterion or yield surface for the rock may then be expressed in the general form

$F(\boldsymbol{\sigma}, \dot{\boldsymbol{\sigma}}, \boldsymbol{\varepsilon}, \dot{\boldsymbol{\varepsilon}}, \mathbf{x}, t) = 0 \,.$

Typical constitutive relations for rocks assume that the deformation process is isothermal, the material is isotropic, quasi-linear, and homogenous and material properties do not depend upon position at the start of the deformation process, that there is no viscous effect and therefore no intrinsic time scale, that the failure criterion is rate-independent, and that there is no size effect. However, these assumptions are made only to simplify analysis and should be abandoned if necessary for a particular problem.

### Yield surfaces for rocks

View of Mohr–Coulomb failure surface in 3D space of principal stresses for $c=2, \phi=-20^\circ$

Design of mining and civil structures in rock typically involves a failure criterion that is cohesive-frictional. The failure criterion is used to determine whether a state of stress in the rock will lead to inelastic behavior, including brittle failure. For rocks under high hydrostatic stresses, brittle failure is preceded by plastic deformation and the failure criterion is used to determine the onset of plastic deformation. Typically, perfect plasticity is assumed beyond the yield point. However strain hardening and softening relations with nonlocal inelasticity and damage have also been used. Failure criteria and yield surfaces are also often augmented with a cap to avoid unphysical situations where extreme hydrostatic stress states do not lead to failure or plastic deformation.

View of Drucker–Prager yield surface in 3D space of principal stresses for $c=2, \phi=-20^\circ$

Two widely used yield surfaces/failure criteria for rocks are the Mohr-Coulomb model and the Drucker-Prager model. The Hoek–Brown failure criterion is also used, notwithstanding the serious consistency problem with the model. The defining feature of these models is that tensile failure is predicted at low stresses. On the other hand, as the stress state becomes increasingly compressive, failure and yield requires higher and higher values of stress.

## Plasticity theory

The governing equations, constitutive models, and yield surfaces discussed above are not sufficient if we are to compute the stresses and displacements in a rock body that is undergoing plastic deformation. A additional kinematic assumption is needed, i.e., that the strain in the body can be decomposed additively (or multiplicatively in some cases) into an elastic part and a plastic part. The elastic part of the strain can be computed from a linear elastic constitutive model. However, determination of the plastic part of the strain requires a flow rule and a hardening model.

Typical flow plasticity theories (for small deformation perfect plasticity or hardening plasticity) are developed on the basis on the following requirements:

1. The rock has a linear elastic range.
2. The rock has an elastic limit defined as the stress at which plastic deformation first takes place, i.e., $\sigma = \sigma_0$.
3. Beyond the elastic limit the stress state always remains on the yield surface, i.e., $\sigma = \sigma_y$.
4. Loading is defined as the situation under which increments of stress are greater than zero, i.e., $d\sigma > 0$. If loading takes the stress state to the plastic domain then the increment of plastic strain is always greater than zero, i.e., $d\varepsilon_p > 0$.
5. Unloading is defined as the situation under which increments of stress are less than zero, i.e., $d\sigma < 0$. The material is elastic during unloading and no additional plastic strain is accumulated.
6. The total strain is a linear combination of the elastic and plastic parts, i.e., $d\varepsilon = d\varepsilon_e + d\varepsilon_p$. The plastic part cannot be recovered while the elastic part is fully recoverable.
7. The work done of a loading-unloading cycle is positive or zero, i.e., $d\sigma\,d\varepsilon = d\sigma\,(d\varepsilon_e + d\varepsilon_p) \ge 0$. This is also called the Drucker stability postulate and eliminates the possibility of strain softening behavior.

### Three-dimensional plasticity

The above requirements can be expressed in three dimensions as follows.

• Elasticity (Hooke's law). In the linear elastic regime the stresses and strains in the rock are related by
$\boldsymbol{\sigma} = \mathsf{C}:\boldsymbol{\varepsilon}$
where the stiffness matrix $\mathsf{C}$ is constant.
• Elastic limit (Yield surface). The elastic limit is defined by a yield surface that does not depend on the plastic strain and has the form
$f(\boldsymbol{\sigma}) = 0 \,.$
• Beyond the elastic limit. For strain hardening rocks, the yield surface evolves with increasing plastic strain and the elastic limit changes. The evolving yield surface has the form
$f(\boldsymbol{\sigma}, \boldsymbol{\varepsilon}_p) = 0 \,.$
• Loading. It is not straightforward to translate the condition $d\sigma > 0$ to three dimensions, particularly for rock plasticity which is dependent not only on the deviatoric stress but also on the mean stress. However, during loading $f \ge 0$ and it is assumed that the direction of plastic strain is identical to the normal to the yield surface ($\partial f/\partial\boldsymbol{\sigma}$) and that $d\boldsymbol{\varepsilon}_p:d\boldsymbol{\sigma} \ge 0$, i.e.,
$d\boldsymbol{\sigma}:\frac{\partial f}{\partial \boldsymbol{\sigma}} \ge 0 \,.$
The above equation, when it is equal to zero, indicates a state of neutral loading where the stress state moves along the yield surface without changing the plastic strain.
• Unloading: A similar argument is made for unloading for which situation $f < 0$, the material is in the elastic domain, and
$d\boldsymbol{\sigma}:\frac{\partial f}{\partial \boldsymbol{\sigma}} < 0 \,.$
• Strain decomposition: The additive decomposition of the strain into elastic and plastic parts can be written as
$d\boldsymbol{\varepsilon} = d\boldsymbol{\varepsilon}_e + d\boldsymbol{\varepsilon}_p \,.$
• Stability postulate: The stability postulate is expressed as
$d\boldsymbol{\sigma}:d\boldsymbol{\varepsilon} \ge 0 \,.$

### Flow rule

In metal plasticity, the assumption that the plastic strain increment and deviatoric stress tensor have the same principal directions is encapsulated in a relation called the flow rule. Rock plasticity theories also use a similar concept except that the requirement of pressure-dependence of the yield surface requires a relaxation of the above assumption. Instead, it is typically assumed that the plastic strain increment and the normal to the pressure-dependent yield surface have the same direction, i.e.,

$d\boldsymbol{\varepsilon}_p = d\lambda\,\frac{\partial f}{\partial \boldsymbol{\sigma}}$

where $d\lambda > 0$ is a hardening parameter. This form of the flow rule is called an associated flow rule and the assumption of co-directionality is called the normality condition. The function $f$ is also called a plastic potential.

The above flow rule is easily justified for perfectly plastic deformations for which $d\boldsymbol{\sigma} = 0$ when $d\boldsymbol{\varepsilon}_p > 0$, i.e., the yield surface remains constant under increasing plastic deformation. This implies that the increment of elastic strain is also zero, $d\boldsymbol{\varepsilon}_e = 0$, because of Hooke's law. Therefore,

$d\boldsymbol{\sigma}:\frac{\partial f}{\partial \boldsymbol{\sigma}} = 0 \quad \text{and} \quad d\boldsymbol{\sigma}:d\boldsymbol{\varepsilon}_p = 0 \,.$

Hence, both the normal to the yield surface and the plastic strain tensor are perpendicular to the stress tensor and must have the same direction.

For a work hardening material, the yield surface can expand with increasing stress. We assume Drucker's second stability postulate which states that for an infinitesimal stress cycle this plastic work is positive, i.e,

$d\boldsymbol{\sigma}: d\boldsymbol{\varepsilon}_p \ge 0 \,.$

The above quantity is equal to zero for purely elastic cycles. Examination of the work done over a cycle of plastic loading-unloading can be used to justify the validity of the associated flow rule.[14]

### Consistency condition

The Prager consistency condition is needed to close the set of constitutive equations and to eliminate the unknown parameter $d\lambda$ from the system of equations. The consistency condition states that $df = 0$ at yield because $f(\boldsymbol{\sigma},\boldsymbol{\varepsilon}_p) = 0$, and hence

$df = \frac{\partial f}{\partial \boldsymbol{\sigma}}:d\boldsymbol{\sigma} + \frac{\partial f}{\partial \boldsymbol{\varepsilon}_p}:d\boldsymbol{\varepsilon}_p = 0 \,.$

## Notes

1. ^ Pariseau (1988).
2. ^ Adams and Coker (1910).
3. ^ Rast (1956).
4. ^ Cheatham and Gnirk (1966).
5. ^ Robertson (1955).
6. ^ Handin and Hager (1957,1958,1963.)
7. ^ Paterson (1958).
8. ^ Mogi (1966).
9. ^ Robinson (1959).
10. ^ Schwartz (1964).
11. ^ Griggs, Turner, Heard (1960)
12. ^ Serdengecti and Boozer (1961)
13. ^ The operators in the governing equations are defined as:
\begin{align} \boldsymbol{\nabla} \mathbf{v} &= \sum_{i,j = 1}^3 \frac{\partial v_i}{\partial x_j}\mathbf{e}_i\otimes\mathbf{e}_j = v_{i,j}\mathbf{e}_i\otimes\mathbf{e}_j \\ \boldsymbol{\nabla} \cdot \mathbf{v} & = \sum_{i=1}^3 \frac{\partial v_i}{\partial x_i} = v_{i,i} \\ \boldsymbol{\nabla} \cdot \boldsymbol{S} &= \sum_{i,j=1}^3 \frac{\partial S_{ij}}{\partial x_j}~\mathbf{e}_i = S_{ij,j}~\mathbf{e}_i ~. \end{align}
where $\mathbf{v}$ is a vector field, $\boldsymbol{S}$ is a symmetric second-order tensor field, and $\mathbf{e}_i$ are the components of an orthonormal basis in the current configuration. The inner product is defined as
$\boldsymbol{A}:\boldsymbol{B} = \sum_{i,j=1}^3 A_{ij}~B_{ij} = \operatorname{trace}(\boldsymbol{A}\boldsymbol{B}^T) ~.$
14. ^ Anandarajah (2010).