Jean-Jacques Marigo

Kyrylo Kazymyrenko

A micromechanical inspired model for the coupled to damage elasto-plastic behavior of geomaterials under compression

Keywords: Plasticity, Damage, compression test, geomaterials, dilatancy

We propose an elasto-plastic model coupled with damage for the behavior of geomaterials in compression. The model is based on the properties, shown in [2], of micro-cracked materials when the micro-cracks are closed with a friction between their lips. That leads to a macroscopic model coupling damage and plasticity where the plasticity yield criterion is of the Drucker-Prager type with kinematical hardening. Adopting an associative flow rule for the plasticity and a standard energetic criterion for damage, the properties of such a model are illustrated in a triaxial test with a fixed confining pressure.

Introduction

Many works have been devoted to the modeling of the mechanical behavior of geo-materials 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 [START_REF] Hoek | The Hoek-Brown failure criterion-a 1988 update[END_REF][START_REF] Eberhardt | The Hoek-Brown failure criterion[END_REF][START_REF] Cuvilliez | An elastoviscoplastic constitutive model for geomaterials: Application to hydromechanical modelling of claystone response to drift excavation[END_REF], whereas it is generally admitted that one cannot reproduce the experimental results observed for concrete without considering a damage evolution law [START_REF] Lemaître | Coupled elasto-plasticity and damage constitutive equations[END_REF][START_REF] Lemaître | Mécanique des Matériaux Solides[END_REF][START_REF] Mazars | A model of a unilateral elastic damageable material and its application to concrete[END_REF][START_REF] Ragueneau | Damage model for concrete-like materials coupling cracking and friction, contribution towards structural damping: first uniaxial applications[END_REF][START_REF] Comi | Fracture energy based bi-dissipative damage model for concrete[END_REF][START_REF] Pontiroli | Predicting concrete behaviour from quasi-static loading to hypervelocity impact An overview of the PRM model[END_REF][START_REF] Richard | Isotropic continuum damage mechanics for concrete under cyclic loading: Stiffness recovery, inelastic strains and frictional sliding[END_REF]. 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 [START_REF] Hoek | The Hoek-Brown failure criterion-a 1988 update[END_REF][START_REF] Eberhardt | The Hoek-Brown failure criterion[END_REF] or the cap model [START_REF] Babadjian | Quasi-static evolution in nonassociative plasticity: the cap model[END_REF]), with common feature being the hypothesis of a non associative 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 of or all 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 to 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 micro-mechanical approach. An exception is [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF] 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 [START_REF] Zhu | A micromechanics-based thermodynamic formulation of isotropic damage with unilateral and friction effects[END_REF] where a full 3D approach is proposed. Even if the analysis in [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF] 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 ε, the plastic strain p:

ψ(ε, p) = 1 2 A 0 (ε -p) • (ε -p) + 1 2 A 1 p • p. (1)
In the expression above one uses the notations precised at the end of the present section. In particular, A 0 denotes the stiffness tensor of the sound material whereas A 1 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 A 1 dependence) and the damage state, not the total strain. The damage state is supposed not to enter in the first term of the energy expression [START_REF] Alessi | Gradient damage models coupled with plasticity: variational formulation and main properties[END_REF] and consequently, even if the stress-strain relation still reads as

σ = A 0 (ε -p),
the stress which appears in the plasticity yield criterion deduced from Coulomb friction law is not σ but the tensor X given by X = σ -A 1 p.

In other words, the yield criterion is of the Drucker-Prager type with a kinematical hardening where the back stress A 1 p depends linearly on the plastic strain but also non linearly 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 [START_REF] Richard | Isotropic continuum damage mechanics for concrete under cyclic loading: Stiffness recovery, inelastic strains and frictional sliding[END_REF][START_REF] Lanoye | An isotropic unilateral damage model coupled with frictional sliding for quasi-brittle materials[END_REF] 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 consider 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 micro-mechanical considerations. We finally study the response of a volume element submitted to a triaxial test at fixed a confining pressure. The damage is introduced in Section 3, still on the basis of the micro-mechanical 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 Table 1): The summation convention on repeated indices is implicitly adopted. The vectors and second order tensors are indicated by boldface letters, like n 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 A 0 or A 0ijkl 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 A 0 ε whose ij-component is A 0ijkl ε kl . The inner product between two vectors or two tensors of the same order is indicated by a dot, like a • b which stands for 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.

a i b i or σ • ε for σ ij ε ij . ε tensor

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 p without other internal variables (the damage variable will be introduced later, in Section 3). Therefore the free energy density ψ is given by

