Tension-Compression Damage Model with Consistent Crack Bandwidths for Concrete Materials

*is paper proposes a tension-compression damage model for concrete materials, formulated within the framework of thermodynamics of irreversible processes. *e aim of this work is to solve the following problems: the premature divergence of numerical solutions under general loading conditions due to the conflict of tensile and compressive damage bounding surfaces, which is a result of the application of the spectral decomposition method to distinguish tension and compression, and the unsatisfactory reproduction of distinct tension-compression behaviors of concrete by strain-driven damage models.*e former is solved by the sign of the volumetric deformation, while the latter is solved via two separated dissipationmechanisms. Moreover, of specific interest is an improved solution to the problem of mesh-size dependency using consistent crack bandwidths, which takes into account situations with irregular meshes and arbitrary crack directions in the context of the crack band approach. *e performance of themodel is validated by the well-documented experimental data.*e simplicity and the explicit integration of the constitutive equations render the model well suitable for large-scale computations.


Introduction
Concrete is still one of the most widely used construction materials in civil engineering applications [1].erefore, developing an efficient constitutive model for concrete is of great importance since it contributes to engineering designs, failure analyses, and life-cycle assessments by reducing costs and time [2].Microscopically, the mechanical behavior of concrete materials is mainly governed by the nucleation and growth of microcracks and their coalescence into macrocracks [3], which may result in the failure of a structure.Physically, it is an irreversible thermodynamic process accompanying energy dissipation and changes of the microstructure [4].Such a progressive degradation of material properties can be phenomenologically described by models [5][6][7][8][9][10][11][12][13][14][15][16][17] based on continuum damage mechanics (CDM), which is generally accepted by the scientific and engineering communities [18].
Despite the achievements of CDM over a wide range of materials (e.g., [19][20][21][22]), the rational and accurate CDM modeling of concrete materials still represents a challenging task.ree modeling challenges have been identified in finite element (FE) analyses.
e first challenge is the criterion distinguishing tension and compression, which is used for activating/inactivating damage evolutions in tension and compression.e existing models usually utilize the spectral decomposition of a stress or strain tensor into tensile (positive) and compressive (negative) eigenvalues (e.g., [23,24]), in order to describe different behaviors of concrete under general loading conditions.However, in the use of this method, the following problem would arise: divergence of numerical solutions will occur under biaxial or reverse loading conditions due to the conflict of tensile and compressive damage bounding surfaces induced by Poisson's ratio effect [25].A better solution to this problem is the introduction of a weighted damage parameter (e.g., [5,25]).Moreover, the sign of the volumetric strain ( [26]) is used to distinguish what is tension and what is compression.However, it should be noted that the aforementioned method is not applicable to incompressible materials which have a value of Poisson's ratio around 0.5 (e.g., fresh concrete and concrete under confinement) since the volumetric strain is null.
erefore, the present model is limited to the compressible materials characterized by Poisson's ratio lower than 0.5.
e second challenge is the modeling of tensioncompression behaviors of concrete materials which exhibit large differences in their peak strengths and postpeak stages, as well as the unilateral effect characterized by a total or partial recovery of the degraded stiffness during the tension-compression load reversal.
e classical straindriven damage models incorporating one damage criterion alone may be incapable of predicting such distinct behaviors.Among them (e.g., [27][28][29][30][31][32]), one case can be mentioned as an example: the predicted strength of concrete under biaxial compression is seriously lower than that under uniaxial compression.To solve the problem, this paper proposes two separate damage criteria describing the state of damage under tension and compression, respectively, whereby the distinct behaviors can be reproduced in a physically realistic fashion.Due to the scope of this paper, the unilateral effect will not be discussed.
e third challenge is to avoid the lack of mesh objectivity of numerical results when the modeling of strain-softening is taken into account.At the stage of strain-softening, the concrete does not lose its entire loadcarrying capacity but exhibits a progressively decreasing residual strength.erefore, capturing the strain-softening is important because it is responsible for predicting the behavior of overloaded structures.However, damage models based on conventional continuum theories usually produce mesh-dependent results in FE analyses [6].e mesh dependency can be further classified into mesh-size [33] and mesh-bias [34] dependencies.
e mesh-size dependency, which originates from the loss of ellipticity of the governing partial differential equations for static problems [35], which means that both the computed energy dissipation due to strain-softening and the corresponding force-displacement response are governed by the fineness of the FE mesh [36].It also reflects a physically unrealistic fact that the infinitesimal element size leads to a vanishing energy dissipation during the structural failure [37].e mesh-bias dependency, which is a result of the inadequacy of the spatial discretization procedure [38,39], means that the numerical result depends on the orientation of the FE mesh in the sense that the computed direction of the strain localization band tends to propagate along continuous mesh lines [40], which may give rise to spurious crack trajectories.Although the problem of meshbias dependency is important for damage modeling of concrete materials, this study focuses only on the solution to the mesh-size dependency.
To solve the problem of mesh-size dependency in the framework of continuum mechanics, various types of regularization methods have been developed that can be distinguished by the crack band approach (or fracture energy regularization) [33,41,42], nonlocal models including the integral-type [6,43,44] and the gradient-type [13,45], micropolar method [46], and viscous or rate-dependent [11] methods.On one hand, the nonlocal method seems to be more recent (e.g., [28,30,47]).However, the nonlocal models have their own difficulties such as (i) incorrect crack initiation, ahead of the crack tip [48], and (ii) incorrect shielding effect with nonzero nonlocal interactions across a crack surface [49].Furthermore, nonlocal methods have no potential for large-scale computations of concrete structures since they require small enough elements, e.g., their sizes are at least one-third of the characteristic length of the material [50], to subdivide strain localization bands of the structure, which are usually a priori unknown.On the other hand, the crack band approach is widely used in most of the commercial FE packages and is still an important analytical tool for practical applications (e.g., [51,52]).One goal of the present work is to propose a constitutive model which needs to ensure a high algorithmic efficiency and an easy implementation; therefore, a fracture energy regularization based on the crack band approach [33] is adopted.
One key parameter in the crack band approach is the crack bandwidth, through which the fracture energy (regarded to be a material property) is related to the volumetric (or specific) fracture energy, and it therefore ensures an objectivity of the energy dissipation regardless of the mesh refinement.
e existing estimations of the crack bandwidth can be divided into three categories: (i) methods based on the square root of the element area (for twodimensional elements) or the cubic root of the element volume (for three-dimensional elements) (e.g., [53][54][55]); (ii) the notable methods proposed by Oliver [56]; and (iii) projection methods (e.g., [57][58][59][60]).e methods in the first category, which estimate the crack bandwidth as a constant during the numerical analysis process, are less applicable to situations with irregular meshes or arbitrary crack directions [56].e methods in the last two categories, in which crack bandwidths are varied with respect to crack directions and spatial positions, solve the aforementioned problems.However, they may encounter difficulties such as excessively large estimation [59] or misprediction of the crack bandwidth.is paper therefore develops an improved estimation of the crack bandwidth to circumvent the above difficulties.Moreover, two different crack bandwidths (one for the tensile strain-softening and the other for the compressive strain-softening) are taken into account in the present study, because of which concrete exhibits quite different localization behaviors in tension and compression [61].It should also be noted that the crack band approach can only be regarded as a partial regularization method for mesh-size dependency in the sense that the objectivity of the energy dissipation is ensured [62], and thus the global structural response can be captured correctly [56], but the mathematical description of the rate boundary value problem (BVP) still remains ill-posed [63].
e main purpose of this study is not to describe all the mechanical behaviors of concrete but rather to give methodological efforts to solve the aforementioned difficulties.erefore, the plastic dissipation is not considered and the isotropic damage is assumed for simplification.Accordingly, 2 Advances in Civil Engineering the objectives of this paper are fourfold: (i) to define the criterion distinguishing tension and compression by the sign of the volumetric deformation; (ii) to propose a model to describe distinct tension-compression behaviors of concrete materials under a framework similar to that of Cervera et al. [11] and Comi and Perego [26]; (iii) to develop an improved estimation of the crack bandwidth using the projection method for solving the problem of mesh-size dependency; and (iv) to validate the proposed model by the selected benchmark tests.

