© J.-J. Marigo and K. Kazymyrenko, published by EDP Sciences 2019

1 Introduction

Many works have been devoted to the modeling of the mechanical behavior of geomaterials like concrete, rocks, and soils and many phenomenological models have been proposed to give an account of the main aspects of the observed phenomena, namely, dilatancy and stress softening, while those materials are submitted to triaxial compressions. The majority of these approaches is based on elasto-plastic formalism with some of them having resort to damage coupling. For instance, in soil or rock mechanics, one generally uses elasto-plastic or more generally (thermo)elasto-viscoplastic models without introducing damage variables [1–3], whereas it is generally admitted that one cannot reproduce the experimental results observed for concrete without considering a damage evolution law [4–10]. In all the cases the elasto-plastic models are based on a restricted choice of yield criteria (like Drucker–Prager one or some variants like Hoek–Brown criterion [1,2] or the cap model [11]), with common feature being the hypothesis of a nonassociative flow rule. The reason generally invoked to put aside the normality rule for the plasticity evolution is that otherwise a too important dilatancy effects would be obtained in the evolution of the volumetric strain once the normality rule is used, for instance, for a standard Drucker–Prager-type criterion.

This, evoked earlier, macroscopic elasto-plastic coupled behavior is interpreted on a lower scale as the consequence of the presence of microcracks inside the material. Because of the compression, a part or all of those cracks are closed and a friction between the lips of the cracks can prevent a free sliding of the lips. This impossibility to relax to the initial strain state results at the macro-level in the presence of residual plastic strains. Moreover, if one assumes that the sliding with friction between the lips follows Coulomb law, then a Drucker–Prager-type criterion for the plasticity law is naturally obtained. Finally, the growth of the microcracks is in straight relation to the observed progressive loss of rigidity of the samples and can be reproduced at the macro-level by using a damage model. But if those mechanisms related to the microcracking are usually invoked to build macroscopic models, very few works have been devoted to an authentic micromechanical approach. An exception is reference [12] where such an attempt is made by using Coulomb friction law for the sliding of the crack lips and Griffith law for the crack propagation; see also a more recent paper [13] where a full 3D approach is proposed. Even if the analysis in reference [12] is conducted in a simplified two-dimensional setting where the uncracked material is assumed isotropic, elastic, and homogeneous whereas the microcracks are assumed straight and small enough to neglect their interactions, some general properties are obtained and can be used to construct more general macroscopic models. In particular, the authors show that the elastic energy stored in the microcracked material reads as the following function of the total strain  and the plastic strain :(1)

In equation (1), one uses the notations precised at the end of this section. In particular,  denotes the stiffness tensor of the sound material, whereas  represents a loss of stiffness due to the microcracks. Therefore, the micromechanical approach leads to a term in the elastic energy which corresponds to the energy blocked by the contact with friction of the lips of the cracks. Let us note that this latter term involves only the plastic strain (and eventually the damage state via  dependence) not the total strain. The damage state is supposed not to enter in the first term of the energy equation (1) and, consequently, even if the stress–strain relation still reads asthe stress which appears in the plasticity yield criterion deduced from Coulomb friction law is not  but the tensor  given by

In other words, the yield criterion is of the Drucker–Prager type with a kinematical hardening where the back stress  depends linearly on the plastic strain but also nonlinearly on the damage state. These two particularities (the presence of a kinematical hardening in the plasticity yield criterion and the fact that damage appears only in the blocked elastic energy) can generally not be seen in the various models proposed in the literature; see however references [10,14] where kinematical hardening is also introduced. It turns out that they have fundamental and rather unexpected consequences on the response of the material under triaxial compression tests. In particular, one can account for reasonable contractance and dilatancy effects without leaving an associate flow rule. The goal of the present paper is to study the specificity of the responses of this type of models.

The paper is organized as follows. In Section 2 we con sider the elasto-plastic model without damage. Presenting first the general ingredients of an associative plasticity law, we progressively refine the model by introducing the particularities coming from the micromechanical considerations. We finally study the response of a volume element submitted to a triaxial test at fixed confining pressure. The damage is introduced in Section 3, still on the basis of the micromechanical considerations. Using a standard law for the damage evolution law, we obtain the complete model of elasto-plasticity coupled with damage. After establishing some general properties, we then consider a family of particular models to finally calculate the response of the volume element submitted to the triaxial test at fixed confining pressure.

Throughout the paper, we use the following notations (see also Tab. 1). The summation convention on repeated indices is implicitly adopted. The vectors and second-order tensors are indicated by boldface letters, like  and  for the unit normal vector and the stress tensor. Their components are denoted by italic letters, like n _i and σ _ij . The fourth-order tensors as well as their components are indicated by a sans serif letter, like  or  for the stiffness tensor. Such tensors are considered as linear maps applying on vectors or second-order tensors and the application is denoted without dots, like  whose ij-component is . The inner product between two vectors or two tensors of the same order is indicated by a dot, like  which stands for a _i b _i or  for σ _ij ε _ij .

