A Plastic Damage Model with Stress Triaxiality-Dependent Hardening for Concrete

X.P.Shenand X.C.Wang

1 Introduction

The inelastic failure of concrete-like materials and structures is characterized by the initiation and evolution of cracks and the frictional sliding on the closed crack surfaces.Plastic damage models are the major measures to deal with cracking-related failure analysis,and are widely used by various researchers,see e.g.Lemaitre(1990);Chaboche(1992);Seweryn and Mroz(1998);de Borst,Pamin,Geers(1999).Owing to its simplicity and a reasonable capacity of problem representation,the isotropic damage model is the most popular damage model in the simulation of the failure phenomena of concrete structures and is therefore the choice of this study.Another important reason for the choice of an isotropic damage model in this study is due to that the aim of this investigation is to analyse the cracking process of concrete material of a prestressed concrete structure.Consequently the modelling of pre-peak nonlinearity of stress-strain curve under compression and the stress-triaxiality dependent plastic hardening law and the corresponding damage evolution law are the most important concerns of this study.

One of the plasticity-based isotropic damage models for concrete is the so-called Barcelona model which is reported by Lubliner,Oliver,Oller,Onate(1989),and has recently been adopted by Lee and Fenves(1998)and Nechnech,Meftah and Reynouard(2002).In this model,a holonomic relationship between damage and equivalent plastic strain has been presented,and two damage variables have been designed for tensile damage and compressive damage respectively.It should be noted that most of the existing plastic damage models have not emphasized on the stress-triaxiality dependent plastic hardening law.

The framework of the constitutive model in this study is constructed on the basis of the plasticity-based damage model reported by Saanouni,Forster and Ben Hatira(1994).Although this model was proposed for the simulation of the plastic damage phenomena of metals,it is still an attractive model for the purpose of this study because of its important advantages over the other available damage models:firstly,in this model,the damage evolution is not only closely connected to the increase of plastic strain,but it is also influenced explicitly by the elastic strain;secondly,the damage evolution is coupled with an increase of plastic strain;finally,this model is relatively easy to modify to make it suitable for the simulation of plastic damage phenomena of concrete-like material.The generalized Drucker-Prager criterion introduced in Liebe and Willam(2001)for plastic loading,together with its plastic potential for non-associated plastic flow rule,is referred to here.

The context of the article is organised in the following order:in section 2,the equations of the constitutive model are given for elasto-plasticity coupled with damage in general;a definition of stress triaxiality is proposed and later introduced in the plastic hardening and damage evolution laws.A driver subroutine for validation of the constitutive models is developed with reference to the principle proposed by Hashash,Wotring,Yao,Lee and Fu(2002).In section 3,results of numerical tests of the proposed model are given for some typical loading cases.Some conclusions are given in section 4.

2 Formulation of the proposed model

2.1 Fundamental equations of the plasticity-based damage model

With the ‘Energy Equivalence Principle’,the fundamental relationships of the plasticity-based damage model proposed by Saanouni,Forster and Ben Hatira(1994)are listed in Eqn.(1)as:

where σijand εijare total stress and strain tensors respectively,superscript p stands for plastic and e represents elastic quantities,overhead tilt～represents the quantities for fictitious net materials,is the deviatoric stress tensor,is the inelastic multiplier,D represents the isotropic damage variable,Y is the damage conjugate force,is the elasticity tensor of the intact material,δijis the 2nd order unit tensor,is the sum of the effective principal stresses,is the second invariant of the deviatoric effective stress tensor,F is a plastic damage potential function.

It is seen from the above equations that the damage evolution is designed to be accompanied by plastic strain increase,and its quantity of increment is also dependent on elastic strain tensor via damage conjugate force.

2.2 Specification for Drucker-Prager type plasticity coupled with damage

With reference to the generalized Drucker-Prager criterion introduced in Menétrey and Willam(1995),together with the hardening model introduced in Besson(2001),the plastic damage loading condition is primarily defined in the effective stress space in the following form(stress triaxiality will be introduced later):