Tension-Compression Damage Model with Consistent Crack Bandwidths
2.1.Helmholtz Free Energy Potential.e thermodynamics of irreversible processes provides a general framework for formulating constitutive equations.In other words, the state of the material can be described by the thermodynamic potential.Consequently, the stresses and other associated variables can be derived from this potential [64].Because of the scope of this paper, the present model is developed based on the following assumptions: (i) damage evolutions are assumed to be isotropic and to occur in small strains and (ii) the rate-dependent behavior and the thermal effect are not considered.To take into account distinct tension-compression behaviors of concrete, the postulated thermodynamic potential in terms of Helmholtz free energy ψ is decomposed into volumetric and deviatoric parts ψ v and ψ d in a similar way as in the study of Cervera et al. [11] and Comi and Perego [26] as follows: where ρ is the mass density of the material, ε is the secondorder strain tensor, d t and d c are the damage internal variables that represent the level of material degradation in tension and compression, K 0 and G 0 are the initial bulk and shear moduli of the undamaged material, respectively, H(•) is the Heaviside function defined by H(x) � 1 if x ≥ 0 and H(x) � 0 if x < 0, and tr(•) is the trace of a tensor.e tensile and compressive damage internal variables are monotonously increasing scalars, both of which grow from zero (undamaged material with the intact stiffness) to one (completely damaged material).ε v and ε d are second-order volumetric and deviatoric strain tensors, defined by the following equation: where I v and I d are the fourth-order volumetric and deviatoric projections, respectively, and I and I s are the second-order identity tensor and the fourth-order symmetric identity tensor, respectively.ε + v and ε − v are positive and negative volumetric strain tensors, defined according to

