下载文档原格式
Unit 2 Classification of MaterialsSolid materials have been conveniently grouped into three basic classifications: metals, ceramics, and polymers. This scheme is based primarily on chemical makeup and atomic structure, and most materials fall into one distinct grouping or another, although there are some intermediates. In addition, there are three other groups of important engineering materials —composites, semiconductors, and biomaterials.译文：译文：固体材料被便利的分为三个基本的类型：金属，陶瓷和聚合物。
固体材料被便利的分为三个基本的类型：金属，陶瓷和聚合物。
固体材料被便利的分为三个基本的类型：金属，陶瓷和聚合物。
这个分类是首先基于这个分类是首先基于化学组成和原子结构来分的，化学组成和原子结构来分的，大多数材料落在明显的一个类别里面，大多数材料落在明显的一个类别里面，大多数材料落在明显的一个类别里面，尽管有许多中间品。
尽管有许多中间品。
除此之外，此之外， 有三类其他重要的工程材料－复合材料，半导体材料和生物材料。
有三类其他重要的工程材料－复合材料，半导体材料和生物材料。
Composites consist of combinations of two or more different materials, whereas semiconductors are utilized because of their unusual electrical characteristics; biomaterials are implanted into the human body. A brief explanation of the material types and representative characteristics is offered next.译文：复合材料由两种或者两种以上不同的材料组成，然而半导体由于它们非同寻常的电学性质而得到使用；生物材料被移植进入人类的身体中。
A constitutive model for the anisotropicelastic–plastic deformation of paper and paperboardQingxi S.Xia,Mary C.Boyce *,David M.ParksDepartment of Mechanical Engineering,Centre of Materials Science and Technology,Massachusetts Institute of Technology,77Massachusetts Avenue,Cambridge,MA 02139-4307,USAReceived 20June 2001;received in revised form 19February 2002AbstractA three-dimensional,anisotropic constitutive model is presented to model the in-plane elastic–plastic deformation of paper and paperboard.The proposed initial yield surface is directly constructed from internal state variables comprising the yield strengths measured in various loading directions and the corresponding ratios of plastic strain components.An associated flow rule is used to model the plastic flow of the material.Anisotropic strain hardening of yield strengths is introduced to model the evolution in the yield surface with strain.A procedure for identifying the needed material properties is provided.The constitutive model is found to capture major features of the highly anisotropic elastic–plastic behavior of paper and paperboard.Furthermore,with material properties fitted to experimental data in one set of loading directions,the model predicts the behavior of other loading states well.Ó2002Published by Elsevier Science Ltd.Keywords:Paper;Paperboard;Plasticity;Yield;Yield surface;Constitutive behavior;Anisotropic1.IntroductionPaper and paperboard are two of the most commonly utilized materials in nearly every industry.Paper is formed by draining a suspension of fibers in a fluid through a filter screen to form a sheet of pulp fibers.Paperboard is in general composed of several pulp fiber sheets bonded by starch or adhesive material,and is usually a multi-layered structure.Schematics of typical paper and paperboard macrostructure and micro-structure are shown in Fig.1,which also depicts the common nomenclature for the three orthogonal di-rections of paper and paperboard.‘‘MD’’refers to the machine (rolling)direction,and ‘‘CD’’refers to the cross or transverse direction.The machine and cross directions form the plane of the structure,and ZD refers to the out-of-plane (or through-thickness)direction.Due to the continuous nature of the paper-making process,fibers are primarily oriented in the plane;furthermore,within the plane,fibers are more highly oriented in the MD than the CD.In this paper,to simplify notation,the 1-direction is used to represent the MD,the 2-direction for ZD and the 3-direction forCD.International Journal of Solids and Structures 39(2002)4053–/locate/ijsolstr*Corresponding author.Tel.:+1-617-253-2342;fax:+1-617-258-8742.E-mail address:mcboyce@ (M.C.Boyce).0020-7683/02/$-see front matter Ó2002Published by Elsevier Science Ltd.PII:S 0020-7683(02)00238-XThe preferential fiber orientation results in highly anisotropic mechanical behavior,including aniso-tropic elasticity,initial yielding,strain hardening and tensile failure strength:•Anisotropic elastic constants for paper and paperboard have been measured by several investigators (e.g.,Mann et al.,1981;Castegnade et al.,1989;Persson,1991;Koubaa and Koran,1995).Their data show that the through-thickness moduli are at least two orders of magnitude less than the in-plane mo-duli.In-plane data of Persson (1991)and Stenberg et al.(2001a)show that moduli in the MD are 2–4times greater than those of the CD.•Persson’s (1991)data on paperboard also shows that the initial yield strength ($proportional limit)in the through-thickness direction is two orders of magnitude lower than the in-plane initial yield strength values.Stenberg’s data (Stenberg et al.(2001a,b);Stenberg,2001)on multi-layer paperboard and single-layer pulp shows similar results.Within the plane,these data show that the initial yield strength of paper and paperboard in the MD is typically greater than that in the CD by a factor of 2–4.Stenberg’s data (Stenberg,2001)also show an asymmetry in the initial yield strength for in-plane tension and compres-sion in both the machine and cross directions.•The in-plane tensile stress–strain curves of Persson (1991)and Stenberg (2001)show substantial strain hardening in which the yield strength increased by more than a factor of two after a strain of less than 5%.The in-plane strain hardening is also highly anisotropic.The percentage strain hardening achieved in MD tension is greater than that obtained in CD tension.•The Persson (1991)and Stenberg (2001)data also show that the in-plane tensile failure strength in the MD is greater than that of the CD by a factor of 2–4.•Biaxial failure stress loci (Gunderson,1983;deRuvo et al.,1980;Fellers et al.,1981)show substantially different failure strengths in machine and cross directions.These data also show that failure tends to be dominated by one or the other of these two directions when subjected to combined loading in both di-rections.Some material models have been proposed to describe the mechanical behavior of paperboard.These models fall into roughly three categories:network models,laminate models,and anisotropic models of the yield surface and/or the failure surface.Perkins and Sinha (1992)and Sinha and Perkins (1995)described a micromechanically based network model for the in-plane constitutive behavior of paper.A meso-element was constructed to represent theFig.1.Schematics of paper and paperboard macrostructure and microstructure.4054Q.S.Xia et al./International Journal of Solids and Structures 39(2002)4053–4071Q.S.Xia et al./International Journal of Solids and Structures39(2002)4053–40714055 microstructure of thefibrous paper network.The mechanical response of the meso-element depends on the fiber properties and properties of the inter-fiber bonds.They found the inelastic behavior of the inter-fiber bonds to play a crucial role in the overall in-plane inelastic behavior of paper.Stahl and Cramer(1998)also developed a network model for low densityfibrous work models can incorporate microlevel mechanisms,such as inter-fiber interaction and bonding.While these models begin to elucidate the un-derlying mechanisms of deformation,they do not provide a continuum-level description of paper or paperboard.Page and Schulgasser(1989)and Schulgasser and Page(1988)developed models of paperboard based on classical laminate theory.While this type of model can predict the elastic response well,it was not extended to capture the anisotropic yielding and subsequent strain hardening response.Gunderson(1983),Gunderson et al.(1986),deRuvo et al.(1980)and Fellers et al.(1981)each used the Tsai–Wu quadratic yield condition to model the failure loci they obtained experimentally.The quadratic nature of this type of model has many shortcomings when applied to paper and paperboard.Experimental data showed the biaxial failure locus to be distinctly non-quadratic.Arramon et al.(2000)developed a multidimensional anisotropic strength criterion based on Kelvin modes that captures the non-quadratic failure envelope.They applied the model to form a strength envelope for paperboard by constructing tensile and compressive modal bounds.However,these efforts only acted to studyfinal failure and did not attempt to study initiation of yield or subsequent strain hardening.In this research,a general three-dimensional constitutive model of the anisotropic elastic–plastic be-havior of paper and paperboard is proposed.The initial elastic behavior is modelled to be linear and or-thotropic.The onset of plasticflow is captured by a non-quadratic yield surface.The yield surface is taken to evolve anisotropically with a scalar measure of plastic strain,with plasticflow modelled using an as-sociatedflow rule.The model is detailed in the following sections,and numerical results are compared to experimental data.2.Experimentally observed behavior2.1.Elastic–plastic behavior of TRIPLEX paperboardDepending onfiber type,fiber density and the chemical/mechanical treatment,the elastic–plastic be-havior of different types of paper and paperboard differs in detail.However,general characteristics of the response remain similar.In this contribution,the anisotropic elastic–plastic behavior of paper and paper-board is illustrated using TRIPLEX e1paperboard as an exemplar material.TRIPLEX e is afive-layer paperboard:the three-layer core is constructed from mechanically processed softwood pulp(commonly referred to as‘‘mechanical’’pulp),and two outer layers(sandwiching the core)are constructed from bleached kraft pulp(commonly referred to as‘‘chemical’’pulp).Stenberg et al.(2001a,b)and Stenberg (2001)conducted an extensive experimental investigation documenting the stress–strain behavior of TRIPLEX e.Note that the outer chemical pulp layers are typically stiffer and stronger than the inner mechanical layers;however,these layers cannot be separately produced for individual evaluation.There-fore,the behavior of TRIPLEX e material will be presented in terms of its effective composite behavior.In principle,testing and modeling can be applied separately for each distinct paper lamina.The experimental results of Stenberg et al.(2001a,b)and Stenberg(2001)are reviewed below.1TRIPLEX e is a trademark of STORA-ENSO,Finland and Sweden.2.1.1.In-plane behaviorThe in-plane uniaxial tensile stress–strain curves for the MD,the CD and an orientation 45°from the MD are plotted together in Fig.2.These stress–strain curves clearly depict the anisotropic in-plane elastic,initial yield and strain hardening behavior.There is a factor of 2–3difference in the modulus and initial yield strength between MD and CD.Hardening achieved in MD (flow strength increases by 300%over a strain of 2%)is higher than that in CD (flow strength increases by 200%over a strain of 5%).MD–CD shear properties are deduced from the 45°test result.In-plane lateral strain (CD)vs.axial strain data for MD-tension and similar data for CD-tension were also recorded by Stenberg (2001),corresponding to their in-plane stress strain curves.Upon subtracting the respective elastic strain components,the lateral plastic strains for both the MD and CD tension cases are computed and shown vs.the respective axial plastic strains in Fig.3.These two curves indicate that for both test orientations,the ratio between lateral plastic strain and axial plastic strain is nearly constant until final fracture.This data provides information for later construction of the plastic flowrule.4056Q.S.Xia et al./International Journal of Solids and Structures 39(2002)4053–4071Tensile loading/unloading/reloading data (Persson,1991;Stenberg,2001)show that after various amounts of plastic strain,upon unloading,the elastic tensile modulus is nearly unaffected by plastic strain,consistent with traditional models of elasto-plasticity.Fig.4shows the in-plane compression stress–strain curves for the MD and CD directions.Note that global specimen buckling was constrained in these tests.These data show that compressive yield is an-isotropic.Furthermore,a comparison of Figs.2and 4shows a yield strength difference between tension and compression,with the compressive yield strengths being smaller than those in tension by 65%and 25%,for MD and CD,respectively.The anisotropic in-plane elastic–plastic properties obtained from these tests are summarized in Table 1.2.1.2.Out-of-plane behaviorStenberg et al.(2001a,b)experimentally obtained the out-of-plane stress–strain behavior of paperboard using a modified Arcan design (Arcan et al.,1978).Nominal stress–strain curves were obtained for TRI-PLEX e under various through-thickness loading conditions.Representative ZD tensile and through-thickness shear (ZD–MD)stress–strain curves obtained by Stenberg et al.(2001a,b)are shown in Fig.5.The stress measure is force per unit initial cross-sectional area;the x -ordinate is the nominal strain,defined as the relative normal/shear separation of the top and bottom surfaces of the laminate,divided by the initial laminate thickness.In the tensile curve,at the earliest stage ofdeformation,the stress increases linearly with strain,exhibiting a composite modulus of E 0ZD ¼18:0MPa.The stress–strain relation shows a small amount of pre-peak non-linearity before reaching a peak stress of 0.4MPa.After the peak,the stress–strain curve exhibits pronounced softening.Features similar to those of the through-thickness tensile curve are observed for the through-thickness shear curve.Thecomposite Table 1Experimental results of uniaxial tensile tests (Stenberg,2001)Elastic modulus(GPa)Poisson’s ratio Tensile yield strength (MPa)Plastic strain ratio,d p ?=d p k Compressive yield strength (MPa)MD5.60.3712.0)0.57.3CD2.00.12 6.5)0.133 5.045° 4.18.0Q.S.Xia et al./International Journal of Solids and Structures 39(2002)4053–40714057transverse shear modulus is observed to be G 12¼34MPa,and the peak shear stress is 1.1MPa.The through-thickness ZD–CD shear stress–strain curve shows similar features.Through-thickness tensile/shear testing conducted within a scanning electron microscope (Dunn,2000)on the same material reveals the nucleation of multiple inter-laminar microcracks near the peak stress,followed by their growth and coalescence,resulting in the observed softening.Similar results have been obtained for tests involving combined through-thickness tension and shear (Dunn,2000;Stenberg et al.(2001a,b)).Therefore,the observed peak stress and subsequent softening for through-thickness loading are due to delamination of the paperboard and will not be further considered here.This inelastic through-thickness behavior is modelled in a separate related study (Xia et al.,in preparation).It is also observed that the amount of lateral (in-plane)strain generated during the through-thickness tensile loading is negligible,indicating that Poisson’s ratios m 21and m 23are close to zero.Through-thickness compression loading/unloading stress–strain curves were also obtained by Stenberg (2001).Up to a nominal compressive strain of 3%,the compressive stress increases linearly with strain.With larger strains,the stress starts to increase exponentially with strain.The data also show that only a small amount of permanent deformation is remained after unloading from a peak strain level of more than 20%,indicating non-linear elastic ZD compressive behavior up to moderately large strains.These obser-vations of through-thickness compressive behavior will be incorporated into the modeling work.The anisotropic linear elastic out-of-plane properties are summarized in Table 2.3.Constitutive modelA three-dimensional,finite deformation constitutive model for paper and paperboard pulp layers is proposed.The model will take the in-plane behavior to be elastic–plastic and the out-of-plane behavior to be elastic.Due to the assumption of elastic out-of-plane behavior,the application of the model alone will be limited to predominant in-plane loading.However,when the model is combined with interlaminardeco-Table 2Elastic out-of-plane properties (Stenberg et al.(2001a,b))Elastic modulusE 0ZD (MPa)Poisson’s ratio m 21Poisson’s ratio m 23Shear modulus G 12(MPa)Shear modulus G 23(MPa)Stiffening constant a 18.0)0.0055)0.003534.026.0 5.44058Q.S.Xia et al./International Journal of Solids and Structures 39(2002)4053–4071hesion models(Xia et al.,in preparation),a general-purpose tool is achieved for modeling behavior under significant out-of-plane loading,such as occurs during converting processes and the in-service behavior ofa broad class of paper and paperboard products.3.1.Stress–strain relationshipFirst,the total deformation gradient F at a material point within a lamina is multiplicatively decom-posed into an elastic part and a plastic part:F¼F e F p;ð1Þwhere F p represents the accumulation of inelastic deformation.Although the maximum in-plane strain level in traditional applications of paper sheets is small,we adopt the presentfinite deformation formulation so that the model can be easily applied to applications such as paperboard converting processes,which generally involvefinite rotations of paperboard layers and may exhibit moderately large,but highly lo-calized in-plane strains.The evolution of F p is given by_F p¼L p F p;ð2Þwhere L p is the plastic velocity gradient,which will be defined by theflow rule.The elastic strain is obtained by using the elastic Green strain measure:E e¼12F eT F eÂÀI Ã;ð3Þwhere I is the second-order identity tensor.The second Piola–Kirchoffstress measure,T,relative to the plastically deformed configuration F p,is then calculated using the linear relationship: T¼C E e½ ;ð4Þwhere C is the fourth-order stiffness tensor,which is taken to be orthotropic.Values of the components of C are defined by the orthotropic elastic moduli.To model the through-thickness non-linear elastic com-pressive stress–strain relationship,the through-thickness engineering elastic constant,E ZD,is taken to be an exponential function of the ZD strain under compression as follows:E ZD¼E0ZD expðÀaE e22Þ;ð5Þwhere E0ZD is the initial elastic modulus in ZD,E e22is the ZD elastic Green strain component,and a is aconstant determined byfitting the compressive through-thickness stress–strain curve;its value is listed in Table2.The stiffness tensor C under ZD compression is in turn determined assuming constant Poisson’s ratios.The Cauchy stress,T,is calculated from its relation to the second Piola–Kirchoffstress by T¼ðdet F eÞF eÀ1TF eÀT:ð6Þ3.2.Yield conditionThe through-thickness strengths(tensile,compressive,shear)of paper and paperboard materials are typically two orders of magnitude lower than those observed in the plane.Therefore,the through-thickness stress components play little role in the inelastic deformation and failure of paperboard under in-plane loading.2Furthermore,from investigation of the mechanisms of deformation and failure of paperboard2Under significant in-plane compression or very large through-thickness compressive strain,this may not be strictly true.However, for the present,these scenarios will not be considered further.Q.S.Xia et al./International Journal of Solids and Structures39(2002)4053–40714059under through-thickness loading (when through-thickness compression is not dominant),it is clear that the majority of through-thickness inelastic deformation occurs in the form of interlaminar microcracking and delamination of the discrete pulp layers,as opposed to inelastic through-thickness deformation distributed quasi-homogeneously within laminae.Thus we can assume that only the in-plane stress components will drive the in-plane inelastic deformation of the pulp layers.Additionally,in classical metal plasticity and in plasticity-based models of yielding in polymers,the deviatoric stress is taken to drive yield because the underlying deformation mechanisms are governed by shear (for example,dislocation glide in crystalline metals).However,in the case of paper,there is no evidence that yield and subsequent plastic flow are driven by deviatoric stress.Micromechanically,yielding is governed by various forms of inter-fiber interactions.Based on these considerations,the total stress will be taken to drive the yield condition described in this research.In order to experimentally define an in-plane yield surface for paper,multi-axial data is required.Al-though the anisotropy of the yield surface is well-recognized,due to numerous studies of the uniaxial behavior in different directions,such as that reported in Fig.2,a literature search reveals virtually no data on the initial and evolving multi-axial yield surfaces of paper and paperboard.However,several researchers have obtained biaxial failure surfaces of paperboard under combinations of MD and CD normal stress.Fig.6shows a representative biaxial failure locus (deRuvo et al.,1980).Similar data have been reported in work by Fellers et al.(1981)and Gunderson (1983).One common feature observed in these failure loci is that the data points tend to form prominent sub-groups,with each sub-group lying on a nearly straight line.For example,for low,but non-zero values of MD stress,failure occurs at nearly the CD uniaxial tensile failure strength;similarly,for low,but non-zero values of CD stress,failure occurs at roughly the MD uniaxial tensile failure strength.This observation suggests that the experimental biaxial failure locus can be well captured by a set of straight lines in two-dimensional biaxial stress space and can be generalized to planes in three-dimensional space.Karafillis and Boyce (1993)and Arramon et al.(2000)developed yield surface and failure surface models,respectively,which capture this non-quadratic feature.In stress space,these lines or planes can be defined by their minimum distance to the origin,together with their corre-sponding normal directions.Given that a comprehensive set of experimental data is generally unavailable (and indeed,is a challenging task to obtain)to determine the full surface,we hereby assume that the yield surface exhibits the same characteristic features observed in the failure surface.Therefore,the yield surface is taken to be constructed of N sub-surfaces,where N K is the normal to the K th such surface,defined in the material coordinates formed by MD,CD and ZD.S K is the equivalent strength of the K th sub-surface,defined by the distance from the origin to each sub-surface,following its normal direction.Thus,the following form of yield criterion isproposed:Fig.6.Biaxial failure surfaces (deRuvo et al.,1980).4060Q.S.Xia et al./International Journal of Solids and Structures 39(2002)4053–4071fðT;S KÞ¼X NK¼1vKTÁN KS K!2kÀ1;ð7Þwhere v K is the switching controller withv K ¼1;if TÁNi>0;0;otherwise;&ð8ÞT is the second Piola–Kirchoffstress measure relative to the F p configuration,and2k is an even integer.N K is the normal of the K th sub-surface,defined relative to the material coordinates.A schematic of a four sub-surface system for biaxial loading,with zero in-plane shear stress,is shown in Fig.7.The normals and corresponding sub-surface strengths are illustrated.The parameter2k is taken to be equal to or larger than4,indicating a non-quadratic yield surface.Fig.8shows the effect of different2k values in controlling the shape of the yield surface in the biaxial stress plane for this simplified four sub-surface system.Higher2k-values give rise to sharper corners between adjacent sub-surfaces and reduce the curvature over increasing central portions of each sub-surface.A schematic of a six sub-surface yield surface,with non-zero in-plane shear stress T13,is shown in Fig.9.Thisfigure graphically illustratestheFig.8.Effect of constant2k on the shape of yield surface.Q.S.Xia et al./International Journal of Solids and Structures39(2002)4053–40714061normals and corresponding equivalent strengths of the sub-surfaces.For this yield surface,the six normals are taken to be of the following form:N K¼X3i;j¼1N Kije i e j;K¼I;...;VIð9ÞandN K ji ¼N Kij;K¼I;...;VI;ð10Þwhere e i are the basis vectors for the material coordinates formed by the MD,CD and ZD.Here,the sub-surface normal index K ranges over the six Roman numerals I–VI.Because out-of-plane stress components are assumed to have no effect on the plastic deformation within a lamina,components of N K involving the2-direction(i.e.,the ZD direction)are set to be zero.This results in elastic through-thickness behavior,as proposed.Determination of each non-zero entry of these matrices will be discussed in Section3.3.3.3.Flow ruleThe plasticflow rule is defined asL p¼D p¼_ c K;ð11Þwhere L p is the plasticflow rate and D p is the symmetric part of L p.For paper and paperboard,the in-plane plastic strains(even at failure)are small;therefore we take the skew part of L p to vanish,or W p¼0,as a simplification.K is theflow direction,and_ c is the magnitude of the plastic stretching rate.K is a second-order tensor with unit magnitude:K¼b K=k b K k;ð12Þwherek b K k¼ffiffiffiffiffiffiffiffiffiffiffiffiffib KÁb Kp¼ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiX3i;j¼1b K ij b K ijv uu t:ð13ÞFig.9.Schematic of yield surface for general in-plane loading,showing six sub-surfaces.4062Q.S.Xia et al./International Journal of Solids and Structures39(2002)4053–4071In Fig.3,the in-plane lateral plastic strain vs.axial plastic strain data showed that the ratio between these two plastic strain components is nearly constant for both the MD and CD simple tension cases.These ratios are taken to define the normal directions of the two respective sub-surfaces of the tensile quadrant of the biaxial yield surface,in the absence of shear stress.Thus,the plasticflow of the material is taken to follow an associatedflow rule:b K¼ð14ÞWith the yield condition defined in Eq.(7),b K can be further calculated asb K¼o T ¼X NK¼12k K2kÀ1KvKN KS K;ð15ÞwhereK K¼TÁN KS K:ð16ÞFor the six sub-surface yield surface shown in Fig.9,assuming an associatedflow rule,the sub-surface normals N K,K¼I;...;VI are defined using the corresponding plastic strain ratios.For example,for sub-surface I(one)in Fig.9,the plastic strain ratio from the MD simple tension test is nearly constant at)0.5. The two non-zero components of N I are then determined by solving the following two equations: N I3311¼À0:5ð17Þand(to make a unit normal)ðN I11Þ2þðN I33Þ2¼1;ð18Þwhich gives N I11¼2=ffiffiffi5pand N I33¼À1=ffiffiffi5p.Similarly,the plastic strain ratio from the CD simple tension testis nearly constant atÀ2=15,giving N II11¼À2=ffiffiffiffiffiffiffiffi229pand N II33¼15=ffiffiffiffiffiffiffiffi229p.With appropriate experimentaldata,similar calculations can determine the normals of each of the sub-surfaces.However,currently,there is no experimental data for the plastic strain ratios for compression in either the machine or cross directions. For the four sub-surface biaxial yield surface shown in Fig.7,the normals for the two sub-surfaces in the compressive quadrants,IV and V,are assumed to have normal directions parallel to those of corresponding tensile sub-surfaces,I and II,respectively,as seen in thefigure,but with generally differing strength levels. For the normal of the sub-surface representing yielding under positive pure shear stress(T13¼T31¼0;T ij¼0otherwise),the two non-zero components are N III13and N III31.Due to the symmetry of shear stress,N III 13¼N III31ð19ÞandðN III13Þ2þðN III31Þ2¼1:ð20ÞThus we have N III13¼N III31¼ffiffiffi2p=2,and the normal of the other sub-surface representing yielding undernegative pure shear stresses is taken to be N VI¼ÀN III.In summary,the components of the normals for the six sub-surfaces determined in the material coordinates are given in Table3.3.4.Strain hardening functionsTo capture the anisotropic strain hardening observed for the in-plane behavior,the equivalent strengths, S K,of each sub-surface are taken to evolve with the accumulated equivalent plastic strain,i.e.。
A plasticity and anisotropic damage model forplain concreteUmit Cicekli,George Z.Voyiadjis *,Rashid K.Abu Al-RubDepartment of Civil and Environmental Engineering,Louisiana State University,CEBA 3508-B,Baton Rouge,LA 70803,USAReceived 23April 2006;received in final revised form 29October 2006Available online 15March 2007AbstractA plastic-damage constitutive model for plain concrete is developed in this work.Anisotropic damage with a plasticity yield criterion and a damage criterion are introduced to be able to ade-quately describe the plastic and damage behavior of concrete.Moreover,in order to account for dif-ferent effects under tensile and compressive loadings,two damage criteria are used:one for compression and a second for tension such that the total stress is decomposed into tensile and com-pressive components.Stiffness recovery caused by crack opening/closing is also incorporated.The strain equivalence hypothesis is used in deriving the constitutive equations such that the strains in the effective (undamaged)and damaged configurations are set equal.This leads to a decoupled algo-rithm for the effective stress computation and the damage evolution.It is also shown that the pro-posed constitutive relations comply with the laws of thermodynamics.A detailed numerical algorithm is coded using the user subroutine UMAT and then implemented in the advanced finite element program ABAQUS.The numerical simulations are shown for uniaxial and biaxial tension and compression.The results show very good correlation with the experimental data.Ó2007Elsevier Ltd.All rights reserved.Keywords:Damage mechanics;Isotropic hardening;Anisotropic damage0749-6419/$-see front matter Ó2007Elsevier Ltd.All rights reserved.doi:10.1016/j.ijplas.2007.03.006*Corresponding author.Tel.:+12255788668;fax:+12255789176.E-mail addresses:voyiadjis@ (G.Z.Voyiadjis),rabual1@ (R.K.AbuAl-Rub).International Journal of Plasticity 23(2007)1874–1900U.Cicekli et al./International Journal of Plasticity23(2007)1874–19001875 1.IntroductionConcrete is a widely used material in numerous civil engineering structures.Due to its ability to be cast on site it allows to be used in different shapes in structures:arc,ellipsoid, etc.This increases the demand for use of concrete in structures.Therefore,it is crucial to understand the mechanical behavior of concrete under different loadings such as compres-sion and tension,for uniaxial,biaxial,and triaxial loadings.Moreover,challenges in designing complex concrete structures have prompted the structural engineer to acquire a sound understanding of the mechanical behavior of concrete.One of the most important characteristics of concrete is its low tensile strength,particularly at low-confining pres-sures,which results in tensile cracking at a very low stress compared with compressive stresses.The tensile cracking reduces the stiffness of concrete structural components. Therefore,the use of continuum damage mechanics is necessary to accurately model the degradation in the mechanical properties of concrete.However,the concrete material undergoes also some irreversible(plastic)deformations during unloading such that the continuum damage theories cannot be used alone,particularly at high-confining pressures. Therefore,the nonlinear material behavior of concrete can be attributed to two distinct material mechanical processes:damage(micro-cracks,micro-cavities,nucleation and coa-lescence,decohesions,grain boundary cracks,and cleavage in regions of high stress con-centration)and plasticity,which its mechanism in concrete is not completely understood up-to-date.These two degradation phenomena may be described best by theories of con-tinuum damage mechanics and plasticity.Therefore,a model that accounts for both plas-ticity and damage is necessary.In this work,a coupled plastic-damage model is thus formulated.Plasticity theories have been used successfully in modeling the behavior of metals where the dominant mode of internal rearrangement is the slip process.Although the mathemat-ical theory of plasticity is thoroughly established,its potential usefulness for representing a wide variety of material behavior has not been yet fully explored.There are many research-ers who have used plasticity alone to characterize the concrete behavior(e.g.Chen and Chen,1975;William and Warnke,1975;Bazant,1978;Dragon and Mroz,1979;Schreyer, 1983;Chen and Buyukozturk,1985;Onate et al.,1988;Voyiadjis and Abu-Lebdeh,1994; Karabinis and Kiousis,1994;Este and Willam,1994;Menetrey and Willam,1995;Grassl et al.,2002).The main characteristic of these models is a plasticity yield surface that includes pressure sensitivity,path sensitivity,non-associativeflow rule,and work or strain hardening.However,these works failed to address the degradation of the material stiffness due to micro-cracking.On the other hand,others have used the continuum damage theory alone to model the material nonlinear behavior such that the mechanical effect of the pro-gressive micro-cracking and strain softening are represented by a set of internal state vari-ables which act on the elastic behavior(i.e.decrease of the stiffness)at the macroscopic level (e.g.Loland,1980;Ortiz and Popov,1982;Krajcinovic,1983,1985;Resende and Martin, 1984;Simo and Ju,1987a,b;Mazars and Pijaudier-Cabot,1989;Lubarda et al.,1994). However,there are several facets of concrete behavior(e.g.irreversible deformations, inelastic volumetric expansion in compression,and crack opening/closure effects)that can-not be represented by this method,just as plasticity,by itself,is insufficient.Since both micro-cracking and irreversible deformations are contributing to the nonlinear response of concrete,a constitutive model should address equally the two physically distinct modes of irreversible changes and should satisfy the basic postulates of thermodynamics.1876U.Cicekli et al./International Journal of Plasticity23(2007)1874–1900 Combinations of plasticity and damage are usually based on isotropic hardening com-bined with either isotropic(scalar)or anisotropic(tensor)damage.Isotropic damage is widely used due to its simplicity such that different types of combinations with plasticity models have been proposed in the literature.One type of combination relies on stress-based plasticity formulated in the effective(undamaged)space(e.g.Yazdani and Schreyer,1990; Lee and Fenves,1998;Gatuingt and Pijaudier-Cabot,2002;Jason et al.,2004;Wu et al., 2006),where the effective stress is defined as the average micro-scale stress acting on the undamaged material between micro-defects.Another type is based on stress-based plastic-ity in the nominal(damaged)stress space(e.g.Bazant and Kim,1979;Ortiz,1985;Lubliner et al.,1989;Imran and Pantazopoulu,2001;Ananiev and Ozbolt,2004;Kratzig and Poll-ing,2004;Menzel et al.,2005;Bru¨nig and Ricci,2005),where the nominal stress is defined as the macro-scale stress acting on both damaged and undamaged material.However,it is shown by Abu Al-Rub and Voyiadjis(2004)and Voyiadjis et al.(2003,2004)that coupled plastic-damage models formulated in the effective space are numerically more stable and attractive.On the other hand,for better characterization of the concrete damage behavior, anisotropic damage effects,i.e.different micro-cracking in different directions,should be characterized.However,anisotropic damage in concrete is complex and a combination with plasticity and the application to structural analysis is straightforward(e.g.Yazdani and Schreyer,1990;Abu-Lebdeh and Voyiadjis,1993;Voyiadjis and Kattan,1999;Carol et al.,2001;Hansen et al.,2001),and,therefore,it has been avoided by many authors.Consequently,with inspiration from all the previous works,a coupled anisotropic dam-age and plasticity constitutive model that can be used to predict the concrete distinct behavior in tension and compression is formulated here within the basic principles of ther-modynamics.The proposed model includes important aspects of the concrete nonlinear behavior.The model considers different responses of concrete under tension and compres-sion,the effect of stiffness degradation,and the stiffness recovery due to crack closure dur-ing cyclic loading.The yield criterion that has been proposed by Lubliner et al.(1989)and later modified by Lee and Fenves(1998)is adopted.Pertinent computational aspects con-cerning the algorithmic aspects and numerical implementation of the proposed constitu-tive model in the well-knownfinite element code ABAQUS(2003)are presented.Some numerical applications of the model to experimental tests of concrete specimens under dif-ferent uniaxial and biaxial tension and compression loadings are provided to validate and demonstrate the capability of the proposed model.2.Modeling anisotropic damage in concreteIn the current literature,damage in materials can be represented in many forms such as specific void and crack surfaces,specific crack and void volumes,the spacing between cracks or voids,scalar representation of damage,and general tensorial representation of damage.Generally,the physical interpretation of the damage variable is introduced as the specific damaged surface area(Kachonov,1958),where two cases are considered:iso-tropic(scalar)damage and anisotropic(tensor)damage density of micro-cracks and micro-voids.However,for accurate interpretation of damage in concrete,one should con-sider the anisotropic damage case.This is attributed to the evolution of micro-cracks in concrete whereas damage in metals can be satisfactorily represented by a scalar damage variable(isotropic damage)for evolution of voids.Therefore,for more reliable represen-tation of concrete damage anisotropic damage is considered in this study.The effective(undamaged)configuration is used in this study in formulating the damage constitutive equations.That is,the damaged material is modeled using the constitutive laws of the effective undamaged material in which the Cauchy stress tensor,r ij,can be replaced by the effective stress tensor, r ij(Cordebois and Sidoroff,1979;Murakami and Ohno,1981;Voyiadjis and Kattan,1999):r ij¼M ijkl r klð1Þwhere M ijkl is the fourth-order damage effect tensor that is used to make the stress tensor symmetrical.There are different definitions for the tensor M ijkl that could be used to sym-metrize r ij(see Voyiadjis and Park,1997;Voyiadjis and Kattan,1999).In this work the definition that is presented by Abu Al-Rub and Voyiadjis(2003)is adopted:M ijkl¼2½ðd ijÀu ijÞd klþd ijðd klÀu klÞ À1ð2Þwhere d ij is the Kronecker delta and u ij is the second-order damage tensor whose evolution will be defined later and it takes into consideration different evolution of damage in differ-ent directions.In the subsequence of this paper,the superimposed dash designates a var-iable in the undamaged configuration.The transformation from the effective(undamaged)configuration to the damaged one can be done by utilizing either the strain equivalence or strain energy equivalence hypoth-eses(see Voyiadjis and Kattan,1999).However,in this work the strain equivalence hypothesis is adopted for simplicity,which basically states that the strains in the damaged configuration and the strains in the undamaged(effective)configuration are equal.There-fore,the total strain tensor e ij is set equal to the corresponding effective tensor e ij(i.e.e ij¼ e ijÞ,which can be decomposed into an elastic strain e eij (= e eijÞand a plastic straine p ij(= e p ijÞsuch that:e ij¼e eij þe p ij¼ e eijþ e p ij¼ e ijð3ÞIt is noteworthy that the physical nature of plastic(irreversible)deformations in con-crete is not well-founded until now.Whereas the physical nature of plastic strain in metals is well-understood and can be attributed to the generation and motion of dislocations along slip planes.Therefore,in metals any additional permanent strains due to micro-cracking and void growth can be classified as a damage strain.These damage strains are shown by Abu Al-Rub and Voyiadjis(2003)and Voyiadjis et al.(2003,2004)to be minimal in metals and can be simply neglected.Therefore,the plastic strain in Eq.(3) incorporates all types of irreversible deformations whether they are due to tensile micro-cracking,breaking of internal bonds during shear loading,and/or compressive con-solidation during the collapse of the micro-porous structure of the cement matrix.In the current work,it is assumed that plasticity is due to damage evolution such that damage occurs before any plastic deformations.However,this assumption needs to be validated by conducting microscopic experimental characterization of concrete damage.Using the generalized Hook’s law,the effective stress is given as follows: r ij¼E ijkl e eklð4Þwhere E ijkl is the fourth-order undamaged elastic stiffness tensor.For isotropic linear-elas-tic materials,E ijkl is given byE ijkl¼2GI dijkl þKI ijklð5ÞU.Cicekli et al./International Journal of Plasticity23(2007)1874–19001877where I dijkl ¼I ijklÀ13d ij d kl is the deviatoric part of the fourth-order identity tensorI ijkl¼12ðd ik d jlþd il d jkÞ,and G¼E=2ð1þmÞand K¼E=3ð1À2mÞare the effective shearand bulk moduli,respectively,with E being the Young’s modulus and m is the Poisson’s ratio which are obtained from the stress–strain diagram in the effective configuration.Similarly,in the damaged configuration the stress–strain relationship in Eq.(4)can be expressed by:r ij¼E ijkl e eklð6Þsuch that one can express the elastic strain from Eqs.(4)and(5)by the following relation:e e ij ¼EÀ1ijklr kl¼EÀ1ijklr klð7Þwhere EÀ1ijkl is the inverse(or compliance tensor)of the fourth-order damaged elastic tensorE ijkl,which are a function of the damage variable u ij.By substituting Eq.(1)into Eq.(7),one can express the damaged elasticity tensor E ijkl in terms of the corresponding undamaged elasticity tensor E ijkl by the following relation:E ijkl¼MÀ1ijmnE mnklð8ÞMoreover,combining Eqs.(3)and(7),the total strain e ij can be written in the following form:e ij¼EÀ1ijkl r klþe p ij¼EÀ1ijklr klþe p ijð9ÞBy taking the time derivative of Eq.(3),the rate of the total strain,_e ij,can be written as _e ij¼_e eijþ_e p ijð10Þwhere_e eij and_e p ij are the rate of the elastic and plastic strain tensors,respectively.Analogous to Eq.(9),one can write the following relation in the effective configuration:_e ij¼EÀ1ijkl _ rklþ_e p ijð11ÞHowever,since E ijkl is a function of u ij,a similar relation as Eq.(11)cannot be used. Therefore,by taking the time derivative of Eq.(9),one can write_e ij in the damaged con-figuration as follows:_e ij¼EÀ1ijkl _r klþ_EÀ1ijklr klþ_e p ijð12ÞConcrete has distinct behavior in tension and compression.Therefore,in order to ade-quately characterize the damage in concrete due to tensile,compressive,and/or cyclic loadings the Cauchy stress tensor(nominal or effective)is decomposed into a positive and negative parts using the spectral decomposition technique(e.g.Simo and Ju, 1987a,b;Krajcinovic,1996).Hereafter,the superscripts‘‘+”and‘‘À”designate,respec-tively,tensile and compressive entities.Therefore,r ij and r ij can be decomposed as follows:r ij¼rþij þrÀij; r ij¼ rþijþ rÀijð13Þwhere rþij is the tension part and rÀijis the compression part of the stress state.The stress tensors rþij and rÀijcan be related to r ij byrþkl ¼Pþklpqr pqð14ÞrÀkl ¼½I klpqÀPþijpqr pq¼PÀklpqr pqð15Þ1878U.Cicekli et al./International Journal of Plasticity23(2007)1874–1900such that Pþijkl þPÀijkl¼I ijkl.The fourth-order projection tensors Pþijkland PÀijklare definedas follows:Pþijpq ¼X3k¼1Hð^rðkÞÞnðkÞi nðkÞj nðkÞpnðkÞq;PÀklpq¼I klpqÀPþijpqð16Þwhere Hð^ rðkÞÞdenotes the Heaviside step function computed at k th principal stress^rðkÞof r ij and nðkÞi is the k th corresponding unit principal direction.In the subsequent develop-ment,the superscript hat designates a principal value.Based on the decomposition in Eq.(13),one can assume that the expression in Eq.(1) to be valid for both tension and compression,however,with decoupled damage evolution in tension and compression such that:rþij ¼Mþijklrþkl; rÀij¼MÀijklrÀklð17Þwhere Mþijkl is the tensile damage effect tensor and MÀijklis the corresponding compressivedamage effect tensor which can be expressed using Eq.(2)in a decoupled form as a func-tion of the tensile and compressive damage variables,uþij and uÀij,respectively,as follows:Mþijkl ¼2½ðd ijÀuþijÞd klþd ijðd klÀuþklÞ À1;MÀijkl¼2½ðd ijÀuÀijÞd klþd ijðd klÀuÀklÞ À1ð18ÞNow,by substituting Eq.(17)into Eq.(13)2,one can express the effective stress tensor as the decomposition of the fourth-order damage effect tensor for tension and compression such that:r ij¼Mþijkl rþklþMÀijklrÀklð19ÞBy substituting Eqs.(14)and(15)into Eq.(19)and comparing the result with Eq.(1), one can obtain the following relation for the damage effect tensor such that:M ijpq¼Mþijkl PþklpqþMÀijklPÀklpqð20ÞUsing Eq.(16)2,the above equation can be rewritten as follows:M ijpq¼Mþijkl ÀMÀijklPþklpq þMÀijpqð21ÞOne should notice the following:M ijkl¼Mþijkl þMÀijklð22Þoru ij¼uþij þuÀijð23ÞIt is also noteworthy that the relation in Eq.(21)enhances a coupling between tensileand compressive damage through the fourth-order projection tensor Pþijkl .Moreover,forisotropic damage,Eq.(20)can be written as follows:M ijkl¼Pþijkl1ÀuþþPÀijkl1ÀuÀð24ÞIt can be concluded from the above expression that by adopting the decomposition of the scalar damage variable u into a positive u+part and a negative uÀpart still enhances adamage anisotropy through the spectral decomposition tensors Pþijkl and PÀijkl.However,this anisotropy is weak as compared to the anisotropic damage effect tensor presented in Eq.(21).U.Cicekli et al./International Journal of Plasticity23(2007)1874–190018793.Elasto-plastic-damage modelIn this section,the concrete plasticity yield criterion of Lubliner et al.(1989)which was later modified by Lee and Fenves(1998)is adopted for both monotonic and cyclic load-ings.The phenomenological concrete model of Lubliner et al.(1989)and Lee and Fenves (1998)is formulated based on isotropic(scalar)stiffness degradation.Moreover,this model adopts one loading surface that couples plasticity to isotropic damage through the effective plastic strain.However,in this work the model of Lee and Fenves(1998)is extended for anisotropic damage and by adopting three loading surfaces:one for plastic-ity,one for tensile damage,and one for compressive damage.The plasticity and the com-pressive damage loading surfaces are more dominate in case of shear loading and compressive crushing(i.e.modes II and III cracking)whereas the tensile damage loading surface is dominant in case of mode I cracking.The presentation in the following sections can be used for either isotropic or anisotropic damage since the second-order damage tensor u ij degenerates to the scalar damage vari-able in case of uniaxial loading.3.1.Uniaxial loadingIn the uniaxial loading,the elastic stiffness degradation variables are assumed asincreasing functions of the equivalent plastic strains eþeq and eÀeqwith eþeqbeing the tensileequivalent plastic strain and eÀeq being the compressive equivalent plastic strain.It shouldbe noted that the material behavior is controlled by both plasticity and damage so that, one cannot be considered without the other(see Fig.1).For uniaxial tensile and compressive loading, rþij and rÀijare given as(Lee and Fenves,1998)rþ¼ð1ÀuþÞE eþe¼ð1ÀuþÞEðeþÀeþpÞð25ÞrÀ¼ð1ÀuÀÞE eÀe¼ð1ÀuÀÞEðeÀÀeÀpÞð26ÞThe rate of the equivalent(effective)plastic strains in compression and tension,eÀep and eþep,are,respectively,given as follows in case of uniaxial loading:1880U.Cicekli et al./International Journal of Plasticity23(2007)1874–1900_eþeq ¼_e p11;_eÀeq¼À_e p11ð27Þsuch thateÀeq ¼Z t_eÀeqd t;eþeq¼Z t_eþeqd tð28ÞPropagation of cracks under uniaxial loading is in the transverse direction to the stress direction.Therefore,the nucleation and propagation of cracks cause a reduction of the capacity of the load-carrying area,which causes an increase in the effective stress.This has little effect during compressive loading since cracks run parallel to the loading direc-tion.However,under a large compressive stress which causes crushing of the material,the effective load-carrying area is also considerably reduced.This explains the distinct behav-ior of concrete in tension and compression as shown in Fig.2.It can be noted from Fig.2that during unloading from any point on the strain soften-ing path(i.e.post peak behavior)of the stress–strain curve,the material response seems to be weakened since the elastic stiffness of the material is degraded due to damage evolution. Furthermore,it can be noticed from Fig.2a and b that the degradation of the elastic stiff-ness of the material is much different in tension than in compression,which is more obvi-ous as the plastic strain increases.Therefore,for uniaxial loading,the damage variable can be presented by two independent damage variables u+and uÀ.Moreover,it can be noted that for tensile loading,damage and plasticity are initiated when the equivalent appliedstress reaches the uniaxial tensile strength fþ0as shown in Fig.2a whereas under compres-sive loading,damage is initiated earlier than plasticity.Once the equivalent applied stressreaches fÀ0(i.e.when nonlinear behavior starts)damage is initiated,whereas plasticityoccurs once fÀu is reached.Therefore,generally fþ¼fþufor tensile loading,but this isnot true for compressive loading(i.e.fÀ0¼fÀuÞ.However,one may obtain fÀ%fÀuin caseof ultra-high strength concrete.3.2.Multiaxial loadingThe evolution equations for the hardening variables are extended now to multiaxial loadings.The effective plastic strain for multiaxial loading is given as follows(Lubliner et al.,1989;Lee and Fenves,1998):U.Cicekli et al./International Journal of Plasticity23(2007)1874–19001881_e þeq ¼r ð^ r ij Þ^_e p maxð29Þ_e Àeq ¼Àð1Àr ð^ r ij ÞÞ^_e p min ð30Þwhere ^_e p max and ^_e p min are the maximum and minimum principal values of the plastic strain tensor _e p ij such that ^_e p 1>^_e p 2>^_e p 3where ^_e p max ¼^_e p 1and ^_e p min ¼^_ep 3.Eqs.(29)and (30)can be written in tensor format as follows:_j p i ¼H ij ^_e p jð31Þor equivalently _e þeq 0_e Àeq8><>:9>=>;¼H þ0000000H À264375^_e p 1^_e p 2^_e p 38><>:9>=>;ð32ÞwhereH þ¼r ð^ rij Þð33ÞH À¼Àð1Àr ð^ r ij ÞÞð34ÞThe dimensionless parameter r ð^ rij Þis a weight factor depending on principal stresses and is defined as follows (Lubliner et al.,1989):r ð^ r ij Þ¼P 3k¼1h ^ r k i P k ¼1j ^ r kj ð35Þwhere h i is the Macauley bracket,and presented as h x i ¼1ðj x j þx Þ,k ¼1;2;3.Note that r ð^ rij Þ¼r ð^r ij Þ.Moreover,depending on the value of r ð^r ij Þ,–in case of uniaxial tension ^ r k P 0and r ð^ r ij Þ¼1,–in case of uniaxial compression ^ rk 60and r ð^ r ij Þ¼03.3.Cyclic loadingIt is more difficult to address the concrete damage behavior under cyclic loading;i.e.transition from tension to compression or vise versa such that one would expect that under cyclic loading crack opening and closure may occur and,therefore,it is a challenging task to address such situations especially for anisotropic damage evolution.Experimentally,it is shown that under cyclic loading the material goes through some recovery of the elastic stiffness as the load changes sign during the loading process.This effect becomes more sig-nificant particularly when the load changes sign during the transition from tension to com-pression such that some tensile cracks tend to close and as a result elastic stiffness recovery occurs during compressive loading.However,in case of transition from compression to tension one may thus expect that smaller stiffness recovery or even no recovery at all may occur.This could be attributed to the fast opening of the pre-existing cracks that had formed during the previous tensile loading.These re-opened cracks along with the new cracks formed during the compression will cause further reduction of the elastic stiffness that the body had during the first transition from tension to compression.The1882U.Cicekli et al./International Journal of Plasticity 23(2007)1874–1900consideration of stiffness recovery effect due to crack opening/closing is therefore impor-tant in defining the concrete behavior under cyclic loading.Eq.(21)does not incorporate the elastic stiffness recovery phenomenon as well as it does not incorporate any coupling between tensile damage and compressive damage and,therefore,the formulation of Lee and Fenves(1998)for cyclic loading is extended here for the anisotropic damage case.Lee and Fenves(1998)defined the following isotropic damage relation that couples both tension and compression effects as well as the elastic stiffness recovery during transi-tion from tension to compression loading such that:u¼1Àð1Às uþÞð1ÀuÀÞð36Þwhere sð06s61Þis a function of stress state and is defined as follows: sð^ r ijÞ¼s0þð1Às0Þrð^ r ijÞð37Þwhere06s061is a constant.Any value between zero and one results in partial recovery of the elastic stiffness.Based on Eqs.(36)and(37):(a)when all principal stresses are positive then r=1and s=1such that Eq.(36)becomesu¼1Àð1ÀuþÞð1ÀuÀÞð38Þwhich implies no stiffness recovery during the transition from compression to tension since s is absent.(b)when all principal stresses are negative then r=0and s¼s0such that Eq.(36)becomesu¼1Àð1Às0uþÞð1ÀuÀÞð39Þwhich implies full elastic stiffness recovery when s0¼0and no recovery when s0¼1.In the following two approaches are proposed for extending the Lee and Fenves(1998)model to the anisotropic damage case.Thefirst approach is by multiplying uþij in Eq.(18)1by the stiffness recovery factor s:Mþijkl ¼2½ðd ijÀs uþijÞd klþd ijðd klÀs uþklÞ À1ð40Þsuch that the above expression replaces Mþijkl in Eq.(21)to give the total damage effecttensor.Another approach to enhance coupling between tensile damage and compressive dam-age as well as in order to incorporate the elastic stiffness recovery during cyclic loading for the anisotropic damage case is by rewriting Eq.(36)in a tensor format as follows:u ij¼d ijÀðd ikÀs uþik Þðd jkÀuÀjkÞð41Þwhich can be substituted back into Eq.(2)to get thefinal form of the damage effect tensor, which is shown next.It is noteworthy that in case of full elastic stiffness recovery(i.e.s=0),Eq.(41)reducesto u ij¼uÀij and in case of no stiffness recovery(i.e.s=1),Eq.(41)takes the form ofu ij ¼uÀijþuþikÀuþikuÀjksuch that both uþijand uÀijare coupled.This means that duringthe transition from tension to compression some cracks are closed or partially closed which could result in partial recovery of the material stiffness(i.e.s>0)in the absence U.Cicekli et al./International Journal of Plasticity23(2007)1874–19001883。
