Simplified damage models applied in the numerical analysis of reinforced concrete structures

advantages of these modeling. Initially, it is shortly described the damage model proposed by Mazars. This constitutive model assumes the concrete as isotropic and elastic material, where locally the damage is due to extensions. On the other hand, the damage model proposed by Pituba, the material is assumed as initial elastic isotropic medium presenting anisotropy, plastic strains and bimodular response (distinct elastic responses whether tension or compression stress states prevail) induced by the damage. To take into account for bimodularity two damage tensors governing the rigidity in tension and compression regimes, respectively, are introduced. Damage activation is expressed by two criteria indicating the initial and further evolution of damage. Soon after, the models are used in numerical analyses of the mechanical behavior of reinforced concrete structures. Accordingly with comparison of the obtained responses, considerations about the application of the isotropic and anisotropic damage models are presented for 1D and 2D reinforced concrete structures modeling as well as the potentialities of the simplified versions of damage models applied in situations of structural engineering.


Introduction
In the Continuum Damage Mechanics (CDM), the damage effects are evidenced in the stiffness constitutive tensor.The damage leads to the reduction of several stiffness components, where the damaged material can either keep its isotropic properties or to become anisotropic.For isotropic models the damage affects neither the direction nor the initial number of symmetry planes presented by the material.Thus, it does not matter if the medium is initially isotropic or anisotropic with some degrees.In this case, those initial characteristics are preserved during the damage process.Some constitutive damage models have been proposed assuming the concrete as an isotropic medium (Mazars [1], Mazars [2], Comi [3] and Berthaud [4]).However, in the last decades anisotropic models which can modify both the direction and the number of the material symmetry planes have been proposed (Brünig [5], Pituba [6], Pietruszczak [7], Ibrahimbegovic [8] and Dragon [9]).Besides, another important characteristic presented by many fiberreinforced composite materials is the intrinsic bimodularity, i.e., distinct responses in tension and compression prevailing states.On the other hand, brittle materials, such as concrete, are a kind of composites that can be initially considered isotropic and unimodular.However when they have been damaged, those materials would start to present some degree of anisotropy and bimodularity.Assuming small deformations, a formulation of constitutive laws for either initially isotropic or anisotropic elastic bimodular materials was proposed by Curnier [10].In order to incorporate damage effects, the formulation of Curnier has been extended by Pituba [11].In particular, a constitutive model for concrete has been derived.Accordingly, the material is initially considered as an isotropic continuous medium with anisotropy and bimodularity induced by the damage.On one side the class of anisotropy induced and considered in the model (transversal isotropy) elapses from the assumption that locally the loaded concrete always presents a diffuse oriented damage distribution as appointed by experimental observations (Van Mier [12]).On the other hand, the bimodularity induced by damage is captured by the definition of two damage tensors: one for dominant tension states and another one for dominant compression states.In order to use the Damage Mechanics in practical situations of Structural Engineering, constitutive models presenting a reduced number of parameters with easy identification and one-dimensional version are desirable.On the other hand, such models must present reliable numerical results in order to estimate the mechanical behavior of the structure as accurate as possible.In this work, two examples of this kind of the damage models are used and here, they are called simplified damage models.This work intends to discuss the related problems to the numerical applications of the isotropic and anisotropic damage models in the context of the one and two-dimensional analyses of reinforced concrete structures.Besides, one intends to show the potentialities of an anisotropic damage model recently proposed by Pituba [6].Then, numerical responses supplied by the models are presented and compared in order to evidence the difficulties and advantages when one deals with this kind of modeling.Finally, some conclusions about the employment of the simplified versions of these kinds of damage models are discussed.

Isotropic Damage Model
This model has been proposed by Mazars [1] and the damage is represented by the scalar variable D (with 0 ≤ D ≤ 1) whose evolution occurs when the equivalent extension deformation ε ~ is bigger than a reference value.The plastic deformations evidenced experimentally are not considered.The equivalent extension deformation is given by: (1) where i ε is a principal deformation component, being The damage activation occurs when = ε ~ed0 , being e d0 the defor- mation referred to the maximum stress of an uniaxial tension test.Thus the criterion is given by: (2) Considering the thermodynamics principles, the damage evolution can be expressed by: (3a) As the concrete behaves differently in tension and compression, the damage variable D is obtained by combining properly the variables D T and D C , related to tension and compression, respectively, as follows: (4) where 1 C T