2 The elasto-plastic behavior at fixed damage state

This section describes a pure elasto-plastic behavior law. Even so the damage is not explicitly present in this formulation, the latter could be considered as the model at fixed damage state. The idea is to identify all the relevant coupling terms which influence dilatancy and softening behavior.

2.1 General considerations

An associative model of elasto-plasticity is defined by two potentials:

1. The volume free energy function of the state variables.

2. The dissipation potential or equivalently the convex set in which the thermodynamical forces associated with the internal variables must lie.

Here we will only consider isothermal processes and a plastic behavior with linear kinematical hardening. Accordingly, the state variables are constituted by the strain tensor  and the plastic tensor  without other internal variables (the damage variable will be introduced in Section 3). Therefore, the free energy density ψ is given by(2)where the state function is assumed to be convex and at least continuously differentiable. By differentiation, one obtains the stress tensor  and the tensor  of the thermodynamical forces associated with the plastic strain tensor:(3)

From the mathematical point of view, , , , and  are order 2 symmetrical tensors which can be identified with 3 × 3 symmetrical matrices after the choice of an orthonormal basis. In other words, those tensors will be considered as elements of .

The plasticity (or yield) criterion is defined by giving the (closed with non-empty interior) convex set  of  where the thermodynamical forces  must lie. Here we make the strong assumption that this domain is fixed, that is independent of time and of the plasticity evolution. (However, we will see in the next section that the image of this domain in the stress space evolves with the plastic strain.) In practice, the convex  is characterized by a convex function  so that the plasticity criterion reads as(4)

The interior points of  correspond to the forces  such that  and the points on the boundary of  to the forces  such that .

Since we only consider an associative model, the evolution of the plastic strain follows the normality rule. When the convex set  has a smooth boundary without angular points, the function f is differentiable and the flow rule can read as(5)where  denotes the plastic multiplier. If the boundary of  is not smooth and the normal is not defined at some points, then the normality rule is extended by using the Hill maximal work principle. In such a case, the flow rule is given by(6)

The inequality (6) says that  must belong to the cone of the outer normals to  at the point  of the boundary. That inequality leads to equation (5) when the normal is well defined at .

From the energetic point of view, the energy which is dissipated by plasticity can be defined in terms of the support function of the convex . Specifically, let us introduce the support function  of :(7)

By construction,  is convex and positively homogeneous of degree 1,

Moreover,  vanishes at , is nonnegative provided that , and can take the value +∞ when  is not bounded. Its interpretation as the volume dissipated power by plasticity comes from Clausius–Duhem inequality. Indeed, in isothermal condition, the dissipated power (by volume unit)  is defined byand the second principle of thermodynamics requires that  be nonnegative. Owing to equations (2) and (3),  reads as

By virtue of Hill maximal work principle (6) and the definition (7) of the support function, one gets(8)which means that  plays also the role of the dissipation potential. Furthermore, it is sufficient to have  for Clausius–Duhem inequality to be satisfied, i.e., .

2.2 Choice of the form of the free energy state function

Guided by the micromechanical considerations presented in reference [12] for microcracked materials, we consider that the plastic strain is due to the friction between the lips of the (closed) microcracks. However, we do not assume that  and hence abandon the usual plastic incompressibility condition. Consequently, the plastic strain is a full symmetric second-order tensor and no more a pure deviator. Moreover, the spherical part of  which comes from the normal displacements between the lips of the microcracks will be also governed by an irreversible evolution law. Accordingly, we choose for the free energy the following quadratic function of :(9)

In equation (9) the first term represents the elastic energy and A ₀ denotes the (fourth-order) stiffness tensor of the sound material. The second term in equation (9) represents the blocked energy by friction and the (fourth-order) stiffness tensor A ₁ represents the loss of stiffness due to the presence of microcracks. As one will see later, this loss of stiffness is visible only when the sliding or the opening between the lips of the cracks are active, i.e., when .

The tensor A ₀ is positive definite and the tensor A ₁ is assumed to be nonnegative. In fact, A ₁ is also positive definite as soon as the material is microcracked, as it is shown in reference [12]. But, we will also consider the case when A ₁ vanishes to emphasize the strong influence of the blocked energy on the mechanical behavior of the microcracked material. Consequently,  is a (strictly) convex function of .

To simplify the presentation, we will assume that the material remains isotropic (even when the damage grows). Therefore, both tensors A ₀ and A ₁ have only two independent moduli. Decomposing the strain and the plastic strain tensors into their spherical and deviatoric parts(10)the free energy finally reads as(11)where K₀ > 0 and μ ₀ > 0 are the compressibility and shear moduli of the sound material, whereas K₁ ≥ 0 and μ ₁ ≥ 0 will be called the kinematical hardening moduli. The elastic moduli K₀ and μ ₀ are related to the Young modulus E₀ and the Poisson ratio ν ₀ of the sound material by(12)

By differentiation of equation (11), the stress–strain relation is linear and reads as(13)

