Congreso de Métodos Numéricos en Ingeniería 2005 Granada, 4 a 7 de Julio, 2005 © SEMNI, España 2005

AN ENERGY BASED PLASTIC-DAMAGE MODEL FOR CONCRETE Jian-Ying Wu1, Jie Li 2, Rui Faria3* 1: Department of Civil Engineering, College of Architecture & Civil Engineering South China University of Technology Guangzhou, Post code: 510641, China e-mail: [email protected] 2: Department of Building Engineering Tongji University Shanghai, China 3: Departamento de Engenharia Civil Faculdade de Engenharia da Universidade do Porto Rua Dr. Roberto Frias, 4200-465 Porto, Portugal e-mail: [email protected]; web: http://www.fe.up.pt/labest

Keywords: Concrete Structures, Damage Mechanics, Plastic Deformations Abstract. A concrete devised constitutive model is presented, where two damage variables leading to a fourth-order damage tensor describe the degradation of the macro-mechanical material properties. Based on the decomposition of the effective stress tensor, the elastic free energy is defined. Evolution of the plastic strains and the expression for the elastoplastic free energy are also determined. The damage energy release rates conjugated to the damage variables are thus derived to establish the damage criteria, on which are based the damage evolution laws. Several applications to concrete specimens and structures are presented, demonstrating the model capabilities under different load conditions.

1. INTRODUCTION When designing and analyzing concrete structures under severe loading conditions, material failure should appropriately be taken into account. Comparing to the extremely developed computational hardware, accurate numerical modelling of the nonlinear behaviour of concrete is still far from being settled, which hinders further insight on the real performance of such structures. With plenty of documented experimental tests and research efforts, it is now generally accepted that the observed nonlinearities in concrete behaviour are mainly attributed to the following microstructural changes, which should be addressed in any appropriate constitutive model: (i) the development of microcracks and microvoids and (ii) the plastic flow along some preferred crack lips. Less phenomenological assumptions would be required in micromechanically derived models, but the inherent complexities and numerical difficulties

Jian-Ying Wu, Jie Li and Rui Faria

of the latter restrict their applicability. Therefore, combining theories of continuum damage mechanics (CDM) and plasticity [1], the plastic-damage approach is a natural way to the constitutive modelling of concrete structures. One of the critical issues in applying CDM to concrete is selecting the appropriate damage criteria and damage variables. Several different criteria, such as the equivalent strain-based ones, the stress-based ones, as well as the damage energy release rate-based proposals, are generally adopted [2-6]. Due to the inherent deficiencies of the former two approaches, and the fact that the damage energy release rates (DERR) are the conjugated forces to the corresponding damage variables, DERR-based damage criteria are considered more consistent thermodynamically. However, the elastic DERR-based unilateral damage model of Mazars [4] can not predict the strength enhancement observed in biaxial compression, and the elastoplastic DERR-based one-scalar damage model suggested by Ju [6] is only suitable under very simple situations (like Von Mises plasticity with linear isotropic hardening). Besides, one-scalar damage models are somewhat too limited to describe the change of Poisson’s ratio due to microcracking, and the unilateral effect inherent to concrete behaviour under cyclic loadings. Therefore, at least a fourth-order damage tensor should be employed even for characterizing the isotropic damage [7], although for applications to concrete structures such approach is too complex for being handled by practical engineers. To describe the unilateral effect and the dissimilar behaviour of concrete under tension and compression, some split of variables like the Cauchy stress tensor σ [2,4,8], the strain tensor ε [5-6] or the effective stress tensor σ [9-10] are usually adopted. In a strain-driven algorithm σ is the unknown to be computed from the constitutive model, which leads the decomposition of the Cauchy stress tensor to introduce computational difficulties; conversely, the split of the strain tensor is consensually considered as the less appropriate strategy to reproduce the macromechanical behaviour of concrete. As it will be demonstrated, based on the effective stress space plasticity [6] the split of σ can be easily and efficiently implemented in a strain-driven finite element implementation, and it provides quite realistic predictions of the concrete performance. In this paper a new energy based plastic-damage model for concrete is introduced. Two damage scalars d+ and d−, describing the two basic degradation mechanisms of concrete under tension and compression, are consistently referred to in the derivation of the constitutive model to be presented herein. Based on the decomposition of the effective stress tensor the elastic Helmholtz free energy (HFE) is defined, and within the framework of CDM an elegant plastic-damage constitutive relation with a fourth-order damage tensor and a plastic strain tensor is derived. Regarding the specific format for the effective stress space plasticity, the evolution law for the plastic strains and the explicit expression for the elastoplastic free energy are determined. The DERR that are conjugated to the damage variables are thus derived to establish the damage criteria, from which the evolution laws for the damage variables are obtained in accordance with the normality rule. Several applications of the proposed model to concrete specimens and structures under different load conditions are presented at the end of the paper, demonstrating its validity and capability. 2