ψ = ψ(ε, p), (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 X of the thermodynamical forces associated with the plastic strain tensor:

σ = ∂ ψ ∂ε (ε, p), X = - ∂ ψ ∂p (ε, p). (3)
From the mathematical point of view, ε, p, σ and X are order 2 symmetrical tensors which can be identified with 3 by 3 symmetrical matrices after the choice of an orthonormal basis. In other words, those tensors will be considered as elements of M 3 s . The plasticity (or yield) criterion is defined by giving the (closed with non empty interior) convex set K of M 3 s where the thermodynamical forces X must lie. Here we make the strong assumption that this domain is fix, that is independent of time and of the plasticity evolution. (However, we will see in the next subsection that the image of this domain in the stress space evolves with the plastic strain.) In practise the convex K is characterized by a convex function f : M 3 s → R so that the plasticity criterion reads as

X ∈ K = {X * ∈ M 3 s : f(X * ) ≤ 0}. ( 4
)
The interior points of K correspond to the forces X * such that f(X * ) < 0 and the points on the boundary of K to the forces X * such that f(X * ) = 0. Since we only consider an associative model, the evolution of the plastic strain follows the normality rule. When the convex set K has a smooth boundary without angular points, the function f is differentiable and the flow rule can read as

ṗ = η ∂f ∂X (X) with η = 0 if f(X) < 0 ≥ 0 if f(X) = 0 , (5)
where η denotes the plastic multiplier. If the boundary of K 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

X ∈ K, (X -X * ) • ṗ ≥ 0, ∀X * ∈ K. (6)
The inequality [START_REF] Cuvilliez | An elastoviscoplastic constitutive model for geomaterials: Application to hydromechanical modelling of claystone response to drift excavation[END_REF] says that ṗ must belong to the cone of the outer normals to K at the point X of the boundary. That inequality leads to (5) when the normal is well defined at X. 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 K. Specifically, let us introduce the support function π K of K:

p ∈ M 3 s → π K (p) = sup X * ∈K X * • p ∈ R. ( 7
)
By construction π K is convex and positively homogeneous of degree 1,

π K (λp) = λπ K (p), ∀λ > 0, ∀p ∈ M 3 s .
Moreover π K vanishes at p = 0, is non negative provided that 0 ∈ K and can take the value +∞ when K 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) D is defined by

D = σ • ε - ψ
and the second principle of thermodynamics requires that D be non negative. Owing to (2)-( 3), D reads as

D = X • ṗ.
By virtue of Hill maximal work principle [START_REF] Cuvilliez | An elastoviscoplastic constitutive model for geomaterials: Application to hydromechanical modelling of claystone response to drift excavation[END_REF] and the definition (7) of the support function, one gets

D = π K ( ṗ), (8)
which means that π K plays also the role of the dissipation potential. Furthermore, it is sufficient to have 0 ∈ K for Clausius-Duhem inequality to be satisfied, i.e. D ≥ 0.

Choice of the form of the free energy state function

Guided by the micromechanical considerations presented in [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF] 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 Tr p = 0 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 p 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 (ε, p):

ψ(ε, p) = 1 2 A 0 (ε -p) • (ε -p) + 1 2 A 1 p • p. (9)
In [START_REF] Francfort | The elasto-plastic exquisite corpse: A Suquet legacy[END_REF] the first term represents the elastic energy and A 0 denotes the (fourth order) stiffness tensor of the sound material. The second term in [START_REF] Francfort | The elasto-plastic exquisite corpse: A Suquet legacy[END_REF] represents the blocked energy by friction and the (fourth order) stiffness tensor A 1 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 ṗ = 0.

The tensor A 0 is positive definite and the tensor A 1 is assumed to be non negative. In fact, A 1 is also positive definite as soon as the material is microcracked as it is shown in [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF]. But, we will also consider the case when A 1 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 (ε, p).

To simplify the presentation, we will assume that the material remains isotropic (even when the damage grows). Therefore both tensors A 0 and A 1 have only two independent moduli. Decomposing the strain and the plastic strain tensors into their spherical and deviatoric parts

ε = 1 3 Tr ε I + ε D p = 1 3 Tr p I + p D , (10)
the free energy finally reads as

ψ(ε, p) = 1 2 K 0 (Tr ε -Tr p) 2 + µ 0 (ε D -p D ) • (ε D -p D ) + 1 2 K 1 (Tr p) 2 + µ 1 p D • p D , (11)
where K 0 > 0 and µ 0 > 0 are the compressibility and shear moduli of the sound material whereas K 1 ≥ 0 and µ 1 ≥ 0 will be called the kinematical hardening moduli. The elastic moduli K 0 and µ 0 are related to the Young modulus E 0 and the Poisson ratio ν 0 of the sound material by:

3K 0 = E 0 1 -2ν 0 , 2µ 0 = E 0 1 + ν 0 . (12)
By differentiation of ( 11), the stress-strain relation is linear and reads as

σ = A 0 (ε -p). (13)
Decomposing the stress tensor into its spherical and deviatoric parts,

σ = σ m I + σ D , σ m = 1 3 Tr σ, (14)
the relation [START_REF] Lanoye | An isotropic unilateral damage model coupled with frictional sliding for quasi-brittle materials[END_REF] becomes