Decomposing the stress tensor into its spherical and deviatoric parts,(14)the relation (13) becomes(15)

Still by differentiation, the thermodynamical forces associated with the plastic strains read as(16)and they differ in general from the stresses when A ₁ ≠ 0. Decomposing the forces X into their spherical and deviatoric parts,(17)the relation (16) becomes(18)

2.3 Choice of the plasticity criterion

Drucker–Prager criterion (or some variants like Mohr–Coulomb or Hoek–Brown criteria) is widely used in the modeling of the behavior of geomaterials under compression, see references [1,2,5]. In general, those criteria are used in a nonstandard context, i.e., with a nonassociative flow rule for the plasticity evolution, like in reference [3]. The main usual argument which is evoked to use a nonassociative flaw rule is that one overestimates the dilatancy effect with the normality rule associated with Drucker–Prager criterion. On the contrary, we propose in this paper to keep an associative flaw rule but with a Drucker–Prager criterion which contains a kinematical hardening. One of the main goals is to show how a kinematical hardening can lead to relevant dilatancy effects.

Specifically, in the present context, the Drucker–Prager criterion reads in terms of the thermodynamical forces X as(19)

In equation (19) k ∈ (0, 1) is a dimensionless coefficient associated with the internal friction, τ _c  ≥ 0 represents a critical stress which vanishes in absence of cohesion, and ∥· ∥ denotes the euclidean norm of a second-order tensor (whereas the factor  is introduced for convenience in order to simplify forthcoming expressions),

Consequently, the elastic domainis really a closed convex subset of  (with non-empty interior). In fact,  is a convex cone which has an angular point at  and which is unbounded in the direction of hydrostatic compressions. Furthermore, even if the elastic domain is fixed in the space of the thermodynamical forces, it varies with the plastic strain in the stress space. Indeed, by virtue of equation (16), the elastic domain expressed in term of  becomes the set  given by

Thus,  is subjected to a linear translation when the plastic strain evolves:

Starting from equation (7), the support function  of the convex  is obtained after tedious calculations which are not reproduced here and finally reads as(20)

Since we adopt the normality rule, the evolution of the plastic strain is directly deduced from equation (19). Considering first a regular point of the boundary of , the flow rule reads as follows:(21)

Let us note that at such a point one has . At the angular point, the rate of the plastic strain tensor must belong to the cone of the outer normal to . That leads to(22)

Consequently, the normality rule forces the trace of the plastic strain to only increase with time, a fundamental property to account for the dilatancy effects.

2.4 Response under a triaxial test with a confining pressure

The elasto-plastic behavior of a material governed by the associative Drucker–Prager law presented in the previous sections is illustrated by considering the response of a volume element during a triaxial test with a confining pressure. Assuming that the volume element starts from a natural reference configuration without plastic strain, the test is divided into the following three stages:

1. First, the volume element is submitted to a hydrostatic compression where the pressure is progressively increased from 0 to a final value p ₀> 0.

2. Then, maintaining the lateral pressure to the value p ₀, one compresses in the axial direction z by prescribing the axial strain ε _z with . Accordingly, the response will be purely elastic as long as |ε _z | is small enough so that the plastic yield criterion is not reached.

3. Finally, whe|ε _z | is larger than a critical value, the volume element plastifies.

Let us determine the evolution of the strains, the plastic strains, and the stresses during those different stages with respect to |ε _z | (which could be considered as the loading parameter).¹

2.4.1 Confining stage

During this stage the strain and stress tensors are purely spherical and the plastic strain remains equal to 0. At the end of this stage, the stresses and the strains are given by(23)

Since p = 0, one has  and Ducker-Prager criterion gives

The axial strain ϵ _z decreases from 0 to −p ₀/3K₀, whereas the volumetric strain  decreases from 0 to −p ₀/K₀. During all this stage, one hasand hence a contractance due to the hydrostatic compression.

2.4.2 Elastic stage at fixed confining pressure

During this stage, the plastic strain remains equal to 0:

By symmetry, the stress and the strain tensors are of the following form²:where p ₀ and ε _z are prescribed. The deviatoric parts read aswith σ _z  and ε _T remain to be determined. They are given by equation (15) with the help (12):(24)

One deduces in particular the evolution of the volumetric strain:(25)

This stage stops when the stresses reach the yield criterion, i.e., when . Sinceand since σ _z  + p ₀ < 0, the yield surface is reached when(26)

The corresponding value of ϵ _z is deduced from equation (24):(27)

2.4.3 Plastification stage

Then, if the confining pressure is constant and if one prescribes , the material will plastify. This plastification stage strongly depends on the presence or not of the kinematical hardening. Therefore, we will distinguish between the two situations by considering first the case where the hardening moduli μ ₁ and K₁ vanish, i.e., the case without hardening.