Jian-Ying Wu, Jie Li and Rui Faria

2. PLASTIC-DAMAGE MODEL 2.1. Effective stress tensor and elastic potential The effective stress tensor σ in a damaged material may be assumed to follow the classical elastoplastic format [6,9]

σ = C0 : ε e = C 0 : ( ε − ε p )

(1)

εe = Λ0 : σ

(2)

or equivalently

where C0 denotes the usual fourth-order linear-elastic constitutive tensor ( Λ 0 = C0−1 ), and ε e and ε p are the rank-two elastic and plastic strain tensors. To describe the different nonlinear performances of concrete under tension and compression, σ and its rate tensor σ& are decomposed into their positive and negative components as

σ+ = P+ : σ

σ− = σ − σ+ = P− : σ

(3)

σ& + = Q + : σ&

σ& − = σ& − σ& + = Q − : σ&

(4)

with the fourth-order projection tensors P ± and Q ± expressed as [11] P + = ∑ H (σ i ) ( pii ⊗ pii )

P− = I − P+

(5)

Q− = I − Q+

(6)

i

Q+ = P+ + 2

3

〈σ i 〉 − 〈σ j 〉

i =1, j > i

σi −σ j

∑

pij ⊗ pij

Here I is the fourth-order identity tensor, H (σ i ) denotes the Heaviside function of the ith eigenvalue σ i of σ , the Macaulay brackets 〈⋅〉 are such that 〈 x〉 = ( x + | x | / 2) , and the second-order symmetric tensor pij depends on the unit eigenvectors ni and nj of σ [11] in the form

pij = p ji =

1 ( ni ⊗ n j + n j ⊗ ni ) 2

(7)

The undamaged elastic HFE potential ψ 0e , defined here as the elastic strain energy, can be written as the addition of its positive and negative components (ψ 0e + ,ψ 0e − )

3

Jian-Ying Wu, Jie Li and Rui Faria

1 2

1 2

1 2

ψ 0e (ε e ) = σ : ε e = σ + : ε e + σ − : ε e = ψ 0e + (ε e ) + ψ 0e − (ε e ) 1 2

1 2

1 2

ψ 0e ± (ε e ) = σ ± : ε e = σ : (P ± : Λ 0 ) : σ = ε e : ( C0 : P ± ) : ε e

(8)

where subscript “0” indicates “initial” states, “ e ” refers to “elastic” and “±” denotes “+” or “−“ as appropriate. Considering the assumed two damage mechanisms in tension and compression, an elastic HFE potential with the form

ψ e (ε e , d + , d − ) = ψ e + (ε e , d + ) + ψ e − (ε e , d − )

(9)

ψ e ± (ε e , d ± ) = (1 − d ± )ψ 0e ± (ε e )

(10)

is postulated

2.2. Plastic strain tensor and effective stress space plasticity To compute σ the so-called “effective stress space plasticity” [6] is resorted to establish the evolution laws for the plastic strains [10]

ε& p = λ&∂ σ F p κ& = λ&H F (σ, κ ) ≤ 0, λ& ≥ 0, λ& F (σ, κ ) ≤ 0

(11)

where F and F p are the plastic yield function and the plastic potential, λ& is the plastic flow parameter, κ and H denote the vectorial hardening parameter and function, and ∂ x y = ∂y ∂x . The rate form of the relation between the effective stress and the strain tensors is thus deduced by standard procedures in classical plasticity, leading to σ& = Cep : ε&

(12)

