# Numerical modelling - constitutive models

A Constitutive model introduces and describes the stress-strain behavior in response to applied loads, relevant to post-failure. This model is used to predict plastic and elastic straining of the material and the physical properties of a given material^{[1]}. In pre-failure, most constitutive models have simple linear elastic behavior and complexity occurs in the post-failure region. So the rate of change of the strength of the material is a significant factor and should be considered in modeling the materials^{[2]}.

The behavior of rock materials generally evolves in three possible ways which are related to the unstable physical system (Figure 1). Average quality rock material tends to show a strain softening behavior after failure occurred similar to the one shown in Figure 1c. The stored energy after failure released kinetic energy, which then radiates away from the source and dissipates. Numerical solution may fail to converge in this situation. A very poor quality rock material shows an elastic-perfectly plastic behavior, see Figure 1b. Hard, good quality rock material tends to show an elastic-brittle behavior, in which the strength drops rapidly, once the material is introduced to plastic straining (Figure 2) ^{[3]}.

The path dependency of nonlinear materials should be considered in modelling. The numerical solution should be able to accommodate different loading paths in order to apply the constitutive model properly. Different loading path should be specified in order for a correct solution to implement the constitutive model. For example, if an excavation is made suddenly, then the solution may be influenced by internal effect that caused additional failure of the material. But this may not be seen, if the excavation is made gradually.
The nonlinearity of the stress-strain response is another factor of material characteristic in modeling. As mentioned before the most complexity constitutive models occurred in the post-failure region. This includes the nonlinear dependence of both the elastic stiffness and the strength envelope on the confining stress. In this situation rock mass behaviour changed according to the stress level (tensile, unconfined and confined regimes) after failure
^{[1]}.

The numerical simulation in geomechanics needs to be able to accommodate these three various forms of difficulties. The evolution of a geologic system in a realistic manner should be considered in order to model the stable numerical solution in which kinetic energy could be generated and dissipated. On the other hand, the law of motion should be included to determine the effect of loading path on the constitutive response. Also dynamic solution in nonlinear constitutive models is needed to calculate quantities in sequence order (forces/stresses and velocities/displacements). It means each element in the model is isolated physically from one another during one calculation step ^{[1]}.

Established constitutive models are classified in three groups depending on mechanical constructs:

1- Algebraic constitutive relations, such as linear and nonlinear elastic model
2- Differential constitutive relations, such as plastic and viscoelastic models
3- Integral constitutive relations, such as hereditary viscoelastic models ^{[4]}.

The strains of a material when exposed to a load can basically be divided into two separate parts, Elastic and Plastic strains. In an elastic material strain developed during loading however in a plastic material strains remained after the material has been unloaded. Plastic strains start to develop once the material reaches its yielding limit, which is dictated by the yield criterion. The major problem in practice is firstly knowing which representation is relevant for the particular case being analyzed, and secondly, how to determine the relevant parameters to describe that behavior ^{[2]}.

The basic constitutive models in the form of linear elastic and nonlinear plastic models are described as below. Two well-known examples of nonlinear plastic constitutive models are the Hoek-Brown and the Mohr-Coulomb models. The parameters defining peak strength are different for each, although the elastic moduli are still represented in the same way i.e. Young’s modulus E and Poisson's ratio.

## Null Model Group

A null material model is used to represent material that is removed or excavated that might be changed in a later stage of simulation such as backfilling. The stresses in this model are set to zero and no body forces even gravity act on these models ^{[1]}.

## Elastic Model Group

The models in this group are characterized by reversible deformations upon unloading; the stress-strain laws are linear and path-independent.

### Elastic, isotopic model

The elastic, isotopic model is the main simplest model of material behavior. This model is considered for homogenous, isotropic and continues materials. The stress- strain behavior is linear with no hysteresis on unloading and is expressed by Robert Hooke law with the deformation increased with the applied forces:

σ = E ε where σ is stress, E is Young’s modulus, ε is strain

The definition of the modulus of elasticity E is attributed to Thomas Young which is measured a normalized force density and a normalized deformation ^{[4]}.

### Elastic, transversely isotropic model

The elastic, transversely isotropic model has a plane of elastic symmetry. It means all directions in a plane, normal and parallel, are elastically equivalent. Then the material properties are the same in all directions. The model is built based on the transversely isotropic elastic material related to plain-strain or plane-stress in two dimension xy plane and oriented to the z-axis which is parallel with the plain isotropy ^{[1]}.

## Plastic Model Group

A consequence of the nonlinearity of the stress-strain relations which is caused by failure criteria is categorized in plastic group. The different models are characterized by their yield function, hardening/softening functions and plastic flow. The plastic flow formulation based on basic assumptions that the total strain may be divided into elastic and plastic parts. The flow rule specifies the direction of the plastic strain increment vector as that normal to the potential surface.