Damage Criteria.
To describe the damage due to tension or compression, it first necessitates a criterion that defines tensile/compressive states.As mentioned in Introduction, the method of spectral decomposition may lead to the premature convergence problem in numerical solutions [25].erefore, the sign of the volumetric deformation (relative volume change), which is calculated by trε, is used to distinguish tension and compression in this paper.e similar application can be found in the work of Comi and Perego [26].In the following, the tensile and compressive strain states are thus defined by trε ≥ 0 and trε < 0, respectively.
e damage criterion characterizes the state of damage and represents a bounding surface whose boundary corresponds to the states at which the damage will grow.To clearly specify the damage due to tension or compression, two damage criteria are introduced as follows: where F t and F c are the tensile and compressive damage criteria, ε eq t and ε eq c are the equivalent tensile and compressive strains, and κ t and κ c are the current tensile and compressive damage thresholds, respectively.It is noteworthy that the damage does not develop if F i < 0 and grows if F i � 0.
e damage threshold represents the largest value of the equivalent strain ever reached by the material during the loading history [65][66][67].erefore, κ t and κ c can be further expressed as follows: Similar to the work of Fichant et al. [68], the equivalent strains in equation ( 10) are assumed to be as follows: with Advances in Civil Engineering where E 0 is the initial Young's modulus of the undamaged material, I 1σ is the first invariant of the effective stress tensor, J 2σ is the second invariant of the deviatoric effective stress tensor, σ e is the effective stress tensor, C 0 is the fourth-order elastic stiffness tensor of the undamaged material, α i and β i are the material constants, f t is the uniaxial tensile strength, and f c and f bc are the uniaxial and biaxial compressive strengths of the material, respectively.Following the experimental result from the study of Kupfer et al. [69], the typical value of f bc /f c � 1.16 for concrete is adopted throughout this paper.

Evolution Laws for Internal Variables.
From the macroscopic viewpoint within the framework of thermodynamics, the damage variable can be introduced as an internal state variable to quantify the average material degradation [22].
In order to capture distinct tension-compression behaviors of concrete, two types of damage evolution laws, which relate the corresponding damage internal variable to the associated damage threshold, are defined in the following general form: In this paper, explicit damage evolution laws are further postulated as follows: where κ 0 t � f t /E 0 and κ 0 c � cf c /E 0 are the initial tensile and compressive damage thresholds, respectively, with c being the material parameter.A t is the nonnegative material parameter that controls the shape of the tensile stress-strain relation, as shown in Figure 1.It is noteworthy that if A t � 1, then the damage evolution law in tension represents a linear softening.κ u t and κ u c are the ultimate tensile and compressive damage thresholds, respectively.According to the crack band approach [33], the ultimate damage threshold depends on the fracture energy and the crack bandwidth.It will be further discussed later.
At any material point, damage can evolve only if the current state reaches the boundary of the elastic domain, characterized by F i (ε eq i , κ i ) � 0. To reflect the loading, unloading, and reloading conditions more generally, the evolution of damage internal variables is governed by the following Kuhn-Tucker relations:

Mechanical Dissipation and Constitutive Laws.
e failure of concrete is related to the degradation of physical properties of the material.During this damage process, the energy dissipation is always nonnegative and accompanied by an increase of entropy [66].It hence implies an irreversible process of thermodynamics.erefore, the mechanical dissipation ϕ d under the isothermal condition needs to satisfy the second principle of thermodynamics, which can be expressed by the Clausius-Duhem inequality as follows: where σ is the second-order stress tensor.Substitution of equation (1) into equation ( 19) yields the following equation: with the tensile and compressive damage energy release rates defined as follows: where C is the fourth-order secant stiffness tensor of the damaged material.
Because of the fact that the inequality of equation ( 20) must be satisfied by the arbitrary change of the strain tensor ε, the coefficient of the first term of the inequality, therefore, must be zero.is requirement yields the following constitutive laws: 4

Advances in Civil Engineering
Using equation ( 22), the Clausius-Duhem inequality expressed by equation ( 20) is eventually written in the following form: erefore, the second principle of thermodynamics is satis ed as long as the damage variables d t and d c are monotonically increasing since Y t and Y c are strictly nonnegative because of their quadratic forms, which involve zC/zd i ≤ 0.

Fracture Energy Regularization with Consistent Crack
Bandwidths.To solve the problem of mesh-size dependency, a fracture energy regularization based on the crack band approach [33] is adopted in the present study.e general idea of the crack band approach is to rescale the softening branch of the stress-strain curve, in order to impose the same energy dissipation per unit area for di erent element sizes [59]. is adjustment consists in relating the fracture energy (regarded to be a material property) to the volumetric (or speci c) fracture energy (Figures 2(c) and 2(d)), which depends on the crack bandwidth.To achieve this goal of regularization, it is nally required to relate the ultimate damage threshold to the fracture energy and the crack bandwidth.
e crack band approach assumes that the failure of material occurs within a damage zone having a limited width, whereas the other parts are unloaded [59] during the process of failure (Figure 2).According to Bažant and Oh [33] and Jirásek and Mar a [43], the same energy dissipation per unit area can be ensured using following relations: where G f t and G f c are the Mode-I and compressive fracture energies de ned as the energy dissipation per unit area of a surface due to cracking or crushing, l (P)  t and l (P) c are the tensile and compressive crack bandwidths, g f t and g f c are the tensile and compressive volumetric fracture energies de ned as the energy dissipation per unit volume, and σ nn t and σ nn c are the tensile and compressive stresses normal and tangential to the surface of cracking, respectively.It is noteworthy that G f t and g f t can be represented by the area under the tractionseparation (tensile stress normal to the crack versus crack opening displacement) curve and under the stress-strain curve, respectively, as shown in Figures 2(c) and 2(d).
Using equation ( 16), substitution of equation ( 22) into equation ( 25) yields the following equation: It is noteworthy that equation ( 27) is derived under the condition of uniaxial loading.Nevertheless, the previous studies (e.g., [20,70,71]) have revealed that it is also an acceptable approximation under general loading conditions.
Similarly, the compressive damage threshold κ u c can be expressed as follows: Generally, the crack bandwidth can be regarded as an FE discretization parameter, which depends on the chosen element shape (e.g., quadrilaterals or triangles), element size, element orientation (or mesh alignment) with respect to the direction of the strain localization band, interpolation function or order of the displacement approximation (e.g., linear or quadratic), and numerical integration scheme (number of integration points) [53,59,60].For the sake of simplicity, only four-node quadrilateral and three-node triangular elements are used throughout the paper.It should also be noted that only low-order elements (i.e., elements with linear or multi-linear interpolation functions) are adopted in this study since the stress elds obtained from high-order elements (i.e., elements with quadratic interpolation functions) are highly disturbed for crack band simulations [59].
Inspired by the works of Červenka et al. [57] and Jirásek and Bauer [59], this paper develops an improved estimation of crack bandwidths in tension and compression by projecting the element onto the directions normal and tangential to the crack direction, respectively.Mathematically, they can be written as follows: with In uence of A t on uniaxial stress-strain relation in tension.
Advances in Civil Engineering 5 where l (P) t and l (P) c are the tensile and compressive crack bandwidths, l max i and l min i (i t or c) are the maximum and minimum projection values, x N is the midpoint of the element edge (for four-node quadrilateral elements) or the nodal point (for three-node triangular elements), N m−n is the number of element midpoints (for four-node quadrilateral elements) or nodes (for three-node triangular elements), x represents the point (e.g., Gaussian point) in the structure where the crack bandwidth is to be evaluated, and n(x) and t(x) are the unit vectors with the direction normal and tangential to the crack direction at the point x, respectively.e geometrical meaning of equation ( 29) is illustrated in Figure 3.
For comparison, the existing estimation formulas of the crack bandwidth mentioned in Introduction are recalled as follows.
e tensile crack bandwidth l (R) t , which is estimated by Rots [53] belonging to the rst category, can be written as follows: , for linear two-dimensional elements,

   (32)