where k is initial shear strength constant,k∞is the strain hardening limit of the fictitious net material,which corresponds to infinite equivalent plastic strain,i.e.λ → ∞,and αFis a material constant designed for pressure-sensitivity properties,parameter b is a model constant which can be determined by fitting experimental phenomena.

The plastic part of the potential,i.e.,is given in the effective stress space as:

where αQis the dilatancy constant for non-associated flow rule if αQ/= αF.

The following form of plastic damage potential function F is adopted in order to have non-associate plastic flow in the effective stress space:

wheres,S,φ are material parameters,Dis damage variable,Yis damage conjugate force.

It is observed that the experimental results of stress-strain curves of concrete-like materials are highly dependent on the stress triaxiality.In fact,phenomena of stress triaxiality dependency of the stress-strain curves exist in engineering for a wide range of materials such as geometerials,ceramics,composites and metals.In order to simulate this kind of phenomena,the stress triaxiality was used by several references(see Alves and Jones,1999;Horstemeyer,Lathrop,Gokhale and Dighe,2000;Borvik,Hopperstad and Berstad,2003;Li,Zhang and Ansari,2002).The forms of expressions of stress triaxiality are different from one to another:Alves and Jones(1999),Borvik et al(2003)and Horstemeyer et al(2000)define their stress triaxiality explicitly,while Li,Zhang and Ansari(2002)implicitly account for stress triaxiality influence in its plastic hardening law by using the stress invariants I1 and J2.

Here,for the convenience of model formulation,together with a reference to the conventional expressions adopted in several references(Sfer,Carol,Gettu and Etse 2002;Etse and Willam 1999),the stress triaxiality is defined as:

The stress triaxiality is introduced into the damage plastic loading condition and the damage plastic potential function in the following form:

where

Consequently plastic strain increment is obtained as:

with

It is obtained with Eqn.(5)that:

The damage evolution law can be derived as:

with

The parameters used in this model are:E,ν,k,αF,αQfor plasticity,ands,S,φ for damage.Stress triaxiality γ is a special variable introduced in this model.

2.3 Constitutive behaviour for a finite displacement increment∆εij

With above constitutive model,the constitutive behaviour can be derived for a known initial stress state()and a given strain increment∆εij.The stress increment can be obtained by making total differential operation over total stress tensor in above Eqn.(1)and a subsequent linearization over the time increment∆t,thus

With the above equations,the following equation is obtained:

The plastic damage multipliercan be determined explicitly with the following consistency condition:

Consequently it is obtained the following equation:

Thus,the following expression for plastic multiplier can be derived:

For the sake of brevity,Eqn.(18)can be re-written in another form as:

with

For a given strain increment∆εij,the stress tensor increment can be obtained by substituting Eqn.(20)into Eqn.(14)such that:

Therefore the algorithmic tangential stiffness tensor can be deduced as

The elastoplastic damage loading condition for a given strain increment∆εijcan be expressed conceptually in the effective stress space as

With Eqn.(2),the following relationship is obtained:

With Eqns.(1)and(2),the tensors and vectors on the right hand side of Eqn.(24)are obtained as

The formulation of Newton-Raphson iteration equation between∆λ and˜fis formed as:

3 Numerical validation at local level

Here a driver subroutine is designed for the purpose of validation of 3-dimensional constitutive model at local level,i.e.for a material point only,with reference to the algorithm proposed in Hashash,Wotring,Yao,Lee,and Fu(2002).Its principle can be illustrated with Fig.1 as:a mixed loading condition is applied with ε11= ε11(t),and σ22= σ33=const,which means that a strain loading will be applied incrementally under a constant stress confinement in the other two directions.The strain loading is applied elastically in direction 11,while the self-equilibrium mechanism at this material point will result in variation of lateral strains nonlinearly(for the sake of damage)in order to keep the lateral confinement constant.The details of the numerical calculations will be given in the following context.

Figure 1:Illustration of the mixed loading condition.

Figure 2:Flow chart of global equilibrium iteration.