= + a a
where D T and D C are given by: e ----= Simplified damage models applied in the numerical analysis of reinforced concrete structures and compression when damaged.This model has been proposed by Pituba [6] and it follows the from the formalism presented in Pituba [11].Moreover, the model respects the principle of energy equivalence between damaged real medium and equivalent continuous medium established in the CDM (Lemaitre [13]).Now, the damage model is shortly described.Initially, for dominant tension states, a damage tensor is given by: (9) The variable D 1 represents the damage in direction orthogonal to the transverse isotropy local plane of the material, while D 4 is representative of the damage due to the sliding movement between the crack faces.The third damage variable, D 5 , is only activated if a previous compression state accompanied by damage has occurred.In the Eq. ( 9), the tensor I is the second-order identity tensor and the tensor A, by definition, Curnier [10], is formed by dyadic product of the unit vector perpendicular to the transverse isotropy plane for himself.On the other hand, for dominant compression states, it is proposed the other damage tensor: (10) where f 1 (D 2 , D 4 , D 5 ) = D 2 -2 f 3 (D 4 , D 5 ) ,f 2 (D 3 ) = D 3 and f 3 (D 4 , D 5 )= 1 -(1-D 4 ) (1-D 5 ).Note that the compression damage tensor introduces two additional scalar variables in its composition: D 2 and D 3 .The variable D 2 (damage perpendicular to the transverse isotropy local plane of the material) reduces the Young's modulus in that direction and in conjunction to D 3 (that represents the damage in the transverse isotropy plane) degrades the Poisson's ratio throughout the perpendicular planes to the one of transverse isotropy.Finally, the constitutive tensors are written as follows: (11a) In Eqs (5) A T and B T are parameters related to uniaxial tension tests while A C and B C are obtained from uniaxial compression tests.To compute the a T and a C values defined in Eq.(4), we have to obtain, initially, the deformations e T and e C associated, respectively, to tension and compression states as follows: (6a) where I is the identity tensor, E the elastic modulus of a non-damaged material, given by: Finally, the constitutive relation can be expressed in terms of the actual deformation tensor as follows:

Anisotropic Damage Model
In this model, it is assumed that the concrete belongs to a category of materials that can be considered initially isotropic and unimodular, however they start to present different behaviours in tension damage induces anisotropy in the concrete.Therefore, it is convenient to separate the damage criteria into two: the first one is only used to indicate damage beginning, or that the material is no longer isotropic and the second one is used for loading and unloading when the material is already considered as transverse isotropic.This second criterion identifies if there is or not evolution of the damage variables.That division is justified by the difference between the complementary elastic strain energies of isotropic and transverse isotropic material.For identifying the damage beginning it is suggested a criterion that compares the complementary elastic strain energy * e W , which is computed locally considering the medium as initially virgin, isotropic and purely elastic, with a certain reference value Y 0T , or Y 0C , obtained from experimental tests of uniaxial tension, or compression, respectively.Accordingly, the criterion for initial activation of damage processes in tension or compression is given by: It is important to note that the damaged medium presents a transverse isotropy plane in correspondence to the current damage level.Then, the complementary elastic energy of the damaged medium is expressed in different forms, depending on whether tension or compression strain states prevail.In the case of dominant tension states (g(e,D T ,D C ) > 0) assuming that direction 1 in the strain space be perpendicular to the transverse isotropy local plane, it can be written: (17) For the damaged medium in dominant compression states, the relationships are similar to the tension case, where the complementary elastic energy is expressed in the following form: (18) in that way the anisotropy and bimodularity induced by damage.Those parameters are given by: The different dyadic products in Eqs. ( 8), ( 10) and ( 11) have the function of allocating the material constants in certain positions of the stiffness constitutive tensors.For more details see Curnier [10] and Pituba [11].Observe that the bimodular character is taken into account by the conditions g (ε, D T , D C ) > 0 or g (ε, D T , D C ) < 0, where g (ε, D T , D C ) is a hypersurface that contains the origin and divides the strain space into a compression and tension sub-domains.A particular form is adopted for the hypersurface in the strain space: a hyperplane g(ε,D) defined by the unit normal N and characterized by its dependence of both the strain and damage states.To simplify the presentation, the hyperplane will be here expressed as the one obtained by enforcing the direction 1 in the strain space to be perpendicular to the transverse isotropy local plane.Thus, the hyperplane is given by: The h(D 1 ) and h(D 2 ) functions are defined, respectively, for the tension and compression cases, assuming for the first one that there was no previous damage of compression affecting the present damage variable D 1 and analogously, for the second one that has not had previous damage of tension affecting variable D 2 .The functions can be written as: (15) As it has already been pointed out, in the model formulation the

