Abstract

A coupled plastic damage model with two damage scalars is proposed to describe the nonlinear features of concrete. The constitutive formulations are developed by assuming that damage can be represented effectively in the material compliance tensor. Damage evolution law and plastic damage coupling are described using the framework of irreversible thermodynamics. The plasticity part is developed without using the effective stress concept. A plastic yield function based on the true stress is adopted with two hardening functions, one for tensile loading history and the other for compressive loading history. To couple the damage to the plasticity, the damage parameters are introduced into the plastic yield function by considering a reduction of the plastic hardening rate. The specific reduction factor is then deduced from the compliance tensor of the damaged material. Finally, the proposed model is applied to plain concrete. Comparison between the experimental data and the numerical simulations shows that the proposed model is able to describe the main features of the mechanical performances observed in concrete material under uniaxial, biaxial, and cyclic loadings.

1. Introduction

The mechanical behavior of concrete is unique, due to the influence of micromechanisms involved in the nucleation and growth of microcracks and plastic flow. This behavior is characterized by several features as follows: the tensile strength, which is different from the compressive strength; the irreversible plastic deformation, which is found as the pressure on the concrete exceeds a certain value; beyond the peak stress, the stress-strain curve displays an unstable region, accompanied by stiffness degradation and strength softening. These features have to be considered in constitutive models of concrete materials. To model these features, several mechanics theories have been used. In general, the damage mechanics theories can be used to model the nucleation and growth of microcracks [1–3], whereas the plasticity theories can be used to model the plastic flow component of the deformation [4].

Because damage mechanics can be used to describe the progressive weakening of solids due to the development of microcracks, and because plasticity theories have been successfully applied to the modeling of slip-type processes, several plastic damage models of concrete have been developed through combining damage mechanics and plasticity theories [5–13]. These coupled plastic damage models (CPDMs) could be formulated in the irreversible thermodynamics framework and can be easily applied to describe the essential nonlinear performances of concrete including the strain softening and the stiffness degradation. One type of CPDMs is based on the concept of effective stress, which was initially proposed by Kachanov [14] for metal creep failure. In these models, the plastic yield function is defined in the effective configuration pertaining to the stresses in the undamaged material [5, 7–9]. Many authors used this approach to couple damage to plasticity. Some of these models, developed by the use of two damage scalars and damage energy release rate-based damage criteria, show excellent performance in reproducing the typical nonlinear behavior of concrete materials under different monotonic and cyclic load conditions. The other type of CPDMs is opposite of the above. In these models, the true stress appears in the plastic process, which clearly couples plasticity to damage [6, 10, 11, 15–17]. However, considering the influence of damage on plastic flow, it is difficult to develop a plastic yield function that can be used to describe the plastic deformation and damage growth of concrete simultaneously. Taqieddin et al. [17] suggested that a damage-sensitive plastic yield function should be introduced in terms of true stress configuration. Even if several types of expressions for the plastic yield function written in terms of the effective stress have been successfully applied to model some of the typical nonlinearities of concrete (such as the volumetric dilation and strength increase under multidimensional compression), they cannot be directly used in the true stress space. A common approach to solve this problem is to introduce a reduction factor of plastic hardening rate into the plastic yield function. The plastic yield function is usually expressed by a function of the stress tensor and plastic hardening function, so the damage parameters are included in the plastic yield function with the introduction of the reduction factor in the plastic hardening function. Several authors applied this approach to develop the CPDM with a true stress space plastic yield function. Shao et al. [15] proposed a single scalar damage framework based on the elastoplastic damage release rate. The influence of damage on the plastic flow is calculated by considering a reduction of the plastic hardening rate. The reduction factor is introduced into the plastic yield function for the coupling between the damage evolution and the plastic flow. Nguyen and Houlsby [10] proposed a double scalar damage-plasticity model for concrete based on thermodynamic principles. In this model, the plasticity part is based on the true stress using a yield function with two hardening functions, one for the tensile loading history and the other for the compressive loading history. Two reduction factors are introduced into the tensile and the compressive hardening functions to consider the influence of the tensile and compressive damage mechanisms on the plastic flow, respectively. Taqieddin et al. [17] proposed a double scalar damage-plasticity model for concrete based on the hypothesis of strain energy equivalence. The model incorporates a plastic yield function written in terms of the true stress. Damage variables are introduced all over the plastic yield function.