σ D = 2µ 0 (ε D -p D ) σ m = K 0 (Tr ε -Tr p) . ( 15
)
Still by differentiation, the thermodynamical forces associated with the plastic strains read as

X = σ -A 1 p (16)
and they differ in general from the stresses when A 1 = 0. Decomposing the forces X into their spherical and deviatoric parts,

X = X m I + X D , X m = 1 3 Tr X, (17)
the relation ( 16) becomes

X D = σ D -2µ 1 p D X m = σ m -K 1 Tr p . ( 18
)

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 geo-materials under compression, see [START_REF] Lemaître | Mécanique des Matériaux Solides[END_REF][START_REF] Hoek | The Hoek-Brown failure criterion-a 1988 update[END_REF][START_REF] Eberhardt | The Hoek-Brown failure criterion[END_REF]. In general those criteria are used in a non standard context, i.e. with a non associative flow rule for the plasticity evolution, like in [START_REF] Cuvilliez | An elastoviscoplastic constitutive model for geomaterials: Application to hydromechanical modelling of claystone response to drift excavation[END_REF]. The main usual argument which is evoked to use a non associative flaw rule is that one overestimates the dilatancy effect with the normality rule associated with Drucker-Prager criterion. On the contrary, we propose in the present 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

f(X) := 1 √ 6 X D + kX m -τ c ≤ 0 . (19)
In [START_REF] Mielke | Evolution of rate-independent systems[END_REF] 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 √ 6 is introduced for convenience in order to simplify forthcoming expressions),

X D = √ X D • X D .
Consequently, the elastic domain

K = {X ∈ M 3 s : f(X) ≤ 0}
is really a closed convex subset of M 3 s (with non empty interior). In fact, K is a convex cone which has an angular point at (X D , X m ) = (0, τ c /k) and which is unbounded in the direction of hydrostatic compressions. Furthermore, even if the elastic domain is fix in the space of the thermodynamical forces, it varies with the plastic strain in the stress space. Indeed, by virtue of ( 16), the elastic domain expressed in term of σ becomes the set K(p) given by

K(p) = σ ∈ M 3 s : 1 √ 6 σ D -2µ 1 p D + k(σ m -K 1 Tr p) ≤ τ c .
Thus K(p) is subjected to a linear translation when the plastic strain evolves:

K(p) = A 1 p + K.
Starting from [START_REF] Eberhardt | The Hoek-Brown failure criterion[END_REF], the support function π K of the convex K is obtained after tedious calculations which are not reproduced here and finally reads as

π K (p) =    τ c k Tr p if Tr p ≥ k √ 6 p D +∞ otherwise . ( 20
)
Since we adopt the normality rule, the evolution of the plastic strain is directly deduced from [START_REF] Mielke | Evolution of rate-independent systems[END_REF]. Considering first a regular point of the boundary of K, the flow rule reads as follows

At X such that f(X) = 0, X = (0, τ c /k),    ṗD = Tr ṗ k √ 6 X D X D Tr ṗ ≥ 0 . ( 21
)
Let us note that at such a point one has Tr ṗ = k √ 6 ṗD . At the angular point, the rate of the plastic strain tensor must belong to the cone of the outer normal to K. That leads to