Simplified damage models applied in the numerical analysis of reinforced concrete structures
Assuming a general situation of damaged medium for dominant tension states, the criterion for the identification of damage increments is represented by the following relationship: where the reference value * T Y 0 is defined by the maximum complementary elastic energy computed throughout the damage process up to the current state.Analogous expressions are valid for dominant compression states.In the loading case, i. e., when , it is necessary to update the values of the scalar damage variables that appear in the D T and D C tensors, considering their evolution laws.Considering just the case of monotonic loading, the evolution laws proposed for the scalar damage variables are resulting of fittings on experimental results and present similar characteristics to those one described in Mazars [1] and Berthaud [4] works.The general form proposed is: where A i , B i and Y 0i are parameters that must be identified.The parameters Y 0i are understood as initial limits for the damage activation, the same ones used in Eq. ( 16).When the damage process is activated, the formulation starts to involve the tensor A that depends on the normal to the transverse isotropy plane.Therefore, it is necessary to establish some rules to identify its location for an actual strain state.Initially, it is established a general criterion for the existence of the transverse isotropy plane.It is proposed that the transverse isotropy due to damage only arises if positive strain rates exist at least in one of the principal directions, Pituba [6].After assuming such proposition as valid, some rules to identify its location must be defined.First of all, considering a strain state in which one of the strain rates is no-null or has sign contrary to the others, the following rule is applied: "In the principal strain space, if two of the three strain rates are extension, shortening or null, the plane defined by them will be the transverse isotropy local plane of the material." The uniaxial tension is an example of this case where the transverse isotropy plane is perpendicular to the tension stress direction.Obviously, it can be suggested criteria based on others formulations, such as for instance, the microplanes theory developed by Bazant [14].

Plain Concrete Beam
In this first numerical application, previously performed by Guello [15], we have considered the concrete beam without any reinforcement bar described in Note that in the one-dimensional analysis the difference between numerical responses of damage models is mainly due to more aggressive reduction of the stiffness in the anisotropic model.While the Mazars' model presents a linear reduction way of the material stiffness (Strain Equivalence -Eq.(8)), the Pituba's model reduces the material stiffness in a quadratic way (Energy Equivalence -Eq.( 11)).In this analysis, it is possible to reproduce a stiffness break about 45 kN.Soon after, a subsequent strength recovery of the structural element is illustrated by the models.Then, the difference between the numerical responses is increased.Finally, note that in 1D and 2D numerical analyses, both models present the same qualitative behaviors.
In order to visualize the damage distribution in the beam, Figure [6] presents the isodamage curves for Mazars´model, Proença [16].Figure [7] shows the damage distribution for 1D and 2D numerical analyses with Pituba´s model considering a loading stage about 50 kN.These curves have been obtained by interpolation of the damage variables along the adopted integration points.
In the context of 1D analysis, when the Proença's response is compared with the Pituba´s model, the damage configuration is very similar, mainly due to damage in tension.Note that Mazars' model presents only one damage variable representing a combination of tension and compression damage processes.For the concrete structures, both models assume that locally the damage is due to extensions.On the other hand, the 2D numerical responses show a more intense damage process in compression given by D 2 and D 3 variables.In fact, the variable D 3 tries to simulate a crushing process of the concrete near of the supports and in upper middle zone of the beam.