3.1 Iteration procedure for the constitutive validation:external equilibrium iteration.

The function of the external equilibrium iteration is:for a known stress and strain state and a set of internal variablesapply a load increment(∆ε1,∆σ2,∆σ3),to find out quasi-elastically the response of(∆σ1,∆ε2,∆ε3,D)by an iterative procedure.

The following quasi-elastic equations(i.e.elastic relationship for a finite time increment∆t)are adopted in the calculation:

where Eijklare components of the elasticity tensor of damaged material which is expressed in Eqn.(22).The principle of the external equilibrium iteration is illustrated in Fig.2.

In Fig.2,UMAT is the constitutive module which checks elastoplastic loading state and makes elastic and/or elastoplastic damage calculations.Subroutine CONSTITUERE will be called in the UMAT.The function of subroutine CONSTITUERE is to carry on the constitutive integration and will be introduced in detail in the following sub-section.The elasto-plastic-damage loading stiffnesswhich also known as algorithmic tangential stiffness,will be updated after every iteration,and will be used in the quasi-elastic calculation of∆ε2and ∆ε3at every first iteration step at each of the loading increments.

3.2 Iteration procedure for the constitutive validation:internal elasto-plastic damage iteration

Byusingthe ‘fixed-point’method described in Chaboche and Cailletaud(1996),for a given finite strain ∆εij,the only unknown in the elastoplastic damage calculation at local level is the plastic-damage multiplier∆λ.

The solution steps adopted in the procedure of CONSTITUERE subroutine are:

·Step 1:Initiate the stress state and state of all the internal variables:

·Step 2:Apply strain increment∆εij,with ∆εij=0 fori/=jobtained from the outer global equilibrium iteration;

·Step 3:Calculate∆λ0with given initial stress state,strain increment and linearized Eqn.(18);

·Step4:With this ∆λ0obtained in step3,calculate consequently the following quantities:

·Step 5:With Eqn.(24)and(30),calculate iteratively the plastic-damage multiplier with the following equations:

·Step 6:Check convergence:if∆(∆λ)≤Tolerance,cease the iteration and continue to the next load increment;

·Otherwise,make

·Return to step 2 to carry on the next iterative calculation up to the maximum iteration limit.

3.3 Numerical examples

In this sub-section,numerical validations of the constitutive model at local level are carried out with the driver subroutine developed here for 3 kind of typical loading cases,i.e.,(1)uniaxial tension;(2)uniaxial compression;and(3)uniaxial compressions under various hydrostatic stress confinements.

With reference to the existing literatures(see Lee and Fenves,1998;Ghavamian and Carol,2003;Etse and Willam,1999),the following values of material parameters are adopted in the calculation.They are:

E=31140MPa,ν=0.2,αF=αQ=0.15,k=2.0MPa,s=1,S=10-16 for tension and 4×10-5 for compression,φ=-1.0,b=88 for tension and 500 for compression,k∞=100MPa for tension and 248MPa for compression.Tolerance=10-20 for internal iteration(i.e.for∆(∆λ))and 10-4 for external iteration(i.e.for constant lateral stress confinement).

3.3.1 Uniaxial tension

The stress-stain behaviour under uniaxial tension of the model is shown in Fig.3.No plastic hardening behaviours is observed in this case.Comparison between the experimental data(Gopalaratnam,Shah,1985)and the numerical results of stress-strain response under uniaxial tension indicates that the prepeak stress-strain behaviour can be predicated very well,while the postpeak behaviour can only be predicted with a reasonable accuracy by the proposed model.Fig.4 shows the numerical results of the response of the lateral strain and volumetric strain.The damage response in Fig.5 shows that damage value asymptotically tends to maximum 1.0 with the increase of strain loading.

3.3.2 Uniaxial compression