Different plasticity models are classified based on the shear yield and potential functions, non-associated flow rules and stress corrections. In this case, only the elastic part of the strain increment can contribute to the stress increment; the latter is corrected by using the plastic flow rule to ensure that the stresses lie on the composite yield function ^{[1]}.

Different plastic models are described below:

### Drucker-Prager Model

This model is mostly used to model soft clays with low friction angles and deals with the plastic deformation of soils. The failure envelope in this model is a Drucker-Prager criterion with tension cutoff. The simplification of the Mohr-Coulomb model where the hexagonal shape of the failure cone was replaced by a simple cone is known as the Drucker-Prager model. The Drucker-Prager criterion is determined whether a material has failed or undergone plastic yield. This model is the smooth version of Mohr-Coulomb yield surface.
The yielding surface of the Drucker-Prager criterion depends on the internal friction angle of the material and its cohesion ^{[5]}. This model is expressed in terms of two generalized stress components, the tangential stress,τ and mean normal stress,σ.
^{[6]}.

### Mohr-Coulomb model

The Mohr-Coulomb model is the most common constitutive model which is described a linear relationship between the shear stress and normal stress (or maximum and minimum principal stresses) of rock or soil at failure in response to applied loads. The combination and generalization of Hooke and Coulomb’s law are formulated in a plasticity framework and is known as Mohr-Coulomb mode. In general stress state, the model’s stress-strain behaves linearly in the elastic range, with two defining parameters from Hooke’s law (Young’s modulus, E and Poisson’s ratio, ν). There are two parameters which define the failure criteria (the friction angle, φ and cohesion, c) and also a parameter to describe the flow rule (dilatancy angle, ψ which comes from the use of non-associated flow rule which is used to model a realistic irreversible change in the volume due to shearing). In the conventional plastic theory, the flow rule is used as the evolution law for plastic strain rates. If the plastic potential function is the same as the yield function, the flow rule is called the associated flow rule and if it is different, it is called the non-associated flow rule ^{[1]}.
Only two strength parameters (normal and shear stress) described the plastic behavior in the Mohr-Coulomb model. These are used to represent shear failure and tensile yield function in soils and rocks at different effective stresses. Many geotechnical analysis methods and programs require use of this strength model due to the simplicity of this model. The Mohr-coulomb criterion states that the shear stress, τ, and the normal stress, σn, is related in a plane in a continuum, which satisfies the yielding occurs c is the intercept of the failure envelope with the τ axis and it is called Cohesion, and φ is the stope of the failure envelope which is called friction angle. (Figure 4).

τ= c + σn tan φ

In the above expression, tension is taken as positive. The principal stresses σ1, σ2, σ3 are used to implement this model and σ1 is the maximum principal stress and σ3 is the minimum principal stress. ^{[3]}. ^{[7]}.

### Ubiquitous-Joint Model

The ubiquitous-joint model represents the weak plane, such as weathering joints and bedding planes, orientation in an anisotropic plasticity model in a Mohr-Coulomb solid in which yielding is occurring in either along the weakness, the solid or both. The failure criteria depend on the orientation of weak plane, material properties of the solid or weak plane and the stress. The term "ubiquitous" means that the joints do not have a fixed location and may occur at any point in the rock mass ^{[1]}.
The normal and shear stress are calculated at any orientation of the weak plane by giving the joint orientation for general failure then displacement, factor of safety and the stability of the excavation is evaluated by using Mohr-Coulomb failure criteria on the weak plane. Then the numerical modelling is taken to suggest the solution under practical geometries and boundary conditions. In the above process, the joints themselves do not influence the stress/displacement field, they are simply a post-processing tool to determine whether the stresses generated from the finite element analysis are sufficient to cause slip on a plane at a defined orientation ^{[8]}.

Figure 5 illustrates the weak plane existing in a Mohr-Coulomb model and the global (x-y) and local (x’-y’) coordinate frame.

### Strain-Hardening/Softening Model

This model represents the nonlinear behavior of a strain softening or hardening material. Rock material tends to lose some of its strength once plastic straining occurs, which is known as softening. However many metals tend to show an increase in strength during plastic straining (Figure 1 a) which is known as hardening. In the hardening case, once yield occurs, the stress needs to be increased in order to drive the plastic deformation. The yield surface may mostly change size, shape and position. The yield surface remains unchanged in the perfectly plastic case and simply depends on the hardening parameters. In order to employ hardening model two methods are considered: isotropic hardening model (Figure 6) and the kinematic hardening model (Figure 7). The difference between these two models is that kinematic hardening shifts the yield surface from one location in stress space to another although in isotropic hardening model, this position and shape is unaltered ^{[3]}.