1. Case without kinematical hardening. In such a case, one gets  and since the elastic domain is fixed in the stress space, the stresses will remain blocked to the value reached at the end of the elastic stage. In other words, the lateral pressure is equal to p ₀ and the axial stress σ _z is given by equation (26). Since , one getsand hence p takes the form(28)

   Since σ _z  + p ₀ < 0, the flow rule (21) gives

   Since  and , one gets

   Therefore, the volumetric strain rate is given by(29)

   Thus, the normality rule predicts a dilatancy of the material during the plastification stage. In general, this predicted dilatancy is too large by comparison with experimental observations and that leads people to use a nonassociative flow rule. We will see just below that the kinematical hardening can strongly modify these results.

2. Case with kinematical hardening. In such a case, the stresses do not remain blocked. Their evolution is obtained by assuming that the strain tensor remains of the form (28) and that the yield criterion is always satisfied at a non-angular point, i.e.,

   Let us note that this assumption is not restrictive because we can use the general result of uniqueness of the response in the case of an associative law with kinematical hardening. Indeed, if we find a solution, then we are sure that is the good one. (Note that such a uniqueness result is no more ensured in the case of nonassociative flow rule.)

Accordingly, one getsand

Differentiating these relations with respect to time and taking into account equation (18) lead to

The flow rule (21) gives(30)

Combining the four previous relations gives the following relationship between  and :(31)

Differentiating the stress–strain relations (13) gives

Eliminating and using the flow rule give two other relations between , , , and :(32)

Using equation (31) allows us to determine , , and in terms of . We finally get(33) (34)

Therefore, equation (33) shows that the kinematical hardening tends to decrease  and even equation (34) shows that the sign of  depends on the intensity of the hardening. Specifically, if the hardening moduli are small enough (by comparison with the elastic moduli), then  and there is dilatancy. But, if the hardening moduli are large enough, then  and there is con tractance. Moreover, let us note that if μ ₁ = K₁ = 0, then we recover that  and the relation (29). But, if µ ₁ and K₁ are very large by comparison to µ ₀ and K₀, then one gets

which is nothing but the relation (25) obtained in the elastic stage for the volumetric strain change. In conclusion, according to the ratio between the hardening moduli and the elastic moduli, the volumetric strain evolution can lead to a strong dilatancy effect corresponding to perfect plasticity as well as a contractance effect like in elasticity. This last conclusion leads to a foretaste consideration (developed in detail in the next section) for a relevant damage coupling − once the damage affects A ₁ tensor, it directly impacts the dilatancy strength of the given law.

3 The elasto-plastic behavior coupled with damage

In this section we introduce the damage into the elasto-plastic model. As mentioned at the end of Section 2, the coupling could be introduced purely from macroscopic considerations. We made a choice for more physical approach starting from a set of hypothesis which will finally lead to the same coupling. Let us first present a list of assumptions used in that construction before particularizing the model.

3.1 Main assumptions on the damage dependence

Although there exist an infinite number of ways to couple damage with plasticity [5–9], we will make the simplest choices. Specifically, the following are the main assumptions:

1. The damage state is characterized by a scalar variable α which can grow from 0 to 1, α = 0 corresponding to a sound material, and α = 1 corresponding to a full damaged material.

2. The damage state enters in the expression of the free energy, but not in the plasticity criterion. In other words, the convex domain  is independent of α.

3. The dissipated energy by damage is only function of the damage state.

By virtue of these assumptions, the state of the material point is characterized by the triple  and the free energy is still a function of state which now be read as(35)

The function  is still assumed to be continuously differentiable from which one deduces the stresses  and the thermodynamical forces X and Y associated with the plastic strain and the damage. Specifically, one sets(36)

In equation (36), Y can be interpreted as the free energy release rate associated with a growth of the damage at constant strain. Since it is natural to assume that the release of energy is really nonnegative, one adds the following condition of positivity for Y:(37)condition which restricts the dependence of the free energy on α.

Since the plasticity criterion is assumed to be independent of α (hypothesis (2)), the plastic dissipation potential still reads aswhere  is a closed convex subset of . Finally, by hypothesis (3), the damage dissipated energy d reads as(38)where  is a nonnegative differentiable function of α which vanishes at α = 0:

This assumption is also motivated by the micromechanical approach of reference [12]: indeed, if one considers microcracked materials whose microcracking evolution is governed by Griffith's criterion, then dissipated energy by damage corresponds to the Griffith surface energy.

Therefore, the part of the dissipated power due to the evolution of damage simply reads in terms of the damage state and the damage rate as

whereas the total dissipated power due to the evolution of both the plastic strain and the damage is given by

It remains to formulate the damage evolution law. As for the elasto-plastic behavior, we adopt a standard law in the sense of the concept of Generalized Standard Materials proposed in reference [15] or [16] in a general context and particularized by [17,18] for damaging materials. Specifically, the evolution of the damage consists in the following three items:

1. The irreversibility condition: Damage can only grow and hence one requires that at any time.

2. The damage criterion: The free energy release rate must remain less or equal to the critical value given by the damage dissipated power:(39)

3. The consistency relation: Damage can only evolve when the free energy release rate is equal to its critical value, condition which can read as

(40)

3.2 The associative elasto-plastic Drucker–Prager model with kinematical hardening coupled with damage