with A being the area of the element.e tensile crack bandwidth l (O) t , which is estimated according to the rst method of Oliver [56], can be expressed as follows: with where N c is the number of corner nodes of the element, N i is the standard shape function for the corner node i belonging to an element with N c corner nodes, ϕ i is the nodal value at the corner node i, and x i and x c are the corner node and the center of the element, respectively.e tensile crack bandwidth l (G) t , which is estimated by Govindjee et al. [58] belonging to the third category, can be given as follows: with e geometrical meanings of the aforementioned tensile crack bandwidths and their di erences with respect to the direction of n(x) are illustrated in Figures 4 and 5, where θ is the angle between the normal to the crack and the global axis x at the point x.29) and (35), respectively) are free from such problems.e discontinuous crack bandwidth obtained from equation ( 33) is caused by the discontinuous nature of ϕ i , whose value is either 1 or 0 as de ned by equation (34).It is also noteworthy that the present method (labelled as l (P) t ) and the method of Govindjee et al. [58] (labelled as l (G)  t ) give exactly the same e further comparative study of these estimation methods will be given in the next section.

Model Verification
To establish the validity of the proposed constitutive model, the computed results were either veri ed theoretically or compared with measured ones from experiments.ese are presented in the following sections.

Analysis of Biaxial Tests of Concrete.
To examine the rationality of the adopted criterion distinguishing tension and compression, and to validate the present model in reproducing the distinct strength features of concrete under biaxial stresses, a series of plain concrete specimens subjected to di erent stress ratios (denoted by α), which were tested by Kupfer et al. [69], were chosen for the numerical analysis.e material properties used for the computation were Young's modulus E 0 3.2 × 10 4 MPa, Poisson's ratio ] 0 0.2, concrete uniaxial tensile strength f t 3. Advances in Civil Engineering mesh (Mesh A was 200.0 mm long × 200.0 mm high) with a four-node plane stress quadrilateral element using a 2 × 2 Gaussian integration scheme [44] was used to represent the plain concrete specimen (Figure 6).
Figure 6 shows the comparison between the experimental and computed results of the strength envelope, which corresponds to the envelope of the maximum stress ever reached in the material during the loading, in a biaxial principal stress space.e results indicate that the predicted biaxial strengths of concrete are very close to the experimental ones, either in the biaxial tensile domain or in the biaxial compressive domain, or even in the tensilecompressive domain.
e numerical results in Figure 6 also demonstrate the rationality of the adopted criterion distinguishing tension and compression.Figure 7 shows the comparison of stress-strain curves of numerical and experimental results [69] under biaxial compression (σ 3 0), according to the following stress ratios (α σ 2 /σ 1 ): α 0, α 0.52, and α 1.It shows that the proposed model can predict well both the strength enhancement of concrete and the corresponding peak strain in biaxial compressive states.