For the models mentioned above, the reduction factor of the plastic hardening rate is equal to the damage variable defined in these models, and there is no distinction in the aspect of choosing a reduction factor between uniaxial and biaxial loadings. The plastic hardening rate is generally defined by the equivalent plastic strain, so the sum of the principal damage values (in three directions: , , and ) seems reasonable to use as the reduction factor of the plastic hardening rate under multiaxial stresses. By introducing such a reduction factor, the widely used plastic yield functions, such as those applied by Lee and Fenves [7], Wu et al. [9], and Lubliner et al. [18], can be used in the true stress space.

In the present work, a coupled plastic damage model is formulated in the framework of thermodynamics. The constitutive formulations are developed by considering an increment in the concrete compliance due to microcrack propagation. The plasticity part is based on the true stress using a yield function with tensile and compressive hardening functions. The plasticity yield function widely used in effective stress space is modified to be applied in this study by considering a reduction in the plastic hardening rate. Finally, the proposed model is applied to plain concrete. Comparison between the experimental data and the numerical simulations shows that the proposed model is able to represent the main features of mechanical performances observed in concrete material under uniaxial, biaxial, and cyclic loadings.

2. General Framework of Coupled Plastic Damage Model

The model presented in this work is thermodynamically consistent and is developed using internal variables to represent the material damage state. The assumption of small strains will be adopted in this work. In the isothermal conditions, the state variables are composed of the total strain tensor , scalar damage variables and , plastic strain , and internal variables for plastic hardening and . The total strain tensor is decomposed into an elastic part and plastic part :

To establish the constitutive law, a thermodynamic potential for the damaged elastoplastic material should be introduced as a function of the internal state variables. Considering that the damage process is coupled with plastic deformation and plastic hardening, the thermodynamic potential can be expressed by where and are the elastic and plastic energy for plastic hardening of damaged material, respectively. is the fourth-order stiffness tensor of the damaged material.

According to the second principle of thermodynamics, any arbitrary irreversible process satisfies the Clausius-Duhem inequality as

Substituting the time derivative of the thermodynamic potential into the inequality yields

Because the inequality must hold for any value of , , , , , and , the constitutive equality can be obtained as follows:

The constitutive relation between the stress and strain tensors can also be expressed utilizing the material’s compliance tensor as where is the fourth-order compliance tensor of the damaged material, and is the matrix inverse function. It is assumed in (6) that the damage can be represented effectively in the material compliance tensor. This assumption is in line with many published papers. Budiansky and O’connell [19] and Horii and Nemat-Nasser [20] indicated that damage induces degradation of the elastic properties by affecting the compliance tensor of the material. Based on this consideration, several researchers proposed an elastic degradation model in which the material stiffness or compliance was adopted as the damage variable. The elastic degradation theory was then extended to the multisurface elastic-damage model [21] and to the combined plastic damage model [6, 15, 16, 22].

Considering that an added compliance tensor is induced by the microcrack propagation, the fourth-order compliance tensor is decomposed as follows [23]: where defines the compliance tensor of the undamaged material and defines the total added compliance tensor due to microcracks. Microcracks in concrete are induced in two modes: the splitting mode and compressive mode [6, 23]. In general, the splitting mode is dominant in tension, whereas the compressive mode is dominant in compression. To incorporate these modes into the formulation, the stress tensor is decomposed into a positive part and negative part : where the eigenvalues of are nonnegative, and the eigenvalues of are all nonpositive. According to Faria et al. [8], Wu et al. [9], and Pelà et al. [24], each part of the stress tensor can be expressed as where the fourth-order projection tensors and are expressed as in which is the fourth-order identity tensor, and are the th eigenvalue and the normalized eigenvector of the stress tensor , respectively, and is the Heaviside function ( for and for ).