Let us particularize now the model that we will use for geo materials in compression.

3.2.1 Choice of the free energy and of the dissipated energy by damage

Let us first introduce the damage variable into the expression (9) of the free energy. The micromechanical approach developed in reference [12] for a microcracked material with contact and friction between the lips of the cracks shows that only the hardening tensor A ₁ depends on the crack state, the tensor A ₀ representing the stiffness of the sound material. Therefore, we will assume that A ₁ only depends on the damage variable and hence the expression of the free energy becomes(41)

Moreover, if we assume that the material is isotropic, then one gets(42)

The positivity of the stiffness tensor A ₀ and of the hardening tensor A ₁ requires that K₀, µ ₀, K₁(α), and μ ₁(α) satisfy the following inequalities:

The energy release rate reads as(43)and hence it depends only on the plastic strain and the damage, and not on the total strain. In order that Y be nonnegative, refer to the condition (37), it is necessary and sufficient that

The damage criterion becomes then

Since the Drucher–Prager yield criterion is assumed to be independent of damage, the model will be complete once the three functions α ↦ K₁ (α), α ↦ μ ₁(α), and  are given for α ∈ [0, 1]. In fact since the choice of the damage variable is arbitrary, it is always possible, after a suitable change of variable, to fix one of the three functions. Here we choose to fix the function by setting(44) where D₁ is a positive constant which has the dimension of a stress (or equivalently an energy by volume unit). Accordingly, the dimensionless damage variable α can be interpreted as the fraction of the volume energy dissipated by damage at the current damage state by comparison to D₁ which represents the volume dissipated energy by damage when the volume element is totally damaged.

3.2.2 The complete model

Finally, the associative elasto-plastic Drucker–Prager model with kinematical hardening coupled with damage consists in the following set of definitions and relations:

The state variables:
The free energy density:
The stress–strain relationships:
The plasticity evolution law:
The plasticity thermodynamical forces:
The plasticity yield criterion:
The plasticity flow rule:
The damage evolution law:
The irreversibility condition:
The damage criterion:
The consistency condition:

3.3 Application to the triaxial test with containing pressure

3.3.1 Some general properties

Let us reconsider the triaxial test by including now the possibility of damage of the volume element. The results obtained in Section 2.4 remain valid as long as the plasticity yield criterion is not reached. In other words, the confining stage and the elastic stage remain unchanged. What happens after depends in particular on the choice of the functions μ ₁(α) and K₁(α). However we can obtain some general properties without specifying those functions provided that they satisfy the following conditions:(45) (46)

Those conditions are inspired by the micromechanical model presented in reference [12]. In particular, the fact that the hardening moduli are infinite when α = 0 ensures that the elastic moduli K₀ and µ ₀ are those of the undamaged material. In the same manner, the fact that the hardening moduli vanish when α = 1 corresponds to the total loss of stiffness when the material is fully damaged. The hypothesis that the second derivatives of K₁ and µ ₁ are positive plays an important role as we will see below. Let us now establish some general properties based on the above assumptions.

1. Damage evolves only with plasticity. Let us first show that damage can grow only when plasticity evolves. Indeed, by virtue of the damage criterion, damage evolves only when Y = D₁, i.e., when

    Differentiating this relation with respect to t leads to(47)from which one deduces that  if , which is the desired result. Let us note that this result is essentially due to the fact that the stiffness tensor A ₀ is assumed to be damage independent.

2. Damage starts at the same time as plasticity. Indeed, by virtue of the hypothesis , the damage criterion at the onset of damage givesand hence can be satisfied only if p = 0.

3. Form of the strain, plastic strain, and stress tensors during the plasticity stage with damage. The previous properties suggest to search a response such that damage and plasticity evolve simultaneously. During this damage with plasticity stage, the strain and stress tensors are assumed to be of the form³ where p ₀ remains constant and ϵ _z is a given function of time which is decreasing from its critical value given in reference (27).

4. Relationship between  and α. Inserting the form of the different tensors into the plasticity flow rule (which does not involve the damage state) gives(48)

   Reporting this relation into the damage criterion gives  in function of α:(49) where R (α) is the following combination of the hardening moduli:(50)

   Let us note that, since the second derivatives of K₁ and µ ₁ are positive,  is a strictly increasing function of α. That means that the required inequality  is satisfied provided that .

5. Relationship between σ _z and α. The plasticity criterion with the definition of the plasticity thermodynamical forces give the following three equations:

   Eliminating X_T , X_z , X_m and using equation (49) allow us to obtain σ _z as the following function of α:(51)where

   One then checks that σ _z is always negative. Its absolute value will be either increasing or decreasing when α grows according to the second derivative of the compliance function S (α) which is negative or positive. The former case, when , corresponds to a stress-hardening response, whereas the second one, when , corresponds to a stress-softening response. Note that the sign of the second derivative of S(α) is not given by the sign of the second derivative of R(α) and hence constitutes an additional choice for the model.