Analysis of a Double-Edge Notched Concrete Specimen under Direct Tension.
To validate the proposed model in capturing the strain-softening, a double-edge notched (DEN) concrete specimen under direct tension, which was investigated experimentally by Hordijk [72], was used for the numerical analysis.e test specimen shown in Figure 8, which had dimensions 125 mm height × 60 mm width × 50 mm thickness with two symmetrical notches of 5 × 5 mm 2 in the middle cross-section, was glued to the loading platens and subjected to an axial tensile load under the deformation control.
e average deformation δ, which denotes the elongation of the specimen, was measured experimentally by using four pairs of extensometers located at the front and rear faces of the specimen.e gauge length of extensometers was 35 mm (Figure 8).e average stress σ avg in the experiment was de ned as the vertical load divided by the net cross-sectional area (50 × 50 mm 2 ).
e material properties of the test specimen used for the numerical analysis were Young's modulus E 0 1.8 × 10 4 MPa, Poisson's ratio ] 0 0.2, concrete uniaxial tensile strength f t 3.5 MPa, A t 6.0, uniaxial

Node of element Gaussian point Crack
Node of element Gaussian point Crack   Advances in Civil Engineering compressive strength f c 34.0 MPa, and Mode-I fracture energy G f t 0.07 N/mm.e DEN concrete specimen was modeled by a 298element mesh with four-node plane stress quadrilateral elements using a 2 × 2 Gaussian integration scheme, as shown in Figure 9.According to the boundary conditions of the experiment, all FE nodes at the lower boundary were assumed to be xed (Figure 9). Figure 10 presents the computed average stress-deformation curve in tension.e numerical result shows a satisfactory agreement with the experimental one, including both the load-carrying capacity and the postpeak behavior.

Analysis of a Bar Subjected to Uniaxial Tension Using
Elements with a Constant Strain Field.To prove the e ciency of the proposed estimation of the crack bandwidth with respect to the FE mesh irregularity, a bar subjected to uniaxial tension was used for the numerical analysis by following the work of Oliver [56].Figure 11 shows the geometry, boundary, and loading conditions of the bar with a thickness of 50 mm.e material properties used for the numerical analysis were Young's modulus E 0 2.0 × 10 4 MPa, Poisson's ratio ] 0 0.0, concrete uniaxial tensile strength f t 2.0 MPa, A t 1.0, uniaxial compressive strength f c 20.0 MPa, and Mode-I fracture energy G f t 0.1 N/mm.According to Oliver [56], the localization of the bar due to strain-softening was induced arti cially by the prescribed imperfection (Figure 11) with a small reduction of the uniaxial tensile strength in shaded elements, as shown in Figure 12.
e bar was discretized by the four-node quadrilateral or three-node triangular plane stress elements, both of which  Advances in Civil Engineering adopted a linear interpolation function (implying a constant strain eld) and a 1-point Gaussian integration scheme.
Figure 12 shows the meshes used for the analysis with di erent mesh irregularities.Due to the fact that the projection methods (e.g., the present method and that of [58]) have solved the di culty of the excessively large estimation of the crack bandwidth arising in the rst method of Oliver [56] shown in the previous section, only the element area-based, Govindjee et al. [58]-based, and the present methods of estimating the crack bandwidth are compared in this section.
Figure 13(c) shows that the load-displacement responses obtained by the present method are almost coincident with the theoretical solution and are hardly detectable in the plots, while the corresponding responses, which were obtained by the crack bandwidth estimated by 2A √ and by Govindjee et al. [58] using equation (35), show an unacceptable mesh dependency with respect to the FE mesh irregularity.Consequently, the present model provides an accurate estimation of the crack bandwidth that yields meshsize independent results with respect to the mesh irregularity.integration) and three-node triangular plane stress elements in the FE analysis, respectively.ree di erent mesh sizes were used for the mesh-size dependency study: coarse, medium, and ne meshes, which had a total of 480, 810, and 1176 elements, respectively (Figure 15).e corresponding minimum element sizes of these three meshes were 19.1 mm, 12.8 mm, and 9.6 mm, respectively.Figure 16 presents the predicted structural responses in terms of the load P at point B versus the crack mouth sliding displacement (CMSD).Figure 17 shows the simulated nal tensile damage distribution, which is denoted by the dark color corresponding to a high damage value close to one.Figure 18 shows the simulated crack patterns at the nal computation stage.
According to the numerical results, it is clear that (i) the global structural responses are well predicted compared to the experimental results (Figure 16); (ii) the proposed model provides an acceptable prediction of the curved macroscopic crack, as compared with the experimentally observed crack trajectory (Figures 17 and 18); and (iii) the simulated crack trajectories obtained with the medium and ne meshes are identical, which indicates that the predicted paths of the crack propagation are almost mesh-size independent.
e predicted thicknesses of the macroscopic crack are inconsistent when using di erent mesh sizes, as shown in Figures 17 and 18.As more clearly shown in Figure 18, the simulated localization band tends to occur in one row of elements.is is due to the inherent de ciency of the crack band approach since the rate BVP still remains ill-posed [63].e similar conclusion to simulations with low-order elements can be found in the study of Jirásek and Bauer [59], who also pointed out that the minimum size of the localization band is the size of one row of integration points (i.e., subelement domains) but not the size of one row of elements in simulations with high-order elements.
As also can be observed in Figures 17 and 18, the predicted crack trajectories are less curved than its experimental counterpart due to the approximation error of the spatial discretization in the standard nite element method.It can be explained by the fact that the localized damage zone cannot be resolved with su cient accuracy because of the minimum element size (19.1 mm) of the coarse mesh being too large, as compared with the thickness of the notch (5 mm).