By introducing the positive/negative projection operators , the total added compliance tensor is decomposed as [25] where and represent the added compliance tensors due to microcracks induced in tension and compression, respectively.

To progress further, evolutionary relations are required for the added compliance tensors and . Based on the previous work of Yazdani and Karnawat [6], Ortiz [23], and Wu and Xu [25], the added compliance tensors are expressed in terms of response tensors and such that where the response tensors and determine the evolution directions of the added compliance tensors and , respectively.

The thermodynamic forces conjugated to the corresponding damage variables are given by where the terms and are the derivatives of the stiffness tensor of the damaged material with respect to the damage variables and , respectively. These terms can be calculated utilizing the material’s compliance tensor as

Substituting (14) into (13) and calling for (5), one obtains

Finally, the rate form of the constitutive equation can be given as

The derivation of the detailed rate equation from (16) requires determining the evolution laws of the damage variables and plastic strains.

2.1. Damage Characterization

The damage criteria are defined as [15] where (tension/compression) represents the threshold of the damage energy release rate at a given value of damage. On the basis of the normality structure in the continuum mechanics, the evolution law for the damage variables can be derived by in which and denote the tensile and compressive damage multipliers, respectively. Under the loading-unloading condition,    and    are expressed according to the Kuhn-Tucker relations:

In the special case of elastic damage loading without plastic flow (), the damage multipliers and can be determined by calling for the damage consistency condition under static loading:

2.2. Plastic Characterization

Irreversible plastic deformation will be formed during the deformation process of concrete. The plastic strain rate can be determined using a plastic yield function with multiple isotropic hardening criteria and a plasticity flow rule. The yield function determines under what conditions the concrete begins to yield and how the yielding of the material evolves as the irreversible deformation accumulates [26]. According to Shao et al. [15] and Salari et al. [27], the plastic yield criterion can be expressed by a function of the stress tensor and the thermodynamic forces and , respectively, conjugated with the internal hardening variables and :

The plastic hardening functions and can be deduced by the derivative of the thermodynamic potential:

The plastic potential is also a function of the stress tensor, the scalar damage variables, and the internal hardening variables. Note that is the yield function in the associated flow rule. At the given function of the plastic potential, the evolution law of the plastic strain is expressed as follows:

In the special case of plastic loading without damage evolution ( and ), the plastic multiplier can be determined from the plastic consistency condition:

2.3. Computational Procedure

The numerical algorithm of the proposed constitutive model is implemented in a finite element code. The elastic predictor-damage and plastic corrector integrating algorithm proposed by Crisfield [28] and then used by Nguyen and Houlsby [10] is used here for calculating the coupled plastic damage model. In this algorithm, the damage and plastic corrector is along the normal at the elastic trial point, which avoids considering the intersection between the predicting increments of elastic stress and the damage surfaces. Furthermore, this algorithm ensures that the plastic and damage consistent conditions are fulfilled at any stage of the loading process. For details on the numerical scheme and the associated algorithmic steps, refer to the published reports [10, 28].

The following notational convention is used in this section: (1) represents at the trial point in time step . (2) For simplicity, the increment is defined as .

The elastic trial stress is defined as follows:

By solving (16), (17), and (21) in terms of the trial stress, the increments of the equivalent plastic strains and , plastic strain , and damage variables and can be obtained:

The expression for the increments of the plastic strain can be obtained by substituting (27b) into (23):

Substituting (28) into (27a) yields the relationship between the stress and the damage variables, written in incremental form as follows: where the fourth-order tensor is the constitutive matrix, which is the tangent with respect to the plasticity and the secant with respect to damage

Substituting (29) into (27c) and (27d) yields the increments of the damage variables and . Back substituting the increments of the damage variables into (29) results in the stress increment .

3. Determination of Specific Functions For Concrete

Specific functions for concrete are proposed based on the general framework of the coupled plastic damage model given in the previous section.