6. Relationships between ϵ _z , ϵ _v , and α. The stress–strain relations allow us to express ϵ _z and ϵ _v in terms of α. Specifically, taking into account equations (49) and (51), one gets(52) (53)

7. The relation between the axial stress σ _z and the axial strain ϵ _z . The relations (51) and (52) give the relation between ϵ _z and σ _z under the form of a curve parametrized by α. Let us note that we are not ensured that ϵ _z is a monotonic function of α. Specifically, when  (stress-hardening behavior), then  when  and hence damage grows when ϵ _z (which is negative) is decreasing. But, when , then the sign of is not guaranteed when , because the term in σ _z is increasing, whereas the term in  is decreasing. Therefore, a snap-back is possible. If such a snap-back exists, the control of the axial strain ϵ _z necessarily implies a discontinuous evolution of the damage. That can even lead to a brutal rupture of the volume element with a sudden jump of α to 1 when ϵ _z reaches its limit value.

8. The relation between the volumetric strain ϵ _v and the axial strain ϵ _z . The relations (51)–(53) give the relation between ϵ _v and σ _z under the form of a curve parametrized by α. Here also one can have a competition between the term in σ _z and the term in . During the hardening phases, one has  and . Consequently, according to the respective weight of the two terms, one will observe either a contractance or a dilatancy. On the contrary, during the softening phases, since  and , one will always observe a dilatancy.

   These different properties are illustrated in the next section by particularizing the model. Let us note that in the triaxial test with a fixed confining pressure, the functions μ ₁(α) and K₁(α) appear only by their com bination R(α) (which involves also the internal friction coefficient k). Therefore, it is sufficient to know the function R(α) and the constants µ ₀, K₀, D₁, k, and τ _c to characterize the model, the confining pressure p ₀ playing the role of a parameter.

3.3.2 Study of a family of models

Let us consider the following family of functions R(α):(54)

This choice of function is governed mainly by our intention to control the asymptotic behavior both at the initiation of damage (α ≈ 0 power law of n) and at the end of microcracking (α ≈ 1 power law of m). One easily checks that such a function α ↦ R(α) is really compatible with the hypotheses (45) and (46). Indeed, R (α) is strictly decreasing from +∞ to 0 when α grows from 0 to 1. Its first and second derivatives are given by

It is also not difficult to verify that  for every α ∈ (0, 1) because m > 1 > n > 0. The first derivative R^′(α) is strictly increasing from −∞ to 0 when α grows from 0 to 1. By virtue of equation (49), one deduces that  grows from 0 to +∞ with the evolution of damage.

The compliance function S (α) and its first and second derivatives are given by See equation above.

The study of the sign of S^′′(α) shows that S^′′(α) is negative in the interval (0, α ₀), positive in the interval (α ₀, 1) with α ₀ given by(55)

Therefore, by virtue of equation (51), the absolute value of the axial stress σ _z increases and one observes a stress-hardening behavior when α grows from 0 to α ₀, then |σ _z | decreases and one observes a stress-softening behavior when α grows from α ₀ to 1. Moreover, since S^′ (0) =  S^′(1) =+ ∞, σ _z starts from the value given by equation (26) and corresponding to the end of the elastic stage, then increases (in absolute value) until its maximal value corresponding to the time at which the damage reaches the value α ₀, and finally decreasing (in absolute value) to come back to the value (26) reached at the end of the elastic stage. The absolute value of the overstress Δσ _z associated with the stress-hardening phase is given by(56)

Let us first study the variations of the axial strain ϵ _z in function of the damage α and let us analyze under which condition a snap-back occurs. Starting from equation (52) and using the fact that ϵ _z is negative, one deduces that the derivative of ϵ _z with respect to α is always negative. Therefore, a snap-back exists if and only if the following condition is satisfied:(57)

Since R (α) is of the form , the left-hand side above can read as(58)where the function α ↦ φ (α) depends only on the exponents m and n of the model. In the neighborhood of α = 1, the different terms behave as follows:

Hence, one gets

Therefore, since φ is continuous, negative in the interval (0, α ₀), positive in the interval (α ₀, 1), and since φ (α ₀) = φ (1) =0, the function φ reaches its upper bound in the interval (α ₀, 1), the upper bound depending only on m and n. Hence, the condition of non-snap-back (57) can read as(59)

In other words, there is no snap-back provided that the hardening modulus R₁ is small enough. Let us note that this condition depends on the parameters E₀, K, m, and n of the model, but is independent of D₁, τ _c , and the confining pressure p ₀.

Let us now study the variations of the volumetric strain ϵ _v as a function of the damage α and let us analyze when one observes a contractance or a dilatancy. We assume that R₁ is small enough and satisfies equation (59) so that there is no snap-back in the response σ _z  − ϵ _z and hence . By virtue of equation (53), one has  and hence con tractance if and only if

With the help of equations (49) and (51), that condition can be expressed in terms of φ (α) defined by equation (58). Specifically, one observes a contractance at the time when the damage reaches the value α if the following condition is satisfied:(60)