In the softening case, the stress-strain curve “turns down”, as in Figure 1 c. The yield surface in this case generally will decrease in size with further straining. Strain hardening/softening of the compression and tension yield surface depends on the minor and major principal plastic strains and same as that is used for the Mohr-Coulomb model. This behavior is controlled by the variation in friction, cohesion and dilation as a function of plastic shear strain, and tension limit as a function of plastic tensile strain. The shear hardening parameter is used to calculate the plastic shear strain and the tensile hardening parameter is measured the accumulated tensile plastic strain ^{[1]}.

### Bilinear Strain-Hardening/Softening Ubiquitous-Joint Model

The bilinear strain-hardening/softening ubiquitous-joint model represents the softening and hardening behavior of material for the matrix of failure envelope and the weak plane incorporated with the ubiquitous-joint model properties (cohesion, friction, dilation, and tensile strength) as the functions of deviatoric and tensile plastic strain. On the other hand, if ubiquitous-joints are introduced, then the variation of material strength properties with mean stress can also be taken into account by using the bilinear option. The softening behaviors for the matrix and the joint are specified which measure the amount of plastic shear and tensile strain. It means in this numerical model, general failure is first detected, and relevant plastic corrections are applied. The new stresses are then analyzed to failure into parallel and perpendicular to the weak plane and updated accordingly for ubiquitous-joint failure ^{[1]}.

This model is mostly used in the post failure of laminated materials that exhibit nonlinear material hardening or softening.

### Double-Yield Model

The double-yield model is intended to represent materials in which there may be significant irreversible compaction in addition to shear yielding, such as hydraulically placed backfill or lightly cemented granular material. This model is developed to characterize the strength and deformation behaviors of coarse granular materials. Two kinds of plastic deformation mechanisms, including the shear and compressive plastic deformation, are taken into account in this mode.

### Hoek-Brown Model

The Hoek-Brown criterion is an empirical non-linear refinement of the Mohr-Coulomb criterion and is specifically designed for rock-like materials. The Hoek-Brown is based on stress conditions in intact rock and rock masses that lead to failure criterion. The failure surface is nonlinear shear-yield function, and is based on the relation between the major and minor effective principal stresses (Figure 8). A plasticity flow rule that varies as a function of the stress level. The Mohr-Coulomb friction and cohesion equivalent is matched to the nonlinear Hoek-Brown strength envelope in particular stress levels. In addition to a failure (or yield) criterion, a “flow rule” is also necessary, in order to provide a relation between the components of strain rate at failure. The flow rule is mostly assumed to be fixed in failure criterion that lead to limiting stress condition that corresponds to the ultimate failure of the material.

## Summery

The most common elastic-plastic constitutive models are summarised in the below table. Representative material and application are also shown in this table.

## References

- ↑
^{1.0}^{1.1}^{1.2}^{1.3}^{1.4}^{1.5}^{1.6}^{1.7}^{1.8}^{1.9}Itasca Consulting group "UDEC, Universal Distinct Element Code 2011" . - ↑
^{2.0}^{2.1}Strafield M. & Cundall P. (1988) "Towards a methodology for rock mechanics modelling International Journal of Rock Mechanics and Mining Science" . - ↑
^{3.0}^{3.1}^{3.2}^{3.3}^{3.4}^{3.5}Smed E.M. & Cundall P. (2012) "Elasto-plasto Strain Hardeninh Mohr-Coulomb Model-Derivation and Implementation, Aalborg, Denmark" . Cite error: Invalid`<ref>`

tag; name "Smed" defined multiple times with different content Cite error: Invalid`<ref>`

tag; name "Smed" defined multiple times with different content Cite error: Invalid`<ref>`

tag; name "Smed" defined multiple times with different content - ↑
^{4.0}^{4.1}Kasper. W.J "Constitutive Models for Engineering Materials, (2002) University of Colorado vol.3" . - ↑ Kelly. P "Solid Mechanics Lecture Notes. New Zealand : The university of Auckland (2014) " .
- ↑
^{6.0}^{6.1}Raisa. E "A study of Geotechnical Constitutive Models using PLAXIS 2D. s.l. : 2012/13 ICE G&S Papers Competition Entry,(2012) " . - ↑ Wikipedia, the free encyclopedia. http://en.wikipedia.org/wiki/Mohr%E2%80%93Coulomb_theory .
- ↑
^{8.0}^{8.1}Salvador I. Drucker-Prager yield criterion application to study the behavior of CFRP confined concrete under compression . s.l. : Department of Construction Engineeering, University of Alicante, Spain, 2010 . - ↑
^{9.0}^{9.1}Gang L. & Hong P. & Harumi K. & Yoshiaki M."Application of Ubiqutios Joints Models in Numerical Modeling, Chinese Journel of Rock Mechanics and Engineering, 2003" . - ↑ PLAXIS "Material Models Manual, 2014" .
- ↑ Itasca "http://icsas.itascacg.com/software/udec/features/materials" .