Reinforced Concrete Bar Structure
The second numerical simulation deals with a reinforced concrete bar structure.Figure [8] shows the geometry of the structure.A reinforcement bar with 10.0 mm is placed in the centre of the section where a force P is applied.This numerical application was performed previously by Mazars [1].Table [2] presents the parameter values of each damage model.These parameters were obtained through procedure similar to the previous example with E c =30000 MPa, where the parameters suggested for the Mazars' model were taken as reference values, see Figure [9].For the steel, E a = 210 GPa has been adopted.Simplified damage models applied in the numerical analysis of reinforced concrete structures In the 1D analysis, a mesh with 20 elements and 10 layers has been considered.On the other hand, in the 2D analyses performed by Guello [15] with Mazars' model, a mesh with 38 x 16 triangle finite elements (6 nodes) has been used in order to represent the concrete, while a mesh with 30 bar elements represents the steel.On the other hand, in this work, a mesh with 100 constant strain quadrilateral (4 nodes) elements divided into 10 layers of 10 elements and placed in whole extension of the structural element has been considered for the 2D analysis with Pituba's model.Note that, in this mesh, one layer represents the reinforcement bar, see In the 1D analysis, once again the stiffness degradation is observed in a more evident way in the anisotropic model by the reasons previ-ously explained.However, that difference is more evident after the adherence loss stage between the reinforcement bar and concrete presented in the Fig. [10] about 28 kN.Note that the 1D models present a stiffness recovery after to total concrete cracking.That new resistance is just owed to the reinforcement bar that presents an elastic behavior up to yielding resistance be achieved.Those observations are in agreement with others works (Mazars [1] and Guello [15]).
On the other hand, in the 2D analysis, the anisotropic model has presented locking problems related to numerical responses due to the distortion of the finite element that should be partly responsible for the excessive stiffness, what does not happen with the Mazars' model (triangular element).Therefore, there was not possibility of more intense stiffness degradation by the anisotropic model.Note that in this analysis the elastic behavior for reinforcement bar has been assumed leading to the increase of difference about numerical responses given by the models, due mainly to the finite elements.Besides, in the 2D analysis the strain localization is an important phenomenon in the behavior of concrete structures.In particular, this structure is mainly tensioned presenting a localized cracking configuration [15] .

Reinforced Concrete Beam
The third numerical application is about a reinforced concrete beam submitted to monotonic loading.This beam was tested by Delalibera [17].The elastic parameters of the concrete are f ck =25MPa and E c =32.3MPa.For the reinforcement has been adopted E a = 205 GPa, yielding stress 590 MPa and ultimate stress 750 MPa.
The geometric characteristics of the beam are given in Figure [11].
The loading is composed by two equal forces applied on the beam.The uniaxial stress tests that were performed by Delalibera [17] and they were used in order to identify the models parameters, see  The numerical responses obtained with isotropic and anisotropic damage models are able to simulate the experimental behavior of the beam, see Figure [13].Both damage models present a strong loss of stiffness about 18 kN trying to evidence a possible damage localization process.In this context, in Delalibera [17] is reported that the first crack occurs about 25 kN.As it can be seen, the numerical responses are quite satisfactory when compared to the experimental results since the beginning of the damage process up to the complete rupture of the beam.On the other hand, it is noted that the numerical responses present a difficulty of convergence evidenced by the increase of the iterations about 115 kN.In this stage, the concrete presents high values for the damage variables in many layers of finite elements located in the medium zone of the beam, therefore the beam stiffness is mainly due to the reinforcement bars.