That can happen only during the phase where φ is negative, that is only during the stress-hardening phases. During the stress-softening phases, since φ > 0, one observes a dilatancy. Therefore, a contractance can happen only at the beginning of the test, as long as α < α ₀. Since the behavior of the first and the second derivatives of R and S in the neighborhood of α = 0 is given byone gets

Consequently, since φ (0) =−∞, one obtains a contractance when α is small enough. On the contrary, as soon as α becomes larger than a critical value (which is necessarily less than α ₀), one observes a dilatancy until the complete rupture of the volume element. Let us note that since tends to 0 when α tends to 1, one gets by virtue of equations (52) and (53)

when α is close to 1. That is in conformity with the result (29) obtained in the case of a perfect plasticity behavior without kinematical hardening, situation to which one tends when α becomes close to 1.

Let us finally study the relation between the axial strain ϵ _z and the axial stress σ _z in the case when the condition (59) of non-snap-back is satisfied. Using the behavior of R^′ and S^′ in the neighborhood of α = 0, one deduces from equations (49) and (51) the following behaviors for  and σ _z :

Thus, the plastic strain is on the order of α ^(1+n)/2 while the variation of the elastic strain is on the order of α ^(1−n)/2. Hence, at the beginning of the plasticity stage, the plastic strain is negligible by comparison with the variation of the elastic strain. Inserting these estimates into equation (52), one deduces that the relation between ϵ _z and σ _z is continuous and continuously differentiable at the transition between the elastic stage and the plasticity stage, the slope being equal to E₀. After the elastic stage, if the axial strain |ε _z | is increased until +∞, then |σ _z | begins by increase, and then passes by a maximum before to decrease to finally take the value that it had at the end of the elastic stage. All those variations are in fact independent of the confining pressure.

3.3.3 Graphs of the response according to the values of the parameters m, n, and R₁

In all the pictures presented in this section, one uses the following values for the Poisson ratio ν ₀ and the parameters τ _c and k of the Drucker–Prager criterion:

One sets(61)

The material constant σ _c fixes the scale of the stresses, whereas the material constant ϵ _c fixes the scale of the strains. The responses depend on the two exponents m and n, and on the two ratios R₁/E₀ and p ₀/σ _c .

• Case where m = 2, n = 1/2. The value α ₀ of the damage from which there is stress-softening is given by equation (55) and hence one finds α ₀ = 0.0883. By virtue of equation (56), the overstress due to the stress-hardening phase depends on R₁/E₀ (but not on the confining pressure p ₀) and is given by(62)

  In order that there is no snap-back in the response σ _z  − ε _z , it is necessary that the condition (59) be satisfied. For the considered values of m and n, one gets  and hence the non-snap-back condition requires that R₁ ≤ 0.275E₀.

  Figure 1 shows the graphs of σ _z and ϵ _v as functions of ϵ _z when there is no confining pressure (p ₀ = 0) and in the case where R₁ = 0.1E₀ (hence there is no snap-back). In that case the value of the overstress is .

  In Figure 2 are plotted the same graphs for different values of R₁/E₀. One can note that the larger the R₁, the larger the overstress, which conforms to equation (62), but also the more rapid is the decrease of the axial stress after the peak. A similar behavior can be seen for the volumetric strain: the larger the R₁, the greater the maximal contraction, but also the more rapid the growing of the dilatancy after the peak.

  When R₁ is greater than 0.275E₀, the graph of σ _z contains a snap-back which induces, during a test where the axial strain is controlled and continuously decreasing, a discontinuity of the evolution of the axial stress and the damage (cf. Fig. 3).

  The graph of the axial stress vs. the axial strain for a given confining pressure during the damage-plasticity phase is simply a translation of that corresponding to the uniaxial compression test, the translation being given by the values of the axial strain and the axial stress at the end of the elastic stage (cf. Fig. 4). This property is due to the linear character of Drucker–Prager criterion. That would not be true anymore if one used a more general plasticity criterion like Hoek–Brown criterion [1,2].

• Influence of the exponents m and n. The shape of the graphs giving the axial stress in terms of the axial strain (during an uniaxial compression test without confining pressure) is always the same whatever the values of m and n provided that they remain inside the admissible intervals, i.e., m > 1 > n > 0. One can see the dependence of the response on m an n in Figure 5.

  We note that as expected the m parameter influence mainly the post-pic stress behavior (for m = 1 the ultimate damage state achieved even for the finite value of ϵ _z ), while n parameter alternate mainly the initial hardening state leaving the large |ε _z | asymptote intact. Both parameters (m and n) contribute (Eqs. (55) and (56)) to the overstress Δσ _z estimation that reads as

  where ϕ (m, n) is only dependent on m and n. At given n, ϕ (m, n) is a decreasing function of m, whereas at given m, ϕ (m, n) is first a decreasing function, then an increasing function of n when n varies from 0 to 1 (cf. Fig. 6).