where Cep is the effective elastoplastic tangent tensor C0  ep C = C0 : ∂ σ F p ) ⊗ ( C0 : ∂ σ F ) ( C0 − ∂ σ F : C0 : ∂ σ F p − ∂ κ F ⋅ H 

if λ& = 0 if λ& > 0

(13)

Here the plastic potential F p is identified with the Drucker-Prager function, expressed as [2,12] F p = α p I1 + 2 J 2

(14)

where I1 is the first invariant of σ , J 2 is the second invariant of s (the deviatoric component

4

Jian-Ying Wu, Jie Li and Rui Faria

of σ ), and α p ≥ 0 is a parameter selected to provide proper dilatancy (within the range of 0.2-0.3 for concrete). Combining Eqs (111) and (14), the evolution law for the plastic strains is obtained

ε& p = λ& ( s s + α p 1)

(15)

s = s : s = 2J 2

(16)

where

is the norm of s , and 1 is the second-order identity tensor. If one introduces the hardening parameters κ + and κ − as “equivalent plastic strains” under uniaxial tension and compression, with the definitions p p κ + = ∫ ε&max ε&max dt

p p κ − = ∫ ε&min ε&min dt

(17)

p p and ε&min are the maximum and minimum eigenvalues of ε& p , the evolution law for where ε&max κ may be postulated as p p κ& = {wκ& + , (1 − w ) κ& − } = {wε&max , − (1 − w ) ε&min } T

T

(18)

with w being a weighting factor expressed as [12] 3

w = ∑ σi i =1

3

∑σ i =1

(19)

i

that renders “1” for pure tension and “0” for pure compression. Calling for Eq (112) H is thus derived [13]

{

H = { H + , H − } = w∂σ i ,max F p , −(1 − w)∂σ i ,min F p T

}

T

(20)

where σ i ,max and σ i ,min are the maximum and minimum values of σ i . Any yield function F appropriate for concrete can be adopted in Eq (113); here a modified function proposed in [12,14] is selected

(

)

F (σ, κ ) = α I1 + 3 J 2 + β 〈σ i ,max 〉 − (1 − α ) c

(21)

with the following definitions for α, β and c

α = (ϑ − 1) ( 2ϑ − 1)

β (κ ) = (1 − α ) f − (κ ) f + (κ ) − (1 + α ) c ( κ ) = f − (κ )

(22)

where ϑ is the ratio between the yield strengths under equibiaxial and uniaxial compression (usually with values in the interval 1.1-1.2). f + (κ ) and f − (κ ) denote evolution stresses 5

Jian-Ying Wu, Jie Li and Rui Faria

(positive) in the effective stress space due to the plastic hardening/softening under uniaxial tension and compression, respectively, with the following simplified linear isotropic definitions f ± ( κ ) = f y± + E p ±κ ±

(23)

where f y+ and f y− are the effective yield strengths under uniaxial tension and compression, approximated as f y+ = ft and f y− = f c , with ft and f c being the uniaxial tensile and compressive strengths of concrete; E p ± are the effective plastic hardening modulus under uniaxial tension and compression, relating to the elastoplastic tangent modulus E ep ± as follows [15] E ep ± =

 E0 E p ± 1 = 1 − p± ± E0 + E  1 + RE

  E0 

(24)

with RE± = E p ± E0 denoting the ratios between E p ± and the initial Young’s modulus E0 . 2.3. Elastoplastic Helmholtz free energy potential Since the contribution to the initial plastic HFE potential ψ 0p due to plastic strains from tension is much smaller than from compression, in the present model ψ 0p + is neglected, and the following definition is adopted [13]  3 p − 1 + − 1 + ν 0  −  3J 2 + − ψ 0p (κ ) = ψ0p− (κ ) =  b0 − α I 3 J I1 I1  1 2  3 E 2 2 0   

(25)

where I1± are the first invariants of σ ± , J 2− is the second invariant of s − , ν0 is the initial Poisson’s coefficient and parameter b0 > 0 is a material property dependent on ϑ . Assuming the plastic flow in concrete to be mainly a consequence of the shear damage mechanism that prevails in compression [16], the plastic HFE potential is supposed to have the form

ψ p (κ , d − ) = (1 − d − )ψ 0p (κ )

(26)

The total elastoplastic HFE potential can be defined as the sum of the elastic component ψ e expressed in Eq (9) and the plastic component ψ p just presented, assumed uncoupled in the following manner

ψ (ε e , κ , d + , d − ) = ψ e (ε e , d + , d − ) +ψ p (κ , d − ) = ψ + (ε e , d + ) +ψ − (ε e , κ , d − ) with the components ψ ± being

6

(27)

Jian-Ying Wu, Jie Li and Rui Faria

ψ + (ε e , d + ) = (1 − d + )ψ 0+ (ε e ) ψ − (ε e , κ , d − ) = (1 − d − )ψ 0− (ε e , κ )

(28)

and the components ψ 0± of the initial elastoplastic HFE taken as

ψ 0+ = ψ 0e + = ( σ + : Λ 0 : σ ) 2

(

ψ 0− = ψ 0e − + ψ 0p − = b0 α I1 + 3 J 2

)

2

(29)

2.4. Constitutive law and damage energy release rates Referring to standard thermodynamics [17], and considering the HFE potentials defined in Eqs (27) and (9), the plastic-damage constitutive law can be established as [10,18] σ=

∂ψ e = ( I − D) : σ ∂ε e

(30)

with the fourth-order damage tensor D being D = d +P+ + d −P−

(31)

It can be clearly seen that Eq (30) is the standard relation in CDM between the Cauchy stress tensor σ and the effective stress tensor σ , according to the “strain equivalence” assumption [16], and Eq (31) is coincident with the conclusion drawn in [7]. Also, for any admissible process the following dissipation condition due to damage evolution has to be fulfilled

∂ψ  ∂ψ  γ& d = −  + d& + + − d& −  = Y + d& + + Y − d& − ≥ 0 ∂d  ∂d 

(32)

where the DERR Y ± , conjugated to the corresponding damage variables, are expressed as

Y+ = −

∂ψ = ψ 0+ = ( σ + : Λ 0 : σ ) 2 ≥ 0 + ∂d

Y− = −

∂ψ = ψ 0− = b0 α I1 + 3J 2 − ∂d

(

)

2

≥0

(33)

These equations demonstrate that the DERR depend on the initial elastoplastic HFE components ψ 0± , and not just on the elastic ψ0e ± ones. Accordingly, it is incorrect to establish the damage criteria taking the elastic DERR alone into consideration, disregarding contribution from plastic strains. 2.5. Damage criteria

With the DERR defined in Eq (33), the state of damage can then be characterized by means of damage criteria with the following form 7

Jian-Ying Wu, Jie Li and Rui Faria

G ± (Y ± , r ± ) = Y ± − r ± ≤ 0

(34)

where r ± are current damage thresholds (energy barriers) controlling the size of the expanding damage surfaces. Here it will be further assumed that r ± = r ± ( µ ± ) , with µ ± being related to the cumulative damages [2], with evolution laws to be defined later. Denoting by r0± = r ± (0) the initial damage thresholds, according to Eq (34) onset of damages implies Y ± = r ± ≥ r0± . Correspondingly, from Eqs (33) the initial damage thresholds are calculated as r0+ = ( f 0+ )

2

( 2 E0 )

r0− = b0 (1 − α ) f 0− 

2

(35)

where f 0+ and f 0− are the stresses (positive) beyond which nonlinearity becomes visible under uniaxial tension and compression, respectively. Under those stress states where the compressive consolidation mechanism can play a significant role, the model predictions when Y − is expressed as in Eq (332) are somewhat conservative, since influence of the third invariant of σ is completely disregarded. To cope with this influence definition of Y − is changed to [18] Y − = b0 α I1 + 3 J 2 − γ〈−σ i ,max 〉   

2

(36)

where

γ = 3 (1 − K c ) ( 2 K c − 1) and K c (denoting the ratio of

(37)

J 2 under the tensile and compressive meridian stress states

for any given value of the hydrostatic pressure I1 ) is assumed as a constant, with a typical value of 2/3 for concrete [14], leading to the value γ = 3 adopted hereinafter. Since influence of coefficient γ disappears in stress states other than triaxial compression (because σ i ,max ≥ 0 ), Eq (36) can be viewed as a minor modification of Eq (332) to improve the predictive capability of the proposed model under stress states other than the plane stress conditions. Under biaxial stress states the damage criteria defined through Eqs (331) and (34-36) lead to the linear domain shown in Figure 1, which agrees rather well with the experimental results from Kupfer et al. [19]. 2.6. Characterization of damage

Applying the normality rule to the damage criteria expressed in Eq (34), the following damage evolution laws can be obtained ± ± ± ± ∂G & d = µ& (r ) = µ& ± (r ± ) ± ∂Y

(38)

This equation states that when µ& ± > 0 further damage takes place, a “loading” scenario that 8

Jian-Ying Wu, Jie Li and Rui Faria

according to Eq (34) is ruled by condition G ± (Y ± , r ± ) = 0

r± = Y ±

⇒

(39)

σ2

Kupfer et al. (1969) Present Model

0.2

f 0−

σ1 0.0 − -1.4 -1.2 -1.0 -0.8 -0.6 -0.4 -0.2 0.0 0.2 f 0 -0.2 -0.4 -0.6

σ1 = σ 2

-0.8 -1.0 -1.2 -1.4

Figure 1. 2D linear domain.

Invoking the damage consistency condition one has G& ± (Y ± , r ± ) = 0

⇒

r& ± = Y& ±

(40)

∂µ ± & ± ∂µ ± ± Y = ± r& = h ± r& ± ± ∂r ∂r

(41)

and ∂r ± r& ± − Y& ± = ± µ& ± − Y& ± = 0 ∂µ

µ& ± =

⇒

Integrating Eqs (38) and (40) one reaches quite clear functional dependences for updating the damage variables and thresholds, that is, d ± = µ ± (r ± )

(42)

  r ± = max r0± , max (Yτ± )  τ∈[ 0, t ]  

(43)

In the present paper the following explicit functions are adopted for d ± = µ ± (r ± ) [9,13]:

d + = µ + (r + ) = 1 −

r0+ r+

   r +    + + +  1 A A exp B 1 − + −   ( )  +    r  0      

−   r −     r0 − − − d = µ (r ) = 1 −  − (1 − A ) + A exp  B  1 − −    r0       r  −

−

−

9

(44)

Jian-Ying Wu, Jie Li and Rui Faria

Eq (441) is capable to describe the tension-stiffening effect due to the interaction between steel rebars and concrete, if any. According to experimental tests [20], in reinforced concrete parameters A+ and B + should be related to the reinforcement ratio ρ s as [13] A+ = 1 − cs ρ s db

(

B + = 270

)

A+ f 0+ E0 ≤ 1000 f 0+ E0

(45)

where db is the rebar diameter (in mm), and cs is a dimension generally taking the value of 75mm. For plain concrete Eq (451) renders A+ = 1, thus Eq (441) reduces to the function first proposed by Oliver et al. in [21], where parameter B + was defined as −1

 G E 1 B =  f +0 2 −  ≥ 0  lch ( f 0 ) 2  +

(46)

with G f being the tensile fracture energy of concrete, and lch is a geometrical characteristic length of the finite element, introduced to minimize dependence of the numerical solution on the mesh refinement. It can be clearly seen that parameter b0 plays no role on Eq (442), so b0 = 1 can be adopted to compute r0− and Y − through Eqs (352) and (36). Parameters A− and B − of Eq (442) may be determined by imposing the numerical 1D σ − ε curve in compression to fit the one obtained from an experimental test. 3. COMPUTATIONAL ASPECTS

For the numerical implementation of the proposed plastic-damage model, the constitutive law in Eq (30) should be differentiated with respect to time, leading to

(

σ& = ( I − ω ) : C0 : ( ε& − ε& p ) − σ + d& + + σ − d& −

)

(47)

with the symmetric fourth-order tensor ω being ω = d +Q + + d −Q −

(48)

According to the concept of operator split [15] Eq (47) can be decomposed into elastic, plastic and damage parts, leading to the inherent numerical algorithm to include elastic-predictor, plastic-corrector and damage-corrector steps, as established in Table 1. It can be seen that during the elastic-predictor and plastic-corrector steps the damage variables are fixed, so the equations in the first two columns, constituting a standard elastoplastic problem in the effective stress space, are decoupled from the ones in the third column, concerning the damage part. Once σ is updated in the elastic-predictor and plastic-corrector steps, the damage variables d ± and the Cauchy stress σ can thus be updated with no more iterations in

10

Jian-Ying Wu, Jie Li and Rui Faria

the damage-corrector step, which is a significant advantage of the effective stress tensor-decomposed approach presented herein. Elastic-predictor

Plastic-corrector

Damage-corrector

ε& = ∇ s u& (t )

ε& = 0

ε& = 0

ε& p = 0

λ&∂ F p if F = F& = 0 ε& p =  σ otherwise 0

ε& p = 0

κ& = 0

κ& = λ&H

κ& = 0

σ& = C0 : ε&

σ& = −C0 : ε& p

σ& = 0

d& ± = 0

d& ± = 0

r& ± h ± if G ± = G& ± = 0 d& ± =  otherwise 0

r& ± = 0

r& ± = 0

Y& ± r& ± =  0

σ& = ( I − ω ) : C0 : ε&

σ& = − ( I − ω ) C0 : ε& p

σ& = − σ + d& + + σ − d& −

(

if G ± = G& ± = 0 otherwise

)

Table 1. Numerical algorithm (operator split).

4. APPLICATIONS

To validate and illustrate the applicability of the proposed model, several numerical examples of concrete specimens and structures under various load conditions are presented in this section. For all cases it is assumed that ν0 = 0.20, ϑ = 1.16 and α p = 0.20. 4.1. Monotonic uniaxial tests

The experimental results from a monotonic 1D tensile test [22] and a monotonic 1D compression one [23] are taken as the first examples to check the adequacy of the numerical model. The material properties used in the simulations are: for the tensile test, E0 = 38GPa, f0+ = 3.4MPa, Gf = 70N/m; for the compression one, E0 = 31.7GPa, f0− = 10.2MPa. Figure 2 compares the predicted σ − ε numerical curves with those obtained from the experimental tests. For both tests predictions from the numerical model agree well with the experimental results, either in the hardening or in the softening branches.

11

Jian-Ying Wu, Jie Li and Rui Faria

Test Result Present Model

4.0 3.5

25

Stress(MPa)

Stress(MPa)

Test Result Present Model

30

3.0 2.5 2.0 1.5

20 15 10

1.0 5

0.5 0.0 0.0

0.5

1.0

1.5

2.0

2.5

3.0

3.5

4.0

0

0

1

Strain(1.0E-4)

2

3

4

5

6

Strain(1.0E-3)

a) Tension [22]