Conclusions
The study concerns to the employment of anisotropic and isotropic damage models to one and two-dimensional concrete structure analyses.In a general way, the numerical results obtained from the damage models presented in this work have been quite satisfactory.The potentialities of the Damage Mechanics when deals with numerical simulation of the non-linear behavior of concrete structures are shown.In particular, the employment of anisotropic models has shown some advantages in 2D analysis when Simplified damage models applied in the numerical analysis of reinforced concrete structures compared to the isotropic ones, such as, the selective stiffness deterioration and evolution of cracking configuration supplying a more realistic numerical response, see Pituba [6].This feature can be more evident in three-dimensional analyses.However, it must be observed that structures with low reinforcement rates can evidence some numerical problems due to plane analysis, see Pituba [6], Pituba [18] and Comi [19].In these cases, the cracking process starts to present a localized distribution limiting the employment of the damage models.In order to overcome numerical problems a non-local version of the anisotropic model can be proposed and implemented in a computational code, for instance, with so-called Generalized Finite Element Method.
On the other hand, it is important to observe that the proposal and parametric identification of evolution laws for damage variables D 4 and D 5 must increase the accuracy of the anisotropic model.In fact, these cracking processes related to shear behavior of the concrete are significant contributions to the released energy.This feature has been studied by Pituba [20] and a theoretical analysis has shown that the anisotropic model has advantage upon constitutive models that use the so called "shear retention factor".
The 1D analysis has shown an efficient and practical employment, without numerical problems and low computational cost.Besides, the parametric identification is simple.In this case, the anisotropic or isotropic damage models could be used in estimative analyses of structures in practical situations, such as: numerical simulation of displacement in cracking concrete beams in order to propose an alternative procedure to the Brazilian Technical Code (Pituba [21]), estimative of ultimate load capacity of frames and beams and collapse configuration of reinforced concrete frames (Pituba [22] and Pituba [23]), included the numerical analyses of the structures submitted to cyclic loading (Pituba [24]) .Finally, this work has demonstrated that simplified damage models are a good alternative to estimate the mechanical behavior of reinforced concrete structures.

.
= dD/dt , i. e., D time derivative; ( ) ε F is written in terms of ε ~ and defined continuous and positive.

are
Lamè constants.The remaining parameters will only exist for no-null damage, evidencing

D
is the elastic fourth order tensor of the non-damaged material.Thus the coefficients T α and C α are obtained by the following expression: , positive and negative parts of the deformations e T and e C defined in Eq. (6); f T,C (s) = * e W -Y 0T,0C < 0 then D T = 0 (i.e., D 1 = D 4 = 0) for dominant tension states or D C = 0 (i.e., D 2 = D 3 = D 5 = 0) for dominant compression states, where the material is linear elastic and isotropic.The reference values Y 0T and Y 0C are model parameters defined by σ 0T and σ 0C are the limit elastic stresses determined in the uniaxial tension and compression regimes.

( 13 ) 11 e
g(e,D T ,D C ) = N(D T ,D C ) .e e = g 1 (D 1 ,D 2 ) e V e + g 2 (D 1 ,D 2 ) ewhere g 1 (D 1 ,D 2 ) = {1+H(D 2 )[H(D 1 )-1]}h(D 1 )+{1+H(D 1 )[H(D 2 )-1]} h(D 2 ) and g 2 (D 1 ,D 2 ) = D 1 +D 2 .The Heaveside functions employed above are given by:(14) The parametric identification of the model is accomplished by uniaxial tension tests in order to obtain A 1 , B 1 and Y 01 = Y 0T , by uniaxial compression tests for the identification of the parameters A 2 , B 2 and Y 02 , and finally by biaxial compression tests in order to obtain A 3 , B 3 and Y 03 = Y 02 = Y 0C .On the other hand, the identification of the parameters for the evolution laws corresponding to the damage variables D 4 and D 5 , which influence the shear concrete behavior, it won't be studied in this version of the model because the experimental tests are not available yet to allow the parameter calibration or, even, the proposition of more realistic evolution laws.

Figure 1 -Figure 2 -
Figure 1 -Geometry details of the plain concrete beam

Figure 8 -Figure 9 -
Figure 8 -Geometry of the reinforced concrete bar structure

Figure 12 -Figure 13 -
Figure 12 -Parametric identification in uniaxial compression and tension tests for reinforced concrete beam

Table 2 -
Parameter values -reinforced concrete bar structure