Conclusions
is paper presents a tension-compression damage model with consistent crack bandwidths developed for the macroscopic description of concrete behaviors.
e model is formulated within the framework of the internal variable theory of thermodynamics so that the constitutive relations are derived without ad hoc assumptions.Based on the numerical and experimental results, the following conclusions can be drawn: rational definition of two equivalent strains, two internal variables, and two separated damage criteria describing distinct damage mechanics under tension and compression.Such features enable the present model to provide a better prediction of the distinct tension-compression behaviors of concrete materials.Furthermore, the adopted method using the sign of the volumetric deformation avoids the premature convergence difficult in the numerical analysis, which is often encountered in existing models using the spectral decomposition of a stress or strain tensor.(ii) e main contribution of this paper is that it provides an improved estimation of the crack bandwidth, which is developed for two-dimensional linear elements by taking into account the element shape, element size, and element orientation with respect to crack directions and spatial positions.Its capability in numerical analyses of irregular meshes and arbitrary crack directions has been validated by the comparative study of existing estimations of the crack bandwidth at the element level.e proposed model incorporating this enhanced description of the crack bandwidth has also been proven to be mesh-size independent at the structural level.e present model is able to provide objective simulation results with respect to the finite element mesh refinements that are very important for the computation of the real engineering structures.(iii) e proposed model has been validated by the welldocumented experimental data of concrete.e good agreement between the predicted and experimental results consequently demonstrates that the present model provides an effective tool for Advances in Civil Engineering modeling concrete behaviors in tension and compression.(iv) e present model has a limited number of material parameters, each of which has a clear physical meaning and can be identified by standard laboratory tests.Moreover, the simplicity and the explicit integration of the constitutive equations make the model well suited for large-scale computations.
e crack bandwidths denoted by l (O) t (according to equation (33)) become very large (Figures 4(b) and 5(b)) under certain circumstances, and they also show a discontinuity (Figures 4(c) and 5(c)), while the crack bandwidths denoted by l (P)

Figure 2 :
Figure 2: Damage zone and fracture energy regularization: (a) damage zone in the concrete specimen under tensile loading; (b) stress-strain response outside the damage zone; (c) traction-separation law in the damage zone; (d) general stress-strain relation for the concrete specimen.

Figure 5 :Figure 6 :
Figure 5: Comparison of estimations of crack bandwidth for a right triangular element: (a) geometrical meanings of crack bandwidths estimated by the present model, Oliver [56], and Govindjee et al. [58]; (b) large crack bandwidth estimated by the rst method of Oliver [56]; (c) relations of tensile crack bandwidth and angle θ.

Figure 9 :
Figure 9: FE mesh for the DEN concrete specimen.

Figure 15 :Figure 16 :
Figure 15: Di erent FE mesh sizes for the SEN concrete beam in the mesh-size dependency study: (a) coarse mesh; (b) medium mesh; (c) ne mesh.

Figure 17 :Figure 18 :
Figure 17: Final tensile damage profiles in the central part of the SEN concrete beam: (a) coarse mesh; (b) medium mesh; (c) fine mesh.