b) Compression [23] Figure 2. Monotonic uniaxial tests.

4.2. Monotonic biaxial tests

The proposed model is also validated with the results under biaxial compression ( σ 3 = 0) reported in Kupfer et al. [19]. The material properties adopted in the simulation are: E0 = 31GPa, f0+ = 3MPa, f0− = 15MPa, Gf = 75N/m. For the tests under the stress ratios σ 2 σ 1 = −1/0, σ 2 σ 1 = −1/−1 and σ 2 σ 1 = −1/−0.52 the predicted stress-strain curves illustrated in Figure 3a agree well with the experimental ones, capturing the overall nonlinear behaviour. σ 2 (MPa)

50

Test Results

40

-1.4

σ 2 σ 1 = −1/ − 0.52

30

σ2

Test Results Present Results Model

Present Model

-1.2 -1.0 -0.8

0.2 0.0 -0.6 -0.4 -0.2 0.0 -0.2

fc

σ1 0.2

fc

-0.4

σ 2 σ 1 = −1/ 0

20

-0.6

σ 2 σ 1 = −1/ − 1

-0.8 -1.0

10

ε2

0 0

1

2

3

4

5

-1.2 -1.4

6

a) Stress-strain curves

b) Biaxial strength envelope

Figure 3. Monotonic biaxial compression tests [19].