The stress-stain behaviour and damage evolution response of the proposed model under uniaxial compression are shown in the following Fig.s.Comparison between the experimental data given by Karsan,Jirsa(1969)and the numerical results of stress-strain response under uniaxial compression in Fig.6 indicates that the prepeak behaviour can be predicated very well,and the postpeak behaviour can be predicted with a reasonable accuracy.The numerical results of the response of the lateral strain and volumetric strain shown in the same Fig.6 indicate the dilatancy property of the model:there is a saturation of dilatancy at which the volumetric strain ceases to increase.The corresponding damage evolution behaviour is presented in Fig.10.

Figure 4:Stress-strain behaviours under uniaxial tension:lateral and volumetric responses.

Figure 5:Damage evolution under uniaxial tension.

Figure 6:Stress-strain behaviours under uniaxial compression.

3.3.3 Compressions with confinement

The stress-strain behaviour of a model for concrete under hydro-static stress confinement is an important aspect:it indicates the pressure-sensitivity behaviour of the model.In the numerical tests performed here,the procedure of loading is given as:the hydro-static confinement,i.e.σmI,is applied before strain loading in direction 11 being carried on.Fig.7 shows the variation of the stress-strain response caused by the confinement of the stress-stain behaviour and damage evolution response,with the other parameters were kept unchanged:with the increment of the stress confinement,the softening phenomena become weaker and weaker.It seems that the stress-triaxiality dependent plastic hardening phenomena are properly simulated by the model proposed here.

In Fig.8 and 9,the stress-strain curves under confinement of-10MPa and-30MPa are given respectively,together with the responses of ε22- σ11and εv-σ11.The dilatancy phenomena become weaker with the increase of confinement pressure.The prepeak nonlinearity of the stress-strain curve seems reasonably simulated.

Figure 7:influence of hydrostatic stress confinement:stress-strain behaviour(σm=0,-10,-30 MPa).

In order to illustrate the proposed model further in-depth,in Fig.11 and 12,the peak strength envelopes obtained numerically with the proposed model are given out.Because of the different parameter values adopted in the calculation for tension and compression,it is seen in Fig.11 that the shape of tensile strength envelope is quite different from the shape of the compressive strength envelope.The axial points in Fig.11 are obtained analytically by using the Dracker-Prager condition directly,because there is no prepeak nonlinearity for tensile case,which is judged by criterionI1≥0.In Fig.12,the envelope of peak-strengths of compression under very high confinement up to-200MPa is shown in order to show the validity of the model for a wide range of hydrostatic stress confinement.

Figure 8:Stress-strain behaviours under compression with-10MPa confinement.

Figure 9:Stress-strain behaviours under compression with-30MPa confinement.

Figure 10::Damage-strain behaviours under compression with various stress confinements.

Figure 11::Peak strength envelope in thespace(with confinements up to-30MPa).

Figure 12::Peak strength envelope in thespace(with confinements up to-200MPa).

4 Conclusions

In this article,a plasticity-based isotropic damage model has been proposed with special concerns being given to the modelling of the stress-triaxiality dependent plastic hardening law and corresponding damage evolution.Numerical validations have been made for the proposed model with a driver subroutine developed in this study.The numerical results of stress-strain behaviour seem reasonably accurate for the uniaxial tension and uniaxial compression cases compared with the experimental data of the references.The phenomena of stress-triaxiality-dependent plastic hardening have been simulated.Owing to the highly nonlinearity existing among the model parameters and the experimental phenomena,it seems necessary to choose values of constitutive parameters by some kind of technique of inverse analysis.

Acknowledgement:Financial support from National Natural Science Foundation of China(contract No.11272116),and support from The Importation and Development of High-Caliber Talent Project of Beijing Municipal Institutions,are gratefully acknowledged.

Alves,M.;Jones,N.(1999):influence of hydrostatic stress on failure of axisymmetric notched specimens.J.Mech.Phys.Solids,vol.47,pp.643-667.

Besson,J.(2001):Mecanique et Ingenierie des Materiaux-Essais mecaniques,Eprouvettes axisymetriques entaillees,Hermes(in French).

