Anisotropic damage model on the effects of damage process due to shearing stress in concrete

This study deals with the proposal of damage variables related to structural behavior of the concrete subjected to shear stress. It is also addressed the laws of evolution and their incorporation into the formulation of a constitutive model for the concrete. Originally, this model considers the anisotropy and bimodularity induced by damage. Initially, the formulation of the damage model is briefly described. The material is assumed to behave as an elastic isotropic medium that start to present anisotropy and bimodular response (distinct elastic responses for prevailing states of tension or compression). Afterwards, we discussed some proposals of damage variables related to the structural behavior of the concrete when subjected to shear stress and their respective laws of evolution. The chosen proposal is then incorporated to the original model and implemented in a finite element code to perform plane analysis. Finally, the numerical responses of the model are compared with the experimental of reinforced concrete beams with different reinforcement rates. The results are presented and discussed aiming to verify the improvements presented by the new version of the model and to show the contribution of the energy dissipation due to shear in the cracked concrete behavior for the analyzed examples.


Introduction
Physical nonlinear response of solids is a macroscopic evidence of irreversible changes in the microstructure. Damage mechanics simulate in equivalent continuous media, the deterioration of materials exclusively owing the microcracking process. Such media are considered intact, but with reduced elastic properties in order to simulate the real medium cracked, being possible to apply the formulations of the Continuous Mechanics. This approach is used for obtaining constitutive relationships through the equivalences of stress, deformation or energy responses.
Whereas the Fracture Mechanics deals with the propagation conditions of macroscopic cracks, the Damage Mechanics studies the effects of microcracks on the material response. The Damage Theory describes locally the phenomena progress which develop between an initial state, relative to a situation of intact material, and at a final state, represented by the formation of a macroscopic crack or, in other words, the rupture of the 'representative volume element' surrounding the considered point.
In the case of concrete, where cracking is the dominant phenomenon in the nonlinear behavior, Damage Mechanics is undoubtedly able to formulate realistic models. Within the Damage Mechanics scope, damage effects are included in the constitutive stiffness tensor. The damage implies a reduction of several stiffness components, and according to the model, the damaged material is capable of being isotropic or anisotropic.
In recent years, many constitutive models for ductile and fragile media have been proposed regarding their characteristics of anisotropy. Among them, Brünig (2004) has developed constitutive models to study the damage process in ductile media; Ibrahimbegovic et al. (2008) made the proposition of a damage-plasticity coupled model to analyze materials with irreversible plastic deformation, changes in elastic response and localized failure, being this last frequently studied with the objective of ensuring the numerical response objectivity HERVÉ et al., 2005;KUCEROVA et al., 2009). Besides, the micromechanics theory has been used for justifying the proposition of constitutive models that deal with cracked media (PICHLER; DORMIEUX, 2009;WELEMANE;CORMERY, 2002;ZHU et al., 2008ZHU et al., , 2009). More recently, given the scope of multi-scale analysis, models based on Damage Mechanics and cohesive fracture are used in analysis on representative volume elements (RVE) and concrete structures (NGUYEN et al., 2012;SKARZYNSKI;TEJCHMAN, 2012).
On the other hand, in the context of small deformations regime, the development of constitutive laws for materials with elastic isotropy or anisotropy which have different behaviors for tension and compression is presented by Curnier et al. (1995). In turn Pituba (2006) modified the formulation of Curnier et al. (1995) formulation in order to consider the damage. From this, the constitutive model for concrete was originated, which can be simulated as a continuous medium initially isotropic with anisotropy induced by damage. The class of induced anisotropy, adopted in the model, is based on the assumption that locally the requested concrete always presents a damage distribution with a well-defined orientation. This assumption, by the way, is justified by the observation of material behavior in experimental tests.
The damage model generated so far did not take into account the concrete behavior when subjected to shear stress. However, this originally proposed model has been tested under several situations evidencing good results (PITUBA, 2006(PITUBA, , 2008(PITUBA, , 2009(PITUBA, , 2010PITUBA;FERNANDES, 2011;PITUBA, 2003). The present study aimed at proposing a form of penalty to the modules of transversal elasticity of the cracked anisotropic medium without considering the friction between the faces of closed microcracks. Besides, we investigated how this proposal interferes with the damage model response originally proposed.

Damage model
Herein, concrete is assumed to be a material which belongs to the category of media initially isotropic which begins to present transversal isotropy and bimodular response induced by damage. The model formulation for concrete is based on the energy equivalence principle and the formalism presented in Pituba (2006).
As seen below, the proposed model is briefly described, beginning with the presentation of damage tensor for predominant states of tension states, whose expression is given by: where: f 1 (D 1 , D 4 , D 5 ) = D 1 -2 f 2 (D 4 , D 5 ); f 2 (D 4 , D 5 ) = 1 -(1-D 4 ) (1-D 5 ). Damage tensor presents two scalar variables in its composition (D 1 and D 4 ) and a third variable of damage D 5 , activated only under previous compression with corresponding damage. The variable D 1 is the damage in the direction perpendicular to the local plane of the material transversal isotropy and D 4 is a representative variable of the damage generated by the shearing between cracks borders belonging to that plane.
In the Eq. (1), tensor I is the identity tensor of second order and tensor A is, by definition, Curnier et al. (1995), formed by the tensorial product of the tensor perpendicular to the transversal isotropy plane. The tensorial product operations between tensors of second order I and A in Eq. (1) and which will be used throughout the formulation are described in Pituba (2006).
There are three scalar variables in its composition: D 2 , D 3 and D 5 , besides D 4 , related to pre-existing tensile effects. The variable D 2 (damage perpendicular to the transverse isotropy local plane) penalizes the elasticity module in this direction, and in conjunction with D 3 (representative of damage on the transverse isotropy local plane) penalizes the Poisson´s rate on planes perpendicular to the transversal isotropy.
Finally, the resultant constitutive tensors are described by: where: The other parameters only exist for non-null damage, demonstrating thus anisotropy and bimodularity induced by damage, and defined as: In Curnier et al. (1995) a hypersurface is defined in the space of tensions or deformations, to be utilized as criterion for the identification of constitutive responses of compression or tension.
In this model it is adopted a particular form to the hypersurface in the space of deformations: a hyperplane g (ε), characterized by its unitary normal N (||N|| = 1). The following equation is applied to the uniaxial version of the model: The Heaveside functions applied to the last relationship are given by: The functions η(D 1 ) and η(D 2 ) are defined, respectively, for cases of tension, assuming no previous damage of compression, and of compression, supposing no previous tensile damage.
Regarding the damage criterion, it is convenient to separate it into criterion to start damage, when material cease to be isotropic; and criterion for loading and unloading, herein considered in a sense of evolution or not of the damage variables, when the material is already transversally isotropic.
The criterion for initial activation of damage in tension or compression is given by: where * e W is the complementary elastic deformation energy considering the medium initially complete, isotropic and purely elastic and is a reference value obtained from uniaxial tests of tension, or compression, respectively, being As for the case of g(ε,D T ,D C ) > 0 and g(ε,D T ,D C ) < 0 , the complementary elastic energies of the damaged medium are given, respectively, by: Regarding a general situation of medium damaged in a predominant tensile regime, the criterion for identification of damage increments is represented by: where the reference value * T Y 0 is defined by the maximum complementary elastic energy determined throughout the damage process until the current state. As for the medium damaged in a predominant regime of compression, analogous relationships can be applied to the case of traction.
Given the loading, i.e. where , it is necessary to update the damage scalar variables values which appear in the tensors D T and D C , considering their laws of evolution.
Focusing on the proportional monotonic loading and uniaxial version of the model, the laws of evolution proposed for the damage scalar variables result from adjustments on the experimental results and present characteristics similar to those found in models of Mazars (1986) and La Borderie et al. (1994). The general form proposed is: where A i , B i and Y 0i are parameters to be identified. Parameters Y 0i are understood as initial limits for damage activation, the same used in Eq. (9). The parametric identification of the model is carried out through experimental tests of uniaxial tensile for obtaining the parameters A 1 , B 1, and Y 01 = Y 0T; of uniaxial compression for the identification of parameters A 2 , B 2 and Y 02 and finally of biaxial compression for obtaining A 3 , B 3 and Y 03 = Y 02 = Y 0C . On the other hand, identification of parameters included in the laws of evolution corresponding to the damage variables D 4 and D 5 , which influence the concrete behavior under shearing, are object of study in this research, where their laws of evolution are proposed. Results of parametric identification through experimental tests corresponding to uni, bi and triaxial states in specimens are described in Proença and Pituba (2003).
In turn, when the damage process is activated, the formulation involves the tensor A which depends on the knowledge of the normal to the transversal isotropy plane. Then, it must be defined rules to locate this normal for a given state of deformation. Pituba (2006) suggested a criterion to locate the transversal isotropy plane based on the signs of the principal strain rates. Initially it is proposed that the transversal isotropy from damage manifests only with positive rates of deformation, at least in one of the main directions. Next, it is defined some rules to identify its location. First, a state of deformation in which one of the deformation rates is non-null or of opposite with others, it is applied the following rule: "In the principal strain space, if two of three strain rates are of elongation, shortening or null, the plane defined by them will be the plane of isotropy".
In this case it is fitted for instance, the uniaxial tension case, resulting in a transversal isotropy plane perpendicular to the direction of tensile stress. However, there are cases not covered by this rule. For instance, the state of plane deformation in which the non-null strains are opposite. For this situation a second rule is applied: "In a state of plane strain, where the principal strain rates in the plane have opposite signs, the transverse isotropy local plane for the material is defined by the directions of the permanently null principal strain and the strain whose rate is positive".
Other particular case occurs when all rates of principal strain are positive. For these, a third rule is valid, where it is assumed that the direction of larger elongation is perpendicular to the transverse isotropy local plane of the material. Obviously, criteria based on other foundations may be suggested, for instance, the theory of microplanes proposed by Bazant and Ozbolt (1990) or models of microplanes used for the case of dynamic loading (WANG et al., 2011).
A very interesting and current theme deals with the deterioration models of the material stiffness, simulating the mechanical behavior of the microstructure, where the anisotropy presented in the macrostructure is a consequence from the microstructural behavior (NGUYEN et al., 2012;SKARZYNSKI;TEJCHMAN, 2012).

Considerations on the concrete shearing
In this section are presented some suggestions of laws of evolution for the variables of damage connected to the structural behavior of concrete subjected to shear stress. Besides, it is regarded here the issue related to the ease of parametric identification. This last, after analyses, motivated the choice for a simple formulation to incorporate the effects of the damage processes related to the concrete shear in the model originally proposed in Pituba (2006). In order to introduce the concrete shear phenomenon, a brief report is given on the cracked concrete behavior when subjected to shear stress. Nogueira et al. (2011) had already stated about the experimental investigations that showed the important contribution of the cracked concrete on the shear strength, in the case of reinforced concrete elements. Two mechanisms are presented in this behavior: the so-called dowel action and the transference of shear stress between the faces of cracks edges.
While the first is specific for reinforced concrete, the second is of major interest in the context of the present study.
Crack surfaces are usually coarse and uneven, generally bypassing large aggregates of concrete. The rough surfaces of the crack can transmit shear stress through friction and aggregates meshing. Application of v shear strength causes a displacement d and, due to the uneven surface causes a separation of crack surfaces or, a dilatation of the medium w. Experimental results showed that the crack initial opening (w 0 ) is the main variable affecting the aggregates meshing action. Large initial openings of cracks imply in lower stiffness and strength to shear. Some approximations to consider the stiffness loss to shear can be found in literature. Among them, the proposal of Hand et al. (1973) modeled the shear stiffness of the concrete (transversal elasticity module) as a reduced constant value in relation to the initial value. However, for Kara and Dundai (2012), the transversal elasticity module decreasing hyperbolically with the normal deformation to the cracked plane (representing the crack opening) is a better approximation in the confrontation with experimental results, being thus reasonable to propose. The major problem found in the proposition of a law of evolution and in the correspondent parametric identification is the difficult to perform and the obtaining of reliable values in experimental tests where concrete is subjected to shear stress.
Other phenomenon involved in the concrete shear is the joint dissipation arising from the evolving damage and the internal friction on the cracks surfaces. In Pituba (2008) is cited that the internal friction is related to the assumption that cracks do not simply open by separation forming faces with smooth surfaces, but tend to combine modes of separation and sliding according to surfaces with certain coarseness. The sliding between crack surfaces originates a strain by sliding, assumed as responsible for a plastic behavior with nonlinear kinematic hardening of the medium damaged. Araújo and Proença (2008), to consider the damage and internal friction, the full stress on a point of the damaged medium is divided into a portion called elastic-damaged stress and other named sliding stress, also dependent on the structural damage level. According to Dragon et al. (2000), an unloading process presents an inelastic behavior, being only a specific effect from damage by microcracking combined with friction of closed microcracks, the last being representative of a specific type of plasticity. In that study, it was developed a model involving the anisotropic damage by microcracking and friction in a series of closed microcracks. The formulation presented employs formalism with internal variable based on some concepts from micromechanics. Nevertheless, some recent studies sought to simulate the stiffness loss of concrete due to shear stress, such as: Xu et al. (2012), Kara and Dundai (2012) and Nogueira et al. (2011).
In the following section two possible proposals are presented to consider the damage effect due to shear in the simulation of the concrete behavior, being one of them chosen in order to meet some requirements to use the damage model originally proposed.

Proposals of laws of evolution for damage variables connected to shear
Concisely, it is reported two forms suggested for the incorporation of the damage effects due to shear stress in the concrete behavior.
PROPOSITION 1: For prevailing tension states, the variable D 4 representative of the damage generated by the shear between the crack edges belonging to the transverse isotropy local plane of the material is dependent on the variable D 1 , which represents the damage in a perpendicular direction to that plane. The same is true for the predominant states of compression, where the involved variables are D 5 and D 3, with the reservation that D 3 (damage in the plane of transverse isotropy) penalizes the stiffness perpendicular to the crack opened by the extension orthogonal to compression applied to the medium.
This proposal can be expressed by: where Ω T and Ω C can be parameters to be identified, or, in a more complex form, it can be proposed a function of scalar value composed of parameters. This second alternative can lead to a higher efficiency of the response obtained with the model, however leading to a higher cost in the parametric identification.
The expressions (13) and (14) consider the indirect association between the penalty of the transversal elasticity module and the strain perpendicular to the crack. This association was reported by Kara and Dudai (2012). For the predominant tension case, the variable D 1 is present in the calculation of strain perpendicular to the crack and, in the case of the predominant states of compression, the variable D 3 is present in the calculation of strain perpendicular to the crack.
PROPOSITION 2: In this case, are proposed independent forms of evolution of damage by shear through laws of evolution for D 4 and D 5 pursuant to Eq. (12), regardless of not being connected with the evolution of other damage variables as in the proposition 1. However, this introduces a considerable number of parameters to be identified and still causes a greater computational cost for the calculation of damage variables, since it would be one more variable to be considered in the iterative process for estimating the damage variables that form an implicit system to the relationships for Y T,C and D T,C . Other inconvenience would be the need for results of experimental tests in specimens of concrete, results illustrated by stress x distortion strain.
In general, both alternatives bring difficulties of execution, such as: How to propose laws of evolution if the experimental tests are difficult to perform? How to relate D 1 and D 4 through Ω T and D 3 and D 5 through Ω C ? The original model has a reduced number of parameters for a concrete anisotropic model, but does it worth to considerably increase this number of parameters to represent the shear?
Examining all these issues, a viable alternative would be the proposition of a simple law of evolution for variables linked to shear stress. The results obtained so far with this model are satisfactory (PITUBA, 2006(PITUBA, , 2008(PITUBA, , 2010, for instance) and the contribution of energy dissipation due to the damage variables D 4 and D 5 can be important, but not essential to the concrete behavior modeling regarding the constitutive model. This leads us to avoid overloading the model formulation due to shear in terms of computational cost, because the iterative procedure to calculate the variables of damage would have to count on more variables. Besides, there would be a higher number of parameters to be identified, because since the beginning of the model conception it was thought about a more practical application, without losing the quality of results. This is proved by the choice of transverse isotropy hypothesis as induced anisotropy and not the orthotropy, which would be more natural, but the number of parameters would be excessive to a small benefit.
Thus, it is considered the plot shear stress (τ) vs. distortion (γ) designed to the concrete behavior ( Figure 1). where G 0 is the transversal elasticity module for the intact material, γ 0 is the corresponding deformation for the peak stress in a test of pure shear, γ u is the last deformation admitted for the material in a pure shear test and α T is the straight inclination in the 'softening' branch. Therefore, in the case of predominant tensile states, the parameters to be identified are γ 0 , γ u and α T .
On the other hand, leveling to Eq. (16) with the relationship stress-deformation of a damaged material, it is obtained the variable of damage value: Therefore, in summarized form, the proposal is given by: being D 5 a constant value given a previous damage in compression.
The difficulty lies in the identification of the parameter γ 0 in a test of pure shear. The parameters γ u and α T can be conveniently chosen to approximate to the softening branch of the concrete. It is proposed that this law must be also valid for the predominant states of compression, i.e., α T = α C , where D 5 is who begins to evolve. However, it may be thought in a milder penalty for G in compression, assigning a higher value to γ u and a lower value to α C .
It is noteworthy that it was opted for simplicity in formulation and parametric identification, because the main dissipative phenomena related to the processes of damage are already considered by the model with the evolution of variables D 1 , D 2 and D 3 .

Results and discussion
A two-dimensional version of the proposed model was implemented in a finite element code for plane analyses. In these analyses, only concrete has physically nonlinear behavior, for steel it is admitted a linear constitutive relationship. Regarding the interaction between these two materials, it was admitted perfect adherence between steel and concrete.
Using the symmetries of loading and geometry, it was discretized only half of the beam, as illustrated in Figure 2.
The modified model was used in the analysis of reinforced concrete beams. The model responses with or without considering penalty for the relative stiffness to concrete shear are confronted in order to obtain a possible improvement in the numerical response given by the damage model.
The beams have two supports, span of 2.40 m, rectangular cross section (12 x 30 cm) and loading consisting of two equal concentrated forces applied on the interspaces thirds. For producing the beams it was used a concrete with E c = 29200 MPa. It was adopted E a = 196000 MPa for the reinforcement steel. The Poisson coefficient adopted was 0.20. More details on the experimental response from every type of beam, collected from tests performed with loading control can be found in Pituba and Fernandes (2011). The beams are described in the Figure 3, where are provided further details on the geometry and reinforcement. The Table 1 lists the values for the parameters related to D 1 , D 2 and D 3 , identified through uniaxial and biaxial stress tests, more details are described in Pituba and Fernandes (2011).   Table 1. Parameters for the model related to D 1 , D 2 and D 3 .

Parameters of traction
Parameters of compression Y 01 = 0.72 x 10 -4 MPa Y 02 = 0.5 x 10 -3 MPa Y 03 = 0.5 x 10 -3 MPa A 1 = 50 A 2 = -0.9 A 3 = -0.6 B 1 = 6700 MPa -1 B 2 = 0.4 MPa -1 B 3 = 70000 MPa -1 Regarding the parameters related to the evolution of the variables of shear damage D 4 and D 5 , aimed at simplicity it was assumed equal values to the behavior under traction and compression. As reference values it was adopted the experimental curve of a direct shear test (Figure 4) as presented by Khaloo and Ahmad (1988) for a concrete with f ck approximately 36 MPa, similar to the value achieved by the concrete used for beam production. The parameters related to the damage process under shear were: α = -90.0 MPa, γ 0 = 2.0 x 10 -6 and γ u = 5.0 x 10 -3 . The figures below present the responses obtained by the model without considering the evolution of damage variables D 4 and D 5 associated with the shear (curves designated by "Model" in the legend) and the corresponding with the addition of damage generated by shear concrete (curves designated by "Model with shear"). Also, it can be observed the instability in the responses for the last stages of loading where the damage location is present pronouncedly. Moreover, the steel is considered elastic, which is not compliant with the experimental response for the last stages of loading. The analyses were performed using meshes consisting of 400 quadrangular finite elements of 4 nodes, and in the case of beam with 7φ10.0 mm a mesh with 600 elements was required. It is noteworthy that only a single layer represents the reinforcement. In general, the contribution of damage related to concrete shear to energy dissipation is not very important, as expected, because the strain energy by bending translated into normal stresses acting in the opening of cracks is the main phenomenon in the concrete, according to the proposed fundament for the damage model. Probably the energy by distortion is not so influent due to the numerical response locking, still perceived, with the use of this type of finite element. However, the numerical responses continued to be satisfactory, where the damage process related to shear is one more phenomenon which contributes to the energy dissipated by the beam. But for the beam with 3#10.0 mm it was evidenced the model incapacity to capture the initial strong break of stiffness given around 12 kN, which interferes with the whole beam response, leading to a more stiff curve than the experimental. In conclusion, a nonlocal formulation of the model with the use of finite element that does not have locking limitations of the response is desirable for improving the proposed model.  Another problem is the difficulty of parametric identification of the model, now related to the damage variables associated with the shear. Despite the use of experimental values of a concrete with characteristics close to those used in the production of tested beams, it was found satisfactory numerical responses with the use of the proposed model.
Finally, by considering anisotropy induced by damage, the model presented a good result regarding the dissipation of energy generated by shear, because selectively penalizes the components of the constitutive tensor corresponding to normal and the tangential stress, which does not occur in an