To further illustrate the capability of the proposed model, the numerical biaxial strength 12

Jian-Ying Wu, Jie Li and Rui Faria

envelope is reproduced in Figure 3b as well, which is found to be almost coincident with the experimental one reported in [19]. As clearly perceptible in Figure 3b, another important attribute of the present model is its ability to predict not only the enhancement of concrete strength under biaxial compression, but also the reduction on the compressive strength induced by orthogonal tensile cracking under tension-compression stress states. 4.3. Monotonic 3D compression test

To check the capability of the proposed model for concrete under compressive confinement, the numerical predictions are compared with the experimental results reported in [24]. To fit the experimental uniaxial compression stress-strain curve the material properties adopted in the simulation are: E0 = 37GPa, f0+ = 4MPa, f0− = 15MPa. The predicted results under two sets of confining stresses, namely σ1 = σ2 = −6.895MPa and σ1 = σ2 = −13.79MPa, are reproduced in Figure 4. It can be clearly seen that the enhancement of strength and ductility due to the compressive confinement, as well as the overall stress-strain experimental curves, are satisfactorily predicted by the numerical model. Test Results Present Model 120

σ 3 (MPa) σ 1 = σ 2 =-13.79MPa

100 80

σ 1 = σ 2 =-6.895MPa

60 40