Referring to the definition of the isotropic compliance tensor, the evolution directions of added compliance tensors can be defined as [25] where is the initial Young’s modulus, and is the initial Poisson’s ratio. The symmetrized outer product “” is defined as .

Generally, the parameters and are the functions of the damage scalars, respectively. Their mathematical expressions can be logarithmic, power, or exponential functions. The specific expressions of and are given as where and are, respectively, the initial damage energy release threshold under tension and compression, and are the parameters controlling the damage evolution rate under tension, and and are the parameters controlling the damage evolution rate under compression.

The plastic yield function adopted in this work was first introduced in the Barcelona model by Lubliner et al. [18] and was later modified by Lee and Fenves [7], Wu et al. [9], and Voyiadjis et al. [11]. It has proven to be excellent in modeling the biaxial strength of concrete. However, the yield function is usually used in effective stress. To apply this function in the CPDM written in terms of true stress, Taqieddin et al. [17] modified the function by introducing damage parameters and achieved better results. Taqieddin et al. showed that using the same model and material parameters, fully coupled plastic damage yield criterion showed a more pronounced effect of damage on the mechanical behavior of the material than did the weakly coupled plastic damage model yield criterion. Considering that damage induces a reduction in the plastic hardening rate, the damage parameters are then introduced into the yield function by the plastic hardening functions of damaged material [ and ]. The yield function used in the present work is expressed as follows: where is the first invariant of the stress tensor . is the second invariant of the nominal deviatoric stress .

The parameter is a dimensionless constant given by Lubliner et al. [18] as follows: where and are the initial equibiaxial and uniaxial compressive yield stresses, respectively. Based on the experimental results, the ratio lies between 1.10 and 1.20; then, as derived by Wu et al. [9], the parameter is between 0.08 and 0.14. In the present work, is chosen as 0.12. The parameter is a parameter expressed from the tensile and compressive plastic hardening functions:

Because the associative flow rule is adopted in the present model, the plastic yield function is also used as the plastic potential to obtain the plastic strain.

According to (22), a specific expression of the plastic energy should be given to obtain the plastic hardening functions. Considering, as previously, a reduction in the plastic hardening rate due to damage, the plastic energy is expressed as follows: where the reduction factors are functions of the damage scalars , and is the plastic hardening energy for undamaged material, which is defined by in which and are the initial uniaxial tensile and compressive yield stresses, respectively. and are the tensile and compressive plastic hardening moduli, respectively. The hardening parameters and are the equivalent plastic strains under tension and compression, respectively, defined as [1, 17] where and are the tensile and compressive equivalent plastic strain rates, respectively.

According to Wu et al. [9], the evolution laws for may be postulated as follows: where and are the maximum and minimum eigenvalues of the plastic strain rate tensor , respectively. and are the maximum and minimum eigenvalues of the stress tensor . is a weight factor, expressed as in which the symbol is the Macaulay bracket, defined as . Note that is equal to one if all of the eigenstresses are positive, and it is equal to zero if all of the eigenstresses are negative.

The plastic tensile and compressive hardening functions of the damaged material are then specified as

It can be observed from (33) and (41) that the damage variables are introduced into the plastic yield function.

The last expressions to be specified are the reduction factors . Defined as the sum of the principal damage values in three directions (, , and ), the reduction factor is obtained from the definition of a particular thermodynamic energy function in the present paper. Based on the previous works of Lemaitre and Desmorat [29], the following particular form of the elastic thermodynamic energy function is used: where is the hydrostatic stress, and , and , , and are the damage values in , , and directions, respectively. The effective damage tensor is given by

Differentiation of the elastic thermodynamic energy yields the strain-stress relations where is the component of the deviatoric tensor.

Substituting (44) into (6), the relation between the second-order damage tensor and the compliance tensor of the damaged material can be given by

By solving (45), the reduction factor can be obtained. For example, considering that concrete is subjected to biaxial tension loading ( and ), the reduction factor can be obtained

Similarly, , corresponding to compression loading, can be calculated