• Possible parameter fitting strategy. Even so our model depends just on three “damage” parameters (R₁, n, m), the common task of parameter choice for the given material could be a rather laborious exercise. Nevertheless, based on the obtained earlier parameter influence properties, one could develop straight experimental data fitting strategy. From the experimental data set characterizing geomaterials, one has usually access to dilatancy () and softening curves (σ _z  (loading)). The softening post-pic behavior would allow us to fit first the m value. Once the m parameter is fixed, one could proceed to R₁ value choice which allow adjusting the strength of dilatancy curve. Finally, the experimental overstress threshold could be fixed by adapting the last n parameter.

4 Conclusion and perspectives

Starting from low-scale mechanical properties of microcracked elastic materials and taking into account the closed cracks friction between its lips, we propose a phenomenological model with an associative elasto-plastic behavior coupled with damage. The main feature of the model which is inspired by the micromechanical considerations is that the free energy contains not only the usual elastic energy but also a stored energy which is due at the microlevel to the friction between the lips of the microcracks. This blocked energy is assumed to depend at the macro-level on the damage state and the plastic strain and induces a kinematical hardening in the plasticity yield criterion. Considering a Drucker–Prager yield criterion and adopting the normality rule for the plasticity flow rule, we have applied such a model to a triaxial test with a fixed confining pressure. It turns out that the model is able to account for the main observed properties of geomaterials in such a situation, like contractance or dilatancy effect on the evolution of volumetric strain and stress-hardening or stress-softening effect for the axial stress.

The next step will be to identify from experimental results the parameters of the proposed model (which can be enriched by considering more complex yield criterion than the Drucker–Prager one). Its validity will be also tested by comparing its predictions with other families of experimental tests like the oedometric test or cyclic tests. An interesting issue is also to use such a model to calculate the response of a true three-dimensional sample (and no more of the volume element). Since the damage evolution is accompanied by stress-softening effect, one can expect that the response of the sample is no more homogeneous (in space). That could require to enhance the damage model by introducing non-local terms like in references [23,24].

Fig. 1 Graphs of the axial stress σ z (left) and of the volumetric strain ϵ v (right) versus the axial strain ϵ z when R1 = 0.1E0 and p 0 = 0 (uniaxial compression test).

Fig. 2 Influence of R1/E0 on the graphs of the axial stress and the volumetric strain versus the axial strain when there is no confining pressure (uniaxial compression test). The value of R1/E0 associated with each graph is indicated in red on the graph.

Fig. 3 When R1/E0 is large enough, the graph of σ z vs. ϵ z contains a snap-back (here R1 = 0.4E0). Consequently, when the axial strain reaches the value corresponding to the limit point A, the response jumps from A to B.

Fig. 4 Graph of the axial stress vs. the axial strain for different values (in red on the graphs) of the confining pressure (in fact of the ratio p 0/σ c ) when R1 = 0.1E0. The first part (the linear part) of the graphs correspond to the hydrostatic and elastic stages.

Fig. 5 Influence of the exponent m (left) and n (right) on the graph of the axial stress vs. the axial strain in a uniaxial compression test (p 0 = 0). The values of m and n are indicated in red on the graph. In the left figure where m varies, one takes n = 1/2, whereas on the right figure where n varies, one takes m = 2. In all cases R1 = 0.1E0.

Fig. 6 Dependence of the overstress Δσ z on the exponents m and n of the model. Left: dependence on m when n = 1/2. Right: dependence on n when m = 2. In all cases R1 = 0.1E0 and k = 0.2.

References

All Tables

All Figures

	Fig. 1 Graphs of the axial stress σ z (left) and of the volumetric strain ϵ v (right) versus the axial strain ϵ z when R1 = 0.1E0 and p 0 = 0 (uniaxial compression test).
In the text

	Fig. 2 Influence of R1/E0 on the graphs of the axial stress and the volumetric strain versus the axial strain when there is no confining pressure (uniaxial compression test). The value of R1/E0 associated with each graph is indicated in red on the graph.
In the text

	Fig. 3 When R1/E0 is large enough, the graph of σ z vs. ϵ z contains a snap-back (here R1 = 0.4E0). Consequently, when the axial strain reaches the value corresponding to the limit point A, the response jumps from A to B.
In the text

	Fig. 4 Graph of the axial stress vs. the axial strain for different values (in red on the graphs) of the confining pressure (in fact of the ratio p 0/σ c ) when R1 = 0.1E0. The first part (the linear part) of the graphs correspond to the hydrostatic and elastic stages.
In the text

	Fig. 5 Influence of the exponent m (left) and n (right) on the graph of the axial stress vs. the axial strain in a uniaxial compression test (p 0 = 0). The values of m and n are indicated in red on the graph. In the left figure where m varies, one takes n = 1/2, whereas on the right figure where n varies, one takes m = 2. In all cases R1 = 0.1E0.
In the text

	Fig. 6 Dependence of the overstress Δσ z on the exponents m and n of the model. Left: dependence on m when n = 1/2. Right: dependence on n when m = 2. In all cases R1 = 0.1E0 and k = 0.2.
In the text