σ 1 = σ 2 =0.0MPa

20 0 0.000

0.002

0.004

0.006

0.008

0.010

0.012

ε3 Figure 4. Monotonic 3D compression tests [24].

4.4. Cyclic uniaxial tests

In Figure 5 the cyclic uniaxial tensile test of [25] and the cyclic compression test of [23] are reproduced numerically to demonstrate the capability of the proposed model under cyclic load conditions. The adopted material properties are: for the tensile test, E0 = 31.7GPa, f0+ = 3.5MPa, Gf = 24N/m; for the compression test, E0 = 30GPa, f0− = 15MPa. As shown in Figure 5, under both tension and compression, the experimentally observed softening and stiffness degrading, as well as the irreversible strains upon unloading, are well reproduced by the proposed constitutive model.

13

Jian-Ying Wu, Jie Li and Rui Faria

Test Result Present Model

4.0

30

3.5

Test Result Present Model

25 20

2.5

Stress(MPa)

Stress (MPa)

3.0

2.0 1.5

15 10

1.0 5

0.5 0.0 0.0000

0.0001

0.0002

0.0003

0.0004

0.0005

0 0.000

0.001

0.002

0.003

0.004

0.005

Strain

Strain

a) Tension [25]

b) Compression [23] Figure 5. Cyclic uniaxial tests.

4.5. Earthquake analysis of Koyna dam

The example of Koyna dam, stroke by the transverse and vertical components of a well-known earthquake recorded in 1967, is analyzed here. The dam-foundation interaction is ignored by assuming a rigid foundation, and the dam-reservoir dynamic interaction is modelled with the simplified added mass technique of Westergaard. The material properties for the concrete are: density ρ = 2643kg/m3, E0 = 31GPa, f0+ = 2.9MPa, Gf = 200N/m and fc = 24.1MPa. A Rayleigh stiffness-proportional damping matrix was assumed to provide a 3% fraction of the critical damping for the first vibration mode of the dam, and the HHT α-method was adopted to integrate the dynamic equation of motion. Horizontal Relative Displacement (mm)

Lee and Fenves (1998) Present Model 50 40 30 20 10 0 0

2

4

6

8

10

-10 -20 -30 -40

Time (s)

a) FE mesh

b) Horizontal relative crest displacement Figure 6. Earthquake analysis of Koyna dam.

14

c) Final pattern of d+

Jian-Ying Wu, Jie Li and Rui Faria