In the special case of uniaxial loadings or combined tension-compression loading (), the reduction factors can be expressed as

4. Numerical Examples and Discussion of Results

In this section, the procedure for the determination of the model parameters is first presented. Then, several numerical examples are provided to investigate the capability of the proposed model in capturing material behavior in both tension and compression under uniaxial and biaxial loadings. The examples are taken from [7], with corresponding experimental data provided by Kupfer et al. [30], Karson and Jirsa [31], Gopalaratnam and Shah [32], and Taylor [33]. To comply with the results of the studies [7, 9], these numerical examples were analyzed using the single quadrilateral finite element shown in Figure 1.

(a)

(b)

(c)

(d)

(a)
(b)
(c)
(d)

Figure 1 

Quadrilateral finite elements under (a) uniaxial tension, (b) biaxial tension, (c) uniaxial compression, and (d) biaxial compression.

4.1. Identification of Model Parameters

The proposed model contains 12 parameters: four elastic constants for the undamaged material (, , , and ), two parameters for plasticity characterization ( and ), and six parameters for damage characterization (, , , and ). All of the parameters can be identified from the uniaxial tension and compression tests.

The initial elastic constants (, , , and ) are determined from the linear part of the typical stress-strain curve in Figure 2.

(a)

(b)

(a)
(b)

Figure 2 

Schematization of typical stress-strain curves of concrete in uniaxial (a) tension and (b) compression tests.

As mentioned before, the damage evolution is coupled with the plastic flow. According to Shao, it is reasonable to assume that the damage initiation occurs at the same time as the plastic deformation. They are identified as the end of the linear part of the stress-strain curves (point in Figure 2). In this case, considering and in the damage criterion and the plastic yield function, the initial damage thresholds and can be determined as follows: where the stress ratio is expressed as .

The two parameters ( and ) characterizing the failure surface and the plastic hardening parameter () can be determined by fitting both curves of stress-strain and stress-plastic strain obtained in the uniaxial tension test. Similarly, the parameters (, , and ) are determined by fitting both curves of stress-strain and stress-plastic strain obtained in the uniaxial compression test.

The effects of the model parameters (, , and ) on the stress-strain response, the stress-plastic strain response, and the damage evolution under uniaxial tension are shown in Figure 3. Each model parameter in turn is varied, while the others are kept fixed, to show the corresponding effect on the behavior of the concrete material. A good agreement exists between the experimental data [33] and the numerical simulations obtained with the associated material parameters given in Table 1. It can also be observed in Figure 3 that not only is the stress-strain curve governed by all of the parameters (, , and ) of the model but also the plastic strain and damage evolution are controlled by all of the parameters. This result is in agreement with that described by the proposed damage energy release rate and the plastic yield criterion in which both the damage scalar and the plastic hardening rate are introduced. These responses observed in Figure 3 show the coupled effect of damage and plasticity on the predicted response.

(a)

(b)

(c)

(a)
(b)
(c)

Figure 3 

Effect of the parameters on the model response in tension: (a) effect of , (b) effect of , and (c) effect of .

To illustrate the effects of the model parameters (, , and ) on the behavior of the proposed model, the stress-strain curve, the stress-plastic strain curve, and the damage evolution under uniaxial compression are plotted in Figure 4 in which the results obtained from varying each parameter while keeping the others fixed are also presented. The model parameters obtained from the experimental data [31] are listed in Table 1. It can be observed in Figure 4 that the numerical predictions, including not only the hardening and softening regimes of the stress-strain curve but also the stress-plastic strain curve, are in good agreement with the experimental data [31]. It can also be observed that all of the responses of the model, including the stress-strain curve, the stress-plastic strain curve, and the damage evolution, are controlled by all of the parameters (, , and ). This finding corresponds to the expressions of the proposed damage energy release rate and the plastic yield criterion in which both the damage scalar and the plastic hardening rate are introduced. These responses observed in Figure 4 again indicate the coupled effect of damage and plasticity on the predicted behavior.