At X = (0, τ c /k), Tr ṗ ≥ k √ 6 ṗD . (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.

Response under a triaxial test with a confining pressure

The elasto-plastic behavior of a material governed by the associative Drucker-Prager presented in the previous subsections 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 an hydrostatic compression where the pressure is progressively increased from 0 to a final value p 0 > 0; 2. Then, maintaining the lateral pressure to the value p 0 , one compresses in the axial direction z by prescribing the axial strain ε z with εz < 0. 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, when |ε 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 in respect of |ε z | (which could be considered as the loading parameter)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

σ = -p 0 I, ε = - p 0 3K 0 I, p = 0. (23)
Since p = 0, one has X = σ and Ducker-Prager criterion gives

f(X) = -kp 0 -τ c < 0.
The axial strain ε z decreases from 0 to -p 0 /3K 0 whereas the volumetric strain ε v = Tr ε decreases from 0 to -p 0 /K 0 . During all this stage, one has εv = 3 εz < 0 and hence a contractance due to the hydrostatic compression.

Elastic stage at fixed confining pressure

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

p = 0.
By symmetry, the stress and the strain tensors are of the following form2 :

X = σ =   -p 0 0 0 0 -p 0 0 0 0 σ z   , ε =   ε T 0 0 0 ε T 0 0 0 ε z   ,
where p 0 and ε z are prescribed. The deviatoric parts read as

σ D = (σ z + p 0 )J, ε D = (ε z -ε T )J with J =   -1/3 0 0 0 -1/3 0 0 0 2/3   , J = 2 3 .
It remains to determine σ z and ε T . They are given by ( 13) with the help of ( 12):

ε T = -ν 0 ε z -(1 -2ν 0 )(1 + ν 0 ) p 0 E 0 , σ z = E 0 ε z -2ν 0 p 0 . ( 24
)
One deduces in particular the evolution of the volumetric strain:

εv = (1 -2ν 0 ) εz < 0 . ( 25
)
This stage stops when the stresses reach the yield criterion, i.e. when f(σ) = 0. Since

f(σ) = 1 3 |σ z + p 0 | + k 3 (σ z -2p 0 ) -τ c
and since σ z + p 0 < 0, the yield surface is reached when

σ z = - 1 + 2k 1 -k p 0 - 3τ c 1 -k . ( 26
)
The corresponding value of ε z is deduced from [START_REF] Zhu | A micromechanics-based thermodynamic formulation of isotropic damage with unilateral and friction effects[END_REF]:

ε z = - 1 -2ν 0 + 2(1 + ν 0 )k 1 -k p 0 E 0 - 3τ c (1 -k)E 0 . (27)

Plastification stage

Then, if the confining pressure is maintaining constant and if one prescribes εz < 0, 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 µ 1 and K 1 vanish, i.e. the case without hardening.

1. Case without kinematical hardening. In such a case, one gets X = σ and since the elastic domain is fix 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 0 and the axial stress σ z is given by (26). Since σ = 0, one gets ṗ = ε and hence p takes the form

p =   p T 0 0 0 p T 0 0 0 p z   , p D = (p z -p T )J, Tr p = p z + 2p T . (28)
Since σ z + p 0 < 0, the flow rule [START_REF] Pontiroli | Predicting concrete behaviour from quasi-static loading to hypervelocity impact An overview of the PRM model[END_REF] gives

ṗT -ṗz = 1 2k ( ṗz + 2 ṗT ) ≥ 0.
Since ṗz = εz et ṗT = εT , one gets

εT = ṗT = - 1 + 2k 2(1 -k) εz > 0.
Therefore, the volumetric strain rate is given by εv

= - 3k 1 -k εz > 0 . ( 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 non associative 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.

f(X) = 0 with X z < X T .
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 non associative flow rule.)

Accordingly, one gets

X D = (X z -X T )J, 3X m = X z + 2X T and X T -X z + 3kX m = 3τ c .
Differentiating these relations with respect to time and taking into account (18) lead to

     ẊT -Ẋz + 3k Ẋm = 0 Ẋz -ẊT = σz -2µ 1 ( ṗz -ṗT ) 3 Ẋm = σz -3K 1 Tr ṗ
The flow rule [START_REF] Pontiroli | Predicting concrete behaviour from quasi-static loading to hypervelocity impact An overview of the PRM model[END_REF] gives

ṗz -ṗT = - Tr ṗ 2k . ( 30
)
Combining the four previous relations gives the following relationship between σz and Tr ṗ:

σz = - µ 1 + 3k 2 K 1 (1 -k)k Tr ṗ . ( 31
)
Differentiating the stress-strain relations [START_REF] Lanoye | An isotropic unilateral damage model coupled with frictional sliding for quasi-brittle materials[END_REF] gives

     σz = 2µ 0 ( εz -εT + ṗT -ṗz ) σz = 3K 0 ( εv -Tr ṗ) εv = εz + 2 εT
.

Eliminating εT and using the flow rule give two other relations between σz , Tr ṗ, εv and εz :

σz = µ 0 3 εz -εv + Tr ṗ k , σz = 3K 0 ( εv -Tr ṗ). ( 32
)
Using (31) allows us to determine σz , Tr ṗ and εv in terms of εz . We finally get

Tr ṗ = - 3k εz 1 -k + 3(µ 1 + 3k 2 K 1 ) (1 -k)E 0 > 0 , (33) εv
= 1 - µ 1 + 3k 2 K 1 k(1 -k)3K 0 Tr ṗ . (34)
Therefore (33) shows that the kinematical hardening tends to decrease Tr ṗ and even (34) shows that the sign of εv depends on the intensity of the hardening. Specifically, if the hardening moduli are small enough (by comparison with the elastic moduli), then εv > 0 and there is dilatancy. But, if the hardening moduli are large enough, then εv < 0 and there is contractance. Moreover let us note that if µ 1 = K 1 = 0, then we recover that εv = Tr ṗ and the relation (29). But, if µ 1 and K 1 are very large by comparison to µ 0 and K 0 , then one gets εv ≈ (1 -2ν 0 ) εz < 0, 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 details in the next section) for a relevant damage coupling -once the damage affects A 1 tensor it then directly impacts the dilatancy strength of the given law.

The elasto-plastic behavior coupled with damage

In this section we introduce the damage into the elasto-plastic model. As mentioned above 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 particularize the model.

Main assumptions on the damage dependence