Figure 6a reproduces the FE mesh adopted for the numerical analyses. The predicted historical evolution of the horizontal displacements at the dam crest is depicted in Figure 6b, which agrees well with the one reported by Lee and Fenves in [12]. Figure 6c reproduces the tensile damages on the dam at the end of the earthquake: at the elevation where the downstream face inclination changes abruptly a damaged region is predicted, in agreement with the major cracks observed after the earthquake and reported in [26]. 4.6. Reinforced concrete shearwalls

The proposed model was also used to simulate the nonlinear performances of thirteen large-scale reinforced concrete shearwalls, tested under various axial and monotonically increasing horizontal forces. The steel configurations, the vertical loadings and other details of the experimental tests were reported in [27, 28]. The material properties of concrete adopted are: E0 = 30.7GPa, f0+ = 3MPa and f0− = 10MPa. In the simulations 8-node plane stress elements are used for discretizing the concrete, whereas 2-node truss elements are adopted for the steel reinforcement. Numerical predictions agree well with the experimental results, as demonstrated in Figure 7, which reproduces the horizontal force-displacement curve registered at the top of two shearwalls (SW11 and SW14). Distributions of the tensile damages at SW11 under horizontal loadings of 122.0kN and 153.1kN are shown in Figure 8. For the thirteen shearwalls the peak loads and displacements computed with the proposed model are compared with the experimental ones in Table 2. It can be clearly seen that the peak forces F and displacements U of all the shearwalls are reproduced rather well by the proposed model: (i) the ratio of the numerically predicted peak loads to the experimental ones has a mean value of 0.998 and a coefficient of variation of 0.042, and (ii) the corresponding values for the displacement ratios are 1.095 and 0.131, respectively. Test Result of SW11 Test Result of SW14 Present Model

Horizontal Load (kN)

300 250 200 150 100 50 0 0

2

4

6

8

Horizental Displacement (mm)

Figure 7. RC shearwall tests [27].

15

10

12

Jian-Ying Wu, Jie Li and Rui Faria

a) Horizontal force: 122.0kN

b) Horizontal force: 153.1kN

Figure 8. Tensile damages in shearwall SW11.

Shearwal l SW11 SW12 SW13 SW14 SW15 SW16 SW17 SW21 SW22 SW23 SW24 SW25 SW26

Ftest

Fmodel

Fmodel

U test

U model

U model

(kN) 260 340 330 265 320 355 247 127 150 180 120 150 123

(kN) 259.6 328.4 352.0 259.6 321.5 354.9 241.8 121.7 148.9 170.2 121.7 164.3 120.1

Ftest

(mm) 9.8 9.8 8.9 11.2 8.8 5.8 10.8 20.6 15.3 13.2 18.1 9.5 20.9

(mm) 10.9 10.0 11.4 10.9 10.0 6.9 10.5 20.0 14.6 14.6 20.0 13.1 21.7

U test

0.998 0.966 1.067 0.980 1.005 1.000 0.979 0.958 0.993 0.946 1.014 1.095 0.976

1.107 1.019 1.284 0.969 1.129 1.197 0.972 0.970 0.954 1.107 1.103 1.383 1.036

Table 2. Peak horizontal forces and displacements.

5. CONCLUSIONS

This paper presents a constitutive model mainly intended for the nonlinear analysis of plain and reinforced concrete structures. Within the framework of continuum damage mechanics, two damage scalars that lead to a fourth-order damage tensor are adopted to describe the degradation of the macro-mechanical properties of concrete. The specific format of the effective stress space plasticity is used to determine the evolution laws for the plastic strains and the plastic Helmholtz free energy. The damage criteria are based on the elastoplastic damage energy release rates, from which the evolution laws for the damage variables are obtained in accordance with the normality rule.

16

Jian-Ying Wu, Jie Li and Rui Faria

The numerical predictions from the proposed model were checked in several applications to concrete specimens and structures, demonstrating that it is able to reproduce the typical nonlinear performances of concrete under monotonic and cyclic load conditions, as well as under 1D, 2D or 3D stress states. Finally, the model was also used for analysing the earthquake performance of Koyna dam, and to predict the ultimate behaviour of several reinforced concrete shearwalls, providing results that agree with the observations. 6. ACKNOWLEDGEMENTS

The support of the National Natural Science Foundation for Outstanding Youth in China (No. 59825105) is gratefully acknowledged. REFERENCES