Borvik,T.;Hopperstad,O.S.;Berstad,T.(2003):On the influence of stress triaxiality and strain rate on the behaviour of a structural stell.Part II.Numerical study.Eur.J.Mech.A/Solids,vol.22,pp.15-32.

Chaboche,J.L.(1992):Damage induced anisotropy:on the difficulties associated with the active/passive unilaterial condition.Int.J.Dama.Mech.,vol.1,pp.148-171.

Chaboche,J.L.;Cailletaud,G.(1996):Integration methods for complex plastic constitutive equations.Comp.Methods.Appl.Mech.Engng.,vol.133,pp.125-155.

de Borst,J.;Pamin,M.G.;Geers,D.(1979):On coupled gradient-dependent plasticity and damage theories with a view to localization analysis.Eur.J.Mech.A/Solids,vol.18,pp.939-962.

Dragon,A.;Mroz,Z.(1995):A continuum model for plastic-brittle behaviour of rock and concrete.Int.J.Mech.Sci.,vol.17,pp.121-137.

Etse,G.;Willam,K.(1999):Failure analysis of elastoviscoplastic material models.ASCEJ.Engng.Mech.,vol.125,no.1,pp.60-69.

Ghavamian,S.;Carol,I.(2003):Benchmarking of concrete cracking constitutive laws:MECA project.In:Computational Modelling of Concrete Structures.N.Bicanic,R.de Borst,H.Mang,G.Meschke(Eds.),Swets&Zeitinger,Lisse,pp.179-187.

Gopalaratnam,V.S.;Shah,S.P.(1985):Softening response of plain concrete in direct tension.ACI J.,vol.82,no.3,pp.310-323.

Hashash,Y.M.A.;Wotring,D.C.;Yao,J.I.C.;Lee,J.S.;Fu,Q.(2002):Visual framework for development and use of constitutive models.Int.J.Numer.Anal.Meth.Geomech.,vol.26,pp.1493-1513.

Horstemeyer,M.F.;Lathrop,J.;Gokhale,A.M.;Dighe,M.(2000):Modeling stress state dependent damage evolution in a cast Al-Si-Mg aluminium alloy.Theo.Appl.Frac.Mech.,vol.33,pp.31-47.

Karsan,I.D.;Jirsa,J.O.(1969):Behaviour of concrete under compressive loading.ASCEJ.Engng Mech,vol.95,no.12,pp.2535-2563.

Lee,J.;Fenves,G.L.(1998):Plastic-damage model for cyclic loading of concrete structures.ASCEJ.Engng Mech.,vol.124,no.8,pp.892-900.

Lemaitre,J.(1990):A Course on Damage Mechanics,2nd ed.,Berlin:Springer.

Li,Q.;Zhang,L.;Ansari,F.(2002):Damage constitutive for high strength con-crete in triaxial cyclic compression.Int.J.Solids Struct.,vol.39,no.15,pp.4013-4025.

Liebe,T.;Willam,K.(2001):Localization properities of generalized Drucker-Prager elastoplasticity.ASCEJ.Engng Mech.,vol.127,no.6,pp.616-619.

Lubliner,J.;Oliver,J.;Oller,S.;Onate,E.(1989):A plastic damage model for concrete.Int.J.Solids&Struct.,vol.25,no.3,pp.299-326.

Nechnech,W.;Meftah,F.;Reynouard,J.M.(2002):An elasto-plastic damage model for plain concrete subjected to high temperatures.Engng.Strut.,vol.24,no.5,pp.597-611.

Saanouni,K.;Forster,C.;Hatira,F.B.(1994):On the anelastic flow with damage.Int.J.Dama.Mech.,vol.3,pp.140-169.

Seweryn,A.;Mroz,Z.(1998):On the criterion of damage evolution for variable multiaxial stress state.Int.J.Solids&Struct.,vol.35,no.14,pp.1589-1616.

Sfer,D.;Carol,I.;Gettu,M.R.;Etse,G.(2002):Study of the behavior of concrete under triaxial compression.ASCEJ.Engng Mech,vol.128,no.2,pp.156-163.