Although there exists an infinite number of manners to couple damage with plasticity [START_REF] Lemaître | Mécanique des Matériaux Solides[END_REF][START_REF] Mazars | A model of a unilateral elastic damageable material and its application to concrete[END_REF][START_REF] Pontiroli | Predicting concrete behaviour from quasi-static loading to hypervelocity impact An overview of the PRM model[END_REF][START_REF] Ragueneau | Damage model for concrete-like materials coupling cracking and friction, contribution towards structural damping: first uniaxial applications[END_REF][START_REF] Comi | Fracture energy based bi-dissipative damage model for concrete[END_REF], we will make the simplest choices. Specifically, the main assumptions are the following ones:

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 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 K is independent of α.

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 (ε, p, α) and the free energy is still a function of state which now be read as

ψ = ψ(ε, p, α) (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

σ = ∂ ψ ∂ε (ε, p, α), X = - ∂ ψ ∂p (ε, p, α), Y = - ∂ ψ ∂α (ε, p, α). (36)
In (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 non negative, one adds the following condition of positivity for Y :

- ∂ ψ ∂α (ε, p, α) ≥ 0, ∀ε, ∀p ∈ M 3 s , ∀α ∈ [0, 1], (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 read as

π K ( ṗ) = sup X * ∈K X * • ṗ,
where K is a closed convex subset of M 3 s . Finally, by hypothesis (3), the damage dissipated energy d reads as

d = D(α), ( 38
)
where D is a non negative differentiable function of α which vanishes at α = 0:

D(0) = 0, D (α) > 0, ∀α ∈ [0, 1].
This assumption is also motivated by the micro-mechanical approach of [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF]: indeed, if one considers microcracked materials whose microcraking 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 ḋ = D (α) α, whereas the total dissipated power due to the evolution of both the plastic strain and the damage is given by

D = π K ( ṗ) + D (α) α.
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 by [START_REF] Halphen | Sur les matériaux standards généralisés[END_REF] or [START_REF] Germain | Continuum Thermodynamics[END_REF] in a general context and particularized by [START_REF] Marigo | Formulation d'une loi d'endommagement d'un matériau élastique[END_REF][START_REF] Marigo | Constitutive relations in Plasticity, Damage and Fracture Mechanics based on a Work Property[END_REF] 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 α ≥ 0 at any time; 2. The damage criterion: the free energy release rate must remain less or equal to the critical value D (α) given by the damage dissipated power:

Y ≤ D (α); (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 (Y -D (α)) α = 0. (40)

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.

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 micro-mechanical approach developed by [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF] for a micro-cracked material with contact and friction between the lips of the cracks shows that only the hardening tensor A 1 depends on the crack state, the tensor A 0 representing the stiffness of the sound material. Therefore we will assume that A 1 only depends on the damage variable and hence the expression of the free energy becomes ψ(ε, p, α)

= 1 2 A 0 (ε -p) • (ε -p) + 1 2 A 1 (α)p • p. (41)
Moreover, if we assume that the material is isotropic, then one gets

ψ(ε, p, α) = 1 2 K 0 (Tr ε -Tr p) 2 + µ 0 (ε D -p D ) • (ε D -p D ) + 1 2 K 1 (α)(Tr p) 2 + µ 1 (α)p D • p D . ( 42
)
The positivity of the stiffness tensor A 0 and of the hardening tensor A 1 requires that K 0 , µ 0 , K 1 (α) and µ 1 (α) satisfy the following inequalities

K 0 > 0, µ 0 > 0, K 1 (α) ≥ 0, µ 1 (α) ≥ 0, ∀α ∈ [0, 1].
The energy release rate reads as

Y = -1 2 K 1 (α)(Tr p) 2 -µ 1 (α)p D • p D (43)
and hence it depends only on the plastic strain and the damage, not on the total strain. In order that Y be non negative, cf the condition (37), it is necessary and sufficient that

K 1 (α) ≤ 0, µ 1 (α) ≤ 0, ∀α ∈ [0, 1].
The damage criterion becomes then

-1 2 K 1 (α)(Tr p) 2 -µ 1 (α)p D • p D ≤ D (α).
Since the Drucher-Prager yield criterion is assumed to be independent of damage, the model will be complete once the three functions α → K 1 (α), α → µ 1 (α) and α → D(α) 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 D by setting

D(α) = D 1 α (44)
where D 1 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 1 which represents the volume dissipated energy by damage when the volume element is totally damaged.

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 :

(ε, p, α) ∈ M 3 s × M 3 s × [0, 1]
The free energy density :

ψ = 1 2 K 0 (Tr ε -Tr p) 2 + µ 0 (ε D -p D ) • (ε D -p D ) + 1 2 K 1 (α)(Tr p) 2 + µ 1 (α)p D • p D
The stress-strain relationships:

σ D = 2µ 0 (ε D -p D ) σ m = K 0 (Tr ε -Tr p)
The plasticity evolution law:

The plasticity thermodynamical forces:

X D = σ D -2µ 1 (α)p D X m = σ m -K 1 (α) Tr p
The plasticity yield criterion : f(X) :

= 1 √ 6 X D + kX m -τ c ≤ 0
The plasticity flow rule :

     ṗD = Tr ṗ k √ 6 X D X D , Tr ṗ ≥ 0 if f(X) = 0, X D = 0 Tr ṗ ≥ k √ 6 ṗD ≥ 0 if f(X) = 0, X D = 0
The damage evolution law:

The irreversibility condition : α ≥ 0

The damage criterion :

Y := -1 2 K 1 (α)(Tr p) 2 -µ 1 (α)p D • p D ≤ D 1
The consistency condition : (Y -D 1 ) α = 0 Remark 3.1. This model can be written in a variational form following the general presentation proposed by [START_REF] Mielke | Evolution of rate-independent systems[END_REF] for rate independent evolution laws. After the choice of the expression of the total energy, sum of the free energy and the dissipated energy, as a function of the state variables, the evolution of the internal variables is deduced from three general physical principles: (i) an irreversibility principle, (ii) a stability criterion, (iii) an energy balance. The interested reader for such an approach could refer to [START_REF] Alessi | Gradient damage models coupled with plasticity: variational formulation and main properties[END_REF][START_REF] Francfort | The taming of plastic slips in von Mises elasto-plasticity[END_REF][START_REF] Francfort | The elasto-plastic exquisite corpse: A Suquet legacy[END_REF] where its application to elasto-plasticity and damage models is developed. That approach can even be used to construct non associative plasticity model, see [START_REF] Babadjian | Quasi-static evolution in nonassociative plasticity: the cap model[END_REF] for an example. Note also that the model has a stress-softening character induced by the damage evolution. That means that the present local model must be regularized (for instance, by introducing gradient of damage terms in the energy) in order that the damage localization in a body be limited and controlled. For the construction of the regularized model, the variational approach should be really useful and even necessary.

Application to the triaxial test with containing pressure

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 remains 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 µ 1 (α) and K 1 (α). However we can obtain some general properties without specifying those functions provided that they satisfy the following conditions:

K 1 (0) = µ 1 (0) = +∞ K 1 (0) = µ 1 (0) = -∞ , K 1 (1) = µ 1 (1) = 0, ( 45
)
∀α ∈ (0, 1) :

     K 1 (α) > 0, µ 1 (α) > 0 K 1 (α) < 0, µ 1 (α) < 0 K 1 (α) > 0, µ 1 (α) > 0 . ( 46
)
Those conditions are inspired by the micro-mechanical model presented in [START_REF] Andrieux | Un modèle de matériau microfissuré pour les bétons et les roches[END_REF]. In particular, the fact that the hardening moduli are infinite when α = 0 ensures that the elastic moduli K 0 and µ 0 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 1 et µ 1 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 1 , i.e. when -1 2 A 1 (α)p • p = D 1 .

Differentiating this relation with respect to t leads to

1 2 A 1 (α)p • p α = -A 1 (α)p • ṗ (47)
from which one deduces that α = 0 if ṗ = 0, which is the desired result. Let us note that this result is essentially due to the fact that the stiffness tensor A 0 is assumed to be damage independent. 2. Damage starts at the same time as plasticity. Indeed, by virtue of the hypothesis K 1 (0) = µ 1 (0) = -∞, the damage criterion at the onset of damage gives

-1 2 K 1 (0)(Tr p) 2 -µ 1 (0)p D • p D ≤ D 1
and 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

3 σ =   -p 0 0 0 0 -p 0 0 0 0 σ z   , ε =   ε T 0 0 0 ε T 0 0 0 ε z   , p =   p T 0 0 0 p T 0 0 0 p z  
where p 0 remains constant and ε z is a given function of time which is decreasing from its critical value given by ( 27). 4. Relationship between Tr p and α. Inserting the form of the different tensors into the plasticity flow rule (which does not involve the damage state) gives

p z -p T = - Tr p 2k . ( 48
)
Reporting this relation into the damage criterion gets Tr p in function of α:

Tr p = 6k 2 D 1 -R (α) ( 49
)
where R(α) is the following combination of the hardening moduli

R(α) := µ 1 (α) + 3k 2 K 1 (α). ( 50
)
Let us note that, since the second derivatives of K 1 and µ 1 are positive, Tr p is a strictly increasing function of α. That means that the required inequality Tr ṗ > 0 is satisfied provided that α > 0. 5. Relationship between σ z and α. The plasticity criterion with the definition of the plasticity thermodynamical forces give the three following equations:

         X T -X z + 3kX m = 3τ c X z -X T = σ z + p 0 + µ 1 (α) Tr p k 3X m = σ z -2p 0 -3K 1 (α) Tr p .
Eliminating X T , X z , X m and using (49) allow us to obtain σ z as the following function of α:

(1 -k)σ z = -3τ c -(1 + 2k)p 0 - 6D 1 S (α) , (51)
where S(α) := 1 R(α) .

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(α) is negative or positive. The former case, when S (α) < 0, corresponds to a stresshardening response, whereas the second one, when S (α) > 0, corresponds to a stresssoftening 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 (49) and ( 51)), one gets:

ε z = σ z E 0 - 1 k -1 Tr p 3 + 2ν 0 p 0 E 0 , (52)
ε v = σ z -2p 0 3K 0 + Tr p . (53)
7. The relation between the axial stress σ z and the axial strain ε z . The relations (51)-( 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 S (α) < 0 (stress-hardening behavior), then εz < 0 when α > 0 and hence damage grows when ε z (which is negative) is decreasing. But, when S (α) > 0, then the sign of εz is not guaranteed when α > 0, because the term in σ z is increasing whereas the term in Tr p 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 Tr p. During the hardening phases, one has σz < 0 and Tr ṗ > 0. 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 σz > 0 and Tr ṗ > 0, one will alway observe a dilatancy.

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

Study of a family of models

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

R(α) = R 1 (1 -α) m α n with R 1 > 0, m > 1, 0 < n < 1 . (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)-( 46). Indeed R(α) is strictly decreasing from +∞ to 0 when α grows from 0 to 1. Its first and second derivatives are given by

       R (α) = -R 1 (m -n)α + n (1 -α) m-1 α n+1 , R (α) = R 1 (m -n)(m -n -1)α 2 + 2n(m -n -1)α + n(n + 1) (1 -α) m-2 α n+2 .
It is also not difficult to verify that R (α) > 0 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 (49), one deduces that Tr p grows from 0 to +∞ with the evolution of damage.

The compliance function S(α) and its first and second derivatives are given by

                 S(α) = 1 R 1 α n (1 -α) m S (α) = 1 R 1 (m -n)α + n α 1-n (1 -α) m+1 S (α) = 1 R 1 (m -n)(m -n + 1)α 2 + 2n(m -n + 1)α -n(1 -n) α 2-n (1 -α) m+2
.

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

α 0 = 1 m -n mn m -n + 1 -n . (55)
Therefore, by virtue of (51), the absolute value of the axial stress σ z increases and one observes a stress-hardening behavior when α grows from 0 to α 0 , then |σ z | decreases and one observes a stress-softening behavior when α grows from α 0 to 1. Moreover, since S (0) = S (1) = +∞, σ z starts from the value given by (26) and corresponding to the end of the elastic stage, then is increasing (in absolute value) until its maximal value corresponding to the time at which the damage reaches the value α 0 , and finally is 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

∆σ z = 6D 1 (1 -k) 2 S (α 0 ) . ( 56
)
Let us first study the variations of the axial strain ε z in function of the damage α and let us analyse under which condition a snap-back occurs. Starting from (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:

S (α) R (α) -R (α) 3/2 S (α) 3/2 ≤ (1 -k) 2 E 0 3 , ∀α ∈ (0, 1). ( 57
)
Since R(α) is of the form R 1 R(α), the left hand side above can read as

S (α) R (α) -R (α) 3/2 S (α) 3/2 =: R 1 ϕ(α) (58)
where the function α → ϕ(α) depends only the exponents m and n of the model. In the neighborhood of α = 1, since the different terms behave as follows

     R (α) ∼ -R 1 m(1 -α) m-1 , R (α) ∼ R 1 m(m -1)(1 -α) m-2 , S (α) ∼ m R 1 (1 -α) -m-1 , S (α) ∼ m(m + 1) R 1 (1 -α) -m-2 one gets ϕ(α) ∼ m + 1 m -1 (1 -α) m .
Therefore, since ϕ is continuous, negative in the interval (0, α 0 ), positive in the interval (α 0 , 1) and since ϕ(α 0 ) = ϕ(1) = 0, the function ϕ reaches its upper bound in the interval (α 0 , 1), the upper bound depending only on m and n. Hence, the condition of non snap-back (57) can read as

R 1 ≤ (1 -k) 2 E 0 3 max α∈[α 0 ,1] ϕ(α) . (59)
In other words, there is no snap-back provided that the hardening modulus R 1 is small enough. Let us note that this condition depends on the parameters E 0 , k, m and n of the model, but is independent of D 1 , τ c , and the confining pressure p 0 .

Let us now study the variations of the volumetric strain ε v as a function of the damage α and let us analysis when one observes a contractance or a dilatancy. We assume that R 1 is small enough and satisfies (59) so that there is no snap-back in the response σ z -ε z and hence εz < 0. By virtue of (53), one has εv < 0 and hence contractance if and only if σz + 3K 0 Tr ṗ > 0.

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

ϕ(α) < -3k(1 -k) K 0 R 1 . ( 60
)
That can happen only during the phase where ϕ is negative, that is only during the stresshardening 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 α < α 0 . Since the behavior of the first and the second derivatives of R and S in the neighborhood of α = 0 is given by

     R (α) ∼ -R 1 nα -n-1 , R (α) ∼ R 1 n(n + 1)α -n-2 , S (α) ∼ n R 1 α n-1 , S (α) ∼ - n(1 -n) R 1 α n-2 , one gets ϕ(α) ∼ - 1 + n 1 -n α -n .
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 α 0 ), one observes a dilatancy until the complete rupture of the volume element. Let us note that since σz tends to 0 when α tends to 1, one gets by virtue of ( 52)-( 53)

εv ≈ Tr ṗ ≈ - 3k 1 -k εz 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 (49) and (51) the following behaviors for Tr p and σ z :

         Tr p = 6k 2 D 1 nR 1 α (1+n)/2 + • • • σ z = - 3τ c 1 -k - (1 + 2k)p 0 1 -k - 6D 1 R 1 (1 -k 2 )n α (1-n)/2 + • • •
Thus the plastic strain is of the order of α (1+n)/2 while the variation of the elastic strain is of 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 (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 0 . After the elastic stage, if the axial strain |ε z | is increased until +∞, then |σ z | begins by increase, 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.

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

In all the pictures presented in the present section, one uses the following values for the Poisson ratio ν 0 and the parameters τ c and k of the Drucker-Prager criterion:

ν 0 = 0.2, τ c = 0, k = 0.2.
One sets

σ c := D 1 E 0 , ε c := D 1 E 0 . (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 1 /E 0 and p 0 /σ c .

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

∆σ z = 1.827 R 1 E 0 σ c . (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 max α∈[α 0 ,1] ϕ(α) = 0.775 and hence the non snap-back condition requires that R 1 ≤ 0.275E 0 .

Are plotted on the figure 1 the graphs of σ z and ε v as functions of ε z when there is no confining pressure (p 0 = 0) and in the case where R 1 = 0.1E 0 (hence there is no snap-back).

In that case the value of the overstress is ∆σ z = 0.578σ c . On the figure 2 same graphs for different values of R 1 /E 0 . One can note that the larger R 1 , the larger the overstress, which is conform to (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 R 1 , the greater the maximal contraction, but also the more rapid the growing of the dilatancy after the peak.

When R 1 is greater than 0.275E 0 , 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 figure 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 Figure 4. This property is due to the linear character of Drucker-Prager criterion. That would not be true any more if one used a more general plasticity criterion like Hoek-Brown criterion [START_REF] Hoek | The Hoek-Brown failure criterion-a 1988 update[END_REF][START_REF] Eberhardt | The Hoek-Brown failure criterion[END_REF]. • Influence of the exponents m and n. The shape of the graphs giving the axial stress in term 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 the 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 harderning state leaving the large |ε z | asymptote intact. Both parameters (m and n) contributes ( 55)-(56) to the overstress ∆σ z estimation, that reads as :

∆σ z = φ(m, n) 1 -k R 1 E 0 σ c
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 Figure 6. • Possible parameter fitting strategy. Even so our model depends just on three "damage" parameters (R 1 , 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 (Tr ε(loading)) and softening curves (σ z (loading)). The softening post-pic behavior would allow us to fit first the m value. Once the m parameter fixed one could proceed to R 1 value choice which allow to adjust the strength of dilatancy curve. Finally the experimental overstress threshold could be fixed by adapting the last n parameter.

Conclusion and Perspectives

Starting from low-scale mechanical properties of micro-cracked elastic materials 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 micro-mechanical considerations is that the free energy contains not only the usual elastic energy but also a stored energy which is due at the micro-level to the friction between the lips of the micro-cracks. 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 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 [START_REF] Comi | A non-local model with tension and compression damage mechanisms[END_REF][START_REF] Pham | Approche variationnelle de l'endommagement: II. Les modèles à gradient[END_REF].

Figure 1 :

 1 Figure 1: Graphs of the axial stress σz (left) and of the volumetric strain εv (right) versus the axial strain εz when R1 = 0.1 E0 and p0 = 0 (uniaxial compression test).

Figure 2 :Figure 3 :Figure 4 :

 234 Figure2: 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.

Figure 5 :

 5 Figure 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 (p0 = 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.

Figure 6 :

 6 Figure 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.

Table 1 :

 1 Table of the notation used throughout the paper for main mechanical quantities

3 2. The elasto-plastic behavior at fixed damage state

Note that, contrarily to a frequent use in Civil Engineering, we keep the convention used in Mechanics of Continuous Media for the sign of the strains and the stresses: the traction are positive and the compression negative.

Throughout the paper, all the tensors are represented by matrices in the basis (ex, ey, ez) where z denotes the axial direction whereas x et y denote the transversal direction.

Let us note that the uniqueness of the response is not guaranteed in presence of damage without making extra assumptions on the functions µ1(α) and K1(α). This study on the uniqueness will not be made here.