[1] W. F. Chen, Constitutive equations for engineering materials. Plasticity and modeling, Elsevier, Vol. 2, (1994). [2] M. Ortiz, “A constitutive theory for inelastic behaviour of concrete”, Mech. Mat., Vol. 4, pp. 67-93, (1985). [3] C. L. Chow and J. Wang, “An anisotropic theory of continuum damage mechanics for ductile fracture”, Eng. Frac. Mech., Vol. 27, pp. 547-558 (1987). [4] J. Mazars, “A model of unilateral elastic damageable material and its application to concrete”, Proc. RILEM Int. Conf. Fracture Mech. Conc., Lausanne, Switzerland, (1985). [5] J. C. Simo and J. W. Ju, “Strain- and stress-based continuum damage models. I: formulation”, Int. J. Solids Struct., Vol. 23(7), pp. 821-840, (1987). [6] J. W. Ju, “On energy-based coupled elastoplastic damage theories: constitutive modeling and computational aspects”, Int. J. Solids Struct., Vol. 25(7), pp. 803-833, (1989). [7] J. W. Ju, “Isotropic and anisotropic damage variables in continuum damage mechanics”, J. Eng. Mech., ASCE, Vol. 116(12), pp. 2764-2770, (1990). [8] S. Yazdani and H. L. Schreyer, “Combined plasticity and damage mechanics model for plain concrete”, J. Eng. Mech., ASCE, Vol. 116(7), pp. 1435-1450, (1990). [9] R. Faria, J. Oliver and M. Cervera, “A strain-based plastic viscous-damage model for massive concrete structures”, Int. J. Solids Struct., Vol. 35(14), pp. 1533-1558, (1998). [10] J. Y. Wu and J. Li, “A new energy-based elastoplastic damage model for concrete”, Proc. XXI Int. Conf. Theor. Applied Mech., Warsaw, Poland, (2004). [11] R. Faria, J. Oliver and M. Cervera, “On isotropic scalar damage models for the numerical analysis of concrete structures”, CIMNE Monograph, N. 198, (2000). [12] L. Lee and G. L. Fenves, “Plastic-damage model for cyclic loading of concrete structures”, J. Eng. Mech. Div., ASCE, Vol. 124, pp. 892-900, (1998). [13] J. Y. Wu, Damage energy release rate-based elastoplastic damage constitutive model for concrete and its application to nonlinear analysis of structures, PhD Dissertation, Tongji University, Shanghai, China, (2004). [14] J. Lubliner, J. Oliver, S. Oller and E. Oñate, “A plastic-damage model for concrete”, Int.

17

Jian-Ying Wu, Jie Li and Rui Faria

[15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28]

J. Solids Struct., Vol. 25(3), pp. 299-326, (1989). J. C. Simo and T. J. Hughes, Computational inelasticity, Springer-Verlag, New York, (1998). J. Lemaitre, “A continuum damage mechanics model for ductile fracture”, J. Eng. Mat. Techn., Vol. 107, pp. 83-89, (1985). B. D. Coleman and M. E. Gurtin, “Thermodynamics with internal state variables”, J. Chem. Physics, Vol. 47, pp. 597-613, (1967). J. Li and J. Y. Wu, “Energy-based CDM model for nonlinear analysis of confined concrete structures”, Int. Symp. Confined Conc., (2004). H. Kupfer, H. K. Hilsdorf and H. Rusch, “Behaviour of concrete under biaxial stress”. J. Am. Conc. Inst., Vol. 66(8), pp. 656-666, (1969). N. J. Stevens, et al., “Constitutive model for reinforced concrete finite element analysis”, ACI Struct. J., Vol. 88(1), pp. 49-59, (1991). J. Oliver, M. Cervera, S. Oller and J. Lubliner, “Isotropic damage models and smeared crack analysis of concrete”, Proc. 2nd Int. Conf. Comp. Aided Analysis Design Conc. Struct., Zell am See, pp. 945-957, (1990). Q. Y. Zhang, Research on the stochastic damage constitutive of concrete material, PhD Dissertation, Tongji University, Shanghai, China, (2001). I. D. Karson and J. O. Jirsa, “Behaviour of concrete under compressive loadings”, J. Struct. Div., ASCE, Vol. 95(12), pp. 2535-2563, (1969). S. J. Green and S. R. Swanson, “Static constitutive relations for concrete”, AFWL-TR-72-244, U.S. Air Force Weapons Laboratory, (1973). R. L. Taylor, “FEAP: A finite element analysis program for engineering workstation”, Rep. No. UCB/SEMM-92, Dept. Civil Eng., University of California, Berkeley, (1992). A. K. Chopra and P. Chakrabarti, “The Koyna earthquake and the damage to Koyna dam”, Bulletin Seism. Soc. Am., Vol. 63(2), pp. 381–397, (1973). D. Lefas, M. D. Kotsovos and N. N. Ambraseys, “Behaviour of reinforced concrete structural walls: strength, deformation characteristics, and failure mechanism”, ACI Struct. J., Vol. 87(1), pp. 23-31, (1990). F. J. Vecchio, “Finite element modeling of concrete expansion and confinement”, J. Eng. Mech., ASCE, Vol. 118(9), pp. 2390-2406, (1992).

18