Coupled carbonation-rust formation-damage modeling and simulation of steel corrosion in 3D mesoscale reinforced concrete

Cement and Concrete Research 74 (2015) 95–107

Contents lists available at ScienceDirect

Cement and Concrete Research journal homepage: http://ees.elsevier.com/CEMCON/default.asp

Coupled carbonation-rust formation-damage modeling and simulation of steel corrosion in 3D mesoscale reinforced concrete T.T.H. Nguyen a, B. Bary a,⁎, T. de Larrard a,b a b

CEA, DEN, DPC, SECR, Laboratoire d'Etude du Comportement des Bétons et des Argiles, F-91191 Gif-sur-Yvette, France LMDC, UPS — INSA Toulouse, France

a r t i c l e

i n f o

Article history: Received 11 July 2014 Accepted 22 April 2015 Available online 16 May 2015 Keywords: E. Modeling C. Corrosion C. Finite element analysis C. Carbonation B. Microcracking

a b s t r a c t This paper presents a modeling strategy to simulate the corrosion of steel reinforcement in atmospheric environment due to carbonation of concrete. Its principal objectives are to analyze the effects of the progressive formation of corrosion products at the steel/concrete interface on concrete cover cracking. The approach is based on modeling studies carried out independently on carbonation, corrosion and creep. These models are coupled and integrated into a numerical 3D simulation procedure for investigating the behavior of concrete mesostructures. A viscodamage model is used to reproduce both creep and damage behaviors of mortar, and the approach is applied to the simulation of a 3D reinforced concrete mesostructure including explicitly the coarse aggregates. The numerical results highlight the influence of aggregates and the effects of creep on crack initiation and propagation. © 2015 Elsevier Ltd. All rights reserved.

1. Introduction During the service life of reinforced concrete structures, one of the main causes leading to mechanical damage in the medium and long term, is the corrosion of steel reinforcement. This process is inevitable because of the porous structure of the material which allows substances such as oxygen and carbon dioxide penetrating in the matrix concrete (see e.g. [1]). From the phenomenological point of view, the transformation of a portion of the steel section into rusts reduces its effective area and then its reinforcement properties. Furthermore, since corrosion products have a volume of 3 to 7 times greater than the initial volume of the steel lost (see e.g. [2–4]), the corresponding expansion may generate tensile stresses important enough to disrupt the concrete cover and initiate internal cracks that propagate from the steel/concrete interface to the exposed surface. The consecutive physical modifications (related both to mechanical and transport properties) then affect the durability of structure. For structures exposed to the atmospheric environment, corrosion relates mainly to the carbonation of concrete cover. Steel in concrete is protected initially against active corrosion by a film of iron oxides on the surface of the embedded steel which prevents the dissolution of iron. Once carbon dioxide in the atmosphere penetrates in the cementitious matrix and reacts with portlandite Ca(OH)2 and other hydrates as C–S–H, concrete interstitial solution is chemically modified (see e.g. [5]). When this process reaches the surface of steel and drops the concrete pH around reinforcement, the passive film disappears ⁎ Corresponding author. Tel.: +33 1 69 08 23 83; fax: +33 1 69 08 84 41. E-mail address: [email protected] (B. Bary).

http://dx.doi.org/10.1016/j.cemconres.2015.04.008 0008-8846/© 2015 Elsevier Ltd. All rights reserved.

and active corrosion of steel occurs. The kinetics of corrosion can then change over several orders of magnitude compared to the passive state. To take into account the corrosion in the calculation of long-term behavior of existing structures, most numerical models are developed at the macroscopic scale and then calibrated on experimental tests. This approach requires certain data about the elastic properties of the layers of corrosion products, such as the Young modulus E and Poisson ratio ν. However, these data are not yet well characterized, especially the value of the elastic modulus of rust. In the modeling, some authors consider that the corrosion products are mainly composed of water. For instance, in [6] the mechanical properties of water are adopted to represent those of the rust. The bulk modulus and Poisson ratio ν are taken to 2 GPa and 0.49, respectively. In [3] it is assumed that the corrosion products have the same properties as sound steel, and E and ν are equal to 210 GPa and 0.3, respectively. Other numerical studies calibrate the Young modulus by inverse analysis on finite element (FE) simulations such that the calculated time of appearance of scaling or amount of rust required for cracking coincides with that reported in the experimental tests. In [7, 8] this value varies between 0.06 GPa and 0.15 GPa for accelerated corrosion in laboratory, and similar values ranging between 0.1 GPa and 0.5 GPa have been found in [9]. On the other hand, for natural corrosion mechanical tests on steel samples corroded in the concrete of a seaport during a time span of 40 years have been carried out in [10], revealing a modulus between 0.1 GPa and 0.6 GPa. By contrast, in a recent study of archeological artifacts over 600 years [11], the authors obtain Young modulus by microindentation between 70 GPa and 200 GPa. At the microscopic scale, some authors assume that the properties of the corrosion products are those of their components. For instance, [12] propose a modulus

96

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

between 50 GPa and 200 GPa corresponding approximately to that of oxide crystal modulus. Similarly, [13] identified this parameter by oedometer tests completed by Hertz contact law on lepidocrocite and α-Fe2O3, and reported very high values between 300 GPa to 400 GPa. The experimental method also plays an important role: the last referenced study shows that the Young's modulus determined by ultrasonic measurements on lepidocrocite is two orders of magnitude smaller (0.6–5 GPa). Some data found in the literature are reported in Fig. 1. This short review raises naturally some questions about the mechanical properties of rust. Indeed, experimental studies reporting high values of modulus were performed on samples taken from structures under natural corrosion during tens or hundreds of years [11,13]. Corrosion products are in this case well-crystallized layers, and the qualitative characterization shows that rusts contain mainly magnetite, lepidocrocite and goethite. By contrast, for accelerated corrosion tests carried out in laboratory, lower values are obtained [8–10,13]. In [9], the authors estimate that the low-crystallized products (assumed to be FeOx(OH)3 − 2x with x equal to 0 or 1) occupy up to 45–65% of total weight of the corrosion products. During these artificial corrosion tests, it has also been observed that green rusts appear at first, then they transform into classical products such as lepidocrocite and akaganeite [14]. The existence of these unstable intermediate phases and green rusts in corrosion was cited by several studies, see e.g. [15–18]. These products contain iron oxides with valence II/III and oxidize easily to become different iron oxides with higher valence according to environmental conditions. As they progressively transform into other stable rusts, they are rarely observed in studies of archeological corrosion as noticed e.g. in [19]. However, their presence during natural corrosion in medium term seems to play a significant role in the mechanical behavior of corroded reinforced concrete because their mechanical characteristics are several orders of magnitude smaller than the ones of crystallized corrosion products. For this reason, we propose in this study to model the development of corrosion due to carbonation of concrete in atmospheric environment based on a simplified hypothesis on the mesoscopic composition of corrosion products. We assume this latter as a two-layer composite: a first layer with low-crystallized phases containing unstable products and possibly green rust, and a second stable layer composed of wellcrystallized rust, as experimentally observed in e.g. [18]. At first, steel oxidizes and transforms into unstable mobile corrosion products with low mechanical properties. These compounds then diffuse through the stable corrosion layer, re-oxidize in contact with oxygen migrating from the exterior, crystallize and become part of this well-crystallized layer. This latter comprises all the classical corrosion products: magnetite, goethite and lepidocrocite, and is assumed to be located next to the interface of virgin steel and to have a higher elastic modulus. This

Fig. 1. Young modulus of corrosion products for different references.

lamellar arrangement of the corrosion products layers allows describing in a simplified way a progressive evolution of the overall mechanical properties of rusts. This model is integrated into a coupled carbonation–corrosion–creep–damage model, which aims at simulating the response and crack initiation of reinforced concrete structure at the mesoscale. In this approach, the transport of both oxygen and carbon dioxide gases from the exposed surfaces to the steel occurs through the partially saturated connected porosity of the matrix. The kinetics of corrosion takes into account oxygen reduction and oxygen diffusion whose properties are function of water saturation degree [20]. The concrete cracking is reproduced via an isotropic damage variable [21]. To account for the creep phenomenon occurring at mid and long term, a linear viscoelastic model is introduced and coupled to the damage approach [22]. The resulting carbonation – rust formation – mechanical formulation is then implemented into the finite element code Cast3m [23] and applied to the simulation of a 3D concrete numerical sample exhibiting four phases: the matrix mortar, elastic aggregates, corrosion layer and steel rebar. The impact of the parameters introduced in the rust formation model on the time evolution of the corrosion products' thickness is then analyzed, and the time of initial concrete cracking is estimated. 2. Modeling 2.1. Coupled carbonation and drying models We briefly recall in this section the atmospheric carbonation modeling used in our approach. The coupled drying and carbonation of concrete are described by a simplified model, initially proposed in [24] and modified and improved in [25] and [26]. This model assumes that the main phenomena involved in the carbonation are the water transport through the connected porosity, the diffusion of carbon dioxide in the gaseous phase and its dissolution followed by the reaction with portlandite and other hydrated phases to form calcite. The formulation rests on two coupled nonlinear mass balance equations for the water and the carbon dioxide in gaseous phase, respectively, expressed as [26]:
 ∂ðð1−Sr ÞϕP c Þ ¼ ∇ DCO2 ∇ðP c Þ þ W CO2 ∂t

ð1Þ

  ∂ðρl ϕSr Þ ρ ¼ ∇ K ðϕÞ l kr ðSr Þ∇ðP l Þ þ W H2 O : η ∂t

ð2Þ

These equations are governed by the pressure of liquid water Pl and the carbon dioxide pressure P c, respectively. Note that in Eq. (2) it is assumed that the transport process for moisture results only from water pressure gradients, meaning that the diffusion of vapor in gas phase is disregarded. Furthermore, the pressure of the gas phase is neglected with respect to the one of the liquid, such that Pc ≈ − Pl, with Pc the capillary pressure. Another remark is that the CO2 is assumed to behave as an ideal gas, which implies the relation ρc ¼ Mc P c =ðRT Þ, with ρc and Mc the density (kg m−3) and molar mass (kg mol−1) of CO2, respectively, R is the perfect gas constant (8.31 J mol−1 K−1) and T is the temperature assumed constant (T = 293 K). In Eqs. (1)–(2), ϕ is the porosity of concrete, Sr is the saturation degree relative to the total porosity and DCO2 is the effective diffusion coefficient of CO2 (m2 s−1). This last parameter is chosen to evolve as a function of the relative humidity hr via the empirn ical function proposed in [27,28] as DCO2 ðhr Þ ¼ ACO2 ð1−hr Þ CO2 , where ACO2 and nCO2 are adjusting parameters. W CO2 stands for the source term related to the dissolution of CO2 in the liquid phase and precipitation of calcite, ρl is the density of water, K(ϕ) is the intrinsic permeability coefficient depending on the porosity, η is the dynamic viscosity of the water (0.001003 Pa.s), kr(Sr) is the relative permeability to water, and W H2 O is the water source term. The expression of the different

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

functions introduced in Eq. (2), including the desorption isotherm curve specifying the dependence between Sr and Pc (or equivalently Pl), is classically taken from [29]: "

 1 #−m  m 2 pffiffiffiffiffi 1 P l 1−m = ; kr ðSr Þ ¼ Sr 1− 1−Sr m ; P0  3   ϕ 1−ϕ0 2 K ðϕÞ ¼ K 0 1−ϕ ϕ0 

Sr ðP l Þ ¼ 1 þ

ð3Þ

in which P0 and m are parameters to identify, and K0 is the permeability corresponding to the initial porosity ϕ0; K(ϕ) is based on the wellknown Kozeny–Carman relation [30]. The dissolution of carbon dioxide W CO2 is assumed to be proportional to the formation rate of calcite. The precipitation of calcite from portlandite and other main hydrated products is supposed to occur simultaneously but at different rates, portlandite being initially more reactive. However, calcite progressively formed around the portlandite crystals becomes a diffusive obstacle for carbonate and calcium ions, leading to reducing the kinetic constants of the calcite formation reactions (se e.g. [25,26,31]). There exists a direct coupling between carbonation and drying through saturation degree and porosity, which appear in both mass balance equations. Carbonation phenomenon indeed modifies the material microstructure and in particular the matrix porosity, while drying makes the saturation degree evolve. Furthermore, the dissolution of hydrates preceding calcite formation implies a release of water in the porosity, which is taken into account through the source term W H2 O in the mass balance Eq. (2). In the present approach, no diffusion of calcium ion through the water phase contained in the connected porosity is considered, contrary to [24,25]. Instead, the hypothesis that the calcium dissolved from portlandite and the other hydrates is locally transformed into calcite, constituting a local chemical equilibrium, is adopted. For further details, the reader may refer to [26]. 2.2. Steel corrosion modeling 2.2.1. Principles and representation of the corrosion products The simplified formulation describing the steel corrosion processes is inspired by the modeling of cathodic reaction [20]. The process is irreversible and controlled by the concentration of oxygen dissolved at the surface of the virgin steel. Corrosion is assumed to simplify to become active as soon as carbonation process reaches the steel surface. Note that steel potential is not considered in the model, which may lead to some imprecisions regarding the effective corrosion initiation time. As mentioned in the Introduction, the corrosion model is based on the assumption of the existence of two layers of corrosion products having different mechanical properties. In these composite corrosion products, the oxygen reduction occurs at several levels. At first, a portion of initial steel consumes oxygen to form a layer of oxides considered as mobile and unstable rusts with low mechanical properties. Later, this porous rusts in contact with the oxygen contained in the gaseous phase will re-oxidize to become thermodynamically stable rusts, crystallized and compact, with greater mechanical properties. At this stage, this stable layer is no more reactive, and allows oxygen to migrate by diffusion. Several experimental results tend to show that these two layers are in fact difficult to distinguish, as they are alternated or on marbling form [32,33]. The marbling permits an electronical connection between concrete and virgin steel. In our study, we assume that the stable layer, dense with high elastic modulus, is closer to the steel, as observed by e.g. [18,34]. Once this layer is constituted, the newly formed mobile oxides at the steel surface diffuse through the stable layer to reach the unstable one, and are subject to the re-oxidization process. Note that this description of the corrosion products as a two-layer lamellar composite allows its overall properties to vary progressively from the weaker ones at the corrosion initiation, to the greater ones (corresponding to

97

the well-crystallized corrosion products), depending on their respective thickness. We also take into account a TM located between the concrete and the layers of corrosion products. This interphase is supposed to be composed of amorphous iron oxides and hydroxides, and unstable rusts that have migrated and precipitated in the pores of the cement matrix. In this study, the corresponding quantity of iron is evaluated based on experimental profiles of iron in the TM medium [35]. The mechanical and diffusive parameters of this layer are assumed to be the same as the ones of concrete, i.e. the presence of iron precipitates does not affect notably the material. The description of the corroded steel–concrete interface, the different transport mechanisms and the consumption of oxygen is illustrated in Fig. 2. 2.2.2. Oxygen reactive transport in the matrix As explained above, the corrosion rate is assumed to be controlled by the oxygen reduction. The kinetics of formation of corrosion products is estimated via the oxygen consumption, which requires determining the profile of the oxygen into the porous structure. Oxygen is assumed to migrate in the material by diffusion in the gas phase, the following mass balance equation then governs its evolutions:


 ∂ ð1−Sr ÞϕC g ¼ ∇ DO2 ðSr Þ∇C g þ W O2 ∂t

ð4Þ

with Cg the oxygen concentration in the gas phase contained in the pores of the material; W O2 the source term related to the oxygen consumption in the corroded zones at the steel–concrete interface; DO2 the effective diffusion coefficient of oxygen through the porous material, expressed in the same form depending on the relative humidity as the one specifying the diffusion coefficient of CO2 [27]. The two parameters AO2 and nO2 of this expression are selected based on experimental results of [28] for oxygen. The unstable layer is reactive, so it consumes oxygen. The local consumption of oxygen within the conductive oxide layer depends on the specific surface related to the volume developed locally by the porous network. On the other hand, this porous layer is diffusive. The reactive transport of oxygen through this layer is expressed by a simplified one-dimensional mass balance equation (meaning that the diameter of the steel rebar is supposed much greater than the thickness of the corrosion products) involving reaction mechanisms in the form [20]: 2

∂C ∂ C ¼ DR1 2 −kR1 sR1 ϕr C ∂t ∂x

ð5Þ

with C standing for the oxygen concentration dissolved in the liquid phase. The first and second terms of the right-hand side relate to the diffusion and local oxygen consumption, respectively; DR1 and kR1 are the diffusion coefficient of oxygen in the unstable rust layer and the kinetic constant of oxygen reduction; sR1 is the specific surface of the pores relative to the total volume of the porous medium; ϕr is the porosity of the unstable corrosion products layer. This equation is valid for x varying between [0, e1] with e1 the thickness of the unstable layer and the position 0 defining the interface between the TM and the corrosion products (see Fig. 2). The stable layer as defined previously is non-reactive with respect to oxygen. The oxygen transport is then purely diffusive for x ∈ [e1, e1 + e2] with e2 the thickness of the stable layer and DR2 its diffusion coefficient. The mass balance equation then simply reads: 2

∂C ∂ C ¼ DR2 2 : ∂t ∂x

ð6Þ

As in [20], we specialize in this study the two previous mass conservation Eqs. (5) and (6) to the case of steady state. This simplification is justified by the fact that when corrosion becomes active (after about

98

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

Fig. 2. Corrosion layer evolutions. Corrosion at time tn: reduction and diffusion of oxygen in the porous network (top); thickness of corrosion products at time tn + 1 (bottom).

several tens of years, see the simulation Section 3), both saturation degree and oxygen profiles in the concrete material exhibit only very slow variations since the drying processes related to external conditions are also in a quasi steady state. Consequently, the oxygen concentration and moisture content in the corrosion layers also evolve extremely slowly, so that the local terms ∂C/∂t in Eqs. (5) and (6) appear negligible. With this hypothesis, the equations are rewritten as follows: 2

DR1

∂ C −kR1 sR1 ϕr C ¼ 0 ∂x2

ð7Þ

      DR2 1 e1 e2 e1 e1 þ þ cosh sinh sinh k λR1 λR1 λR1 λR1 λR1       γ¼− DR2 e1 e1 e1 þ e2 cosh þ λR1 sinh cosh k λR1 λR1 λR1

In the stable layer defined by x ∈ [e1, e1 + e2], the oxygen concentration solution of Eq. (8) is linear. The oxygen concentration C2 at the surface of the steel (at x = e1 + e2) can then be expressed as a function of C0 and with the help of the continuity condition at x = e1 given by Eq. (9) as follows: 

2

∂ C DR2 2 ¼ 0: ∂x

ð8Þ

The mass balance Eq. (7) related to the reactive unstable layer is a second order differential equation. Considering the boundary conditions at the steel surface (x = e1 + e2 in Fig. 2) specifying that the oxygen does not penetrate into the steel (i.e. the oxygen flux is zero) and the conditions of continuity at the interface of the two layers of corrosion products (x = e1), the solution in the unstable layer for x ∈ [0, e1] takes the form:  C ðxÞ ¼ C 0 cosh

x λR1



  x þ γλR1 C 0 sinh λR1

with C0 the oxygen concentration at the unstable layer/TM interface at x = 0, which is obtained by the FE resolution of the whole equations set on the considered mesostructure; λR1 is a parameter defined by Eq. (10) and γ a coefficient determined from the boundary conditions and expressed by Eq. (11); k is the kinetic constant of oxygen reduction at the steel surface.

λR1

sffiffiffiffiffiffiffiffiffiffiffiffiffi DR1 ¼ kR1 sR1

C2 ¼

        e2 e1 e1 e1 e1 þ e2 γ cosh þ cosh þ γλR1 sinh C0 : sinh λR1 λR1 λR1 λR1 λR1

ð12Þ Finally, at the surface of the steel rebar, the consumption of oxygen is a pure interfacial reaction step. The corresponding reduction of oxygen W ST 02 can be expressed as a linear relation with respect to the local concentration C2, involving the kinetic constant k and the fraction of porosity filled with water as [20]: W ST 02 ¼ kSr ϕr C 2 ¼

ð9Þ

ð10Þ

ð11Þ

RT kSr ϕr C g2 KH

ð13Þ

where Cg2 is the concentration of oxygen in the gaseous phase at x = e1 + e2 and KH is Henry dissolution constant of gaseous oxygen (7.47 × 104 Pa m3mol−1). The quantity of oxygen reduced per second in the unstable layer W R1 02 of rust is obtained by integrating the source term appearing in Eq. (7) along the layer thickness as:

W R1 02 ¼

RT kR1 sR1 ϕr KH

Z

e1 0

C g ðxÞ:

ð14Þ

The sum of the oxygen consumption in the unstable rust and at steel surface defines the source term W O2 appearing in the mass balance

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

Eq. (4). Solving this equation then allows us to get the oxygen profile in the structure. W O2 takes the form: R1 W O2 ¼ W ST 02 þ W 02 :

ð15Þ

2.2.3. Effects of the saturation degree on the oxygen reduction To illustrate the influence of the saturation degree on the oxygen profile within the corrosion products layers (excluding the TM), 1-D simulations are carried out with the following conditions. The thickness of the two layers is supposed equal to 10 μm each; the oxygen diffusion coefficients are taken as DR1 = DR2 = 1 × 10−9 m2 s−1; the other parameters are set to [20] sR1 = 3.7 × 107 m−1, ϕr = 0.1 [36], kR1 = 1 × 10−9 m s−1, k = 1 × 10−5 m s−1 [20]. Fig. 3 (up) shows the profiles of normalized oxygen concentration within both stable and unstable corrosion layers, for the three following different saturation degrees: 0.25, 0.5 and 1. We observe as expected that these profiles are linear in the stable corrosion layer (i.e. for x ∈ [e1, e1 + e2]), while they evolve nonlinearly in the unstable one. We also note that with the considered set of parameters, and given the relatively narrow range of saturation degree variations in practice (between about 0.6 to 0.9), the influence of Sr on the oxygen profiles appears quite moderate. Fig. 3 (down) presents the evolutions of the source term W O2 defined in Eq. (15) as a function of the oxygen diffusion coefficients in both layers DRi = DR1 = DR2, and for three different saturation degrees. We remark that for DRi ranging between about 1 × 10−10 and 1 × 10−8 m2 s−1, W O2 evolves rapidly whatever the saturation degree, while its variations are more progressive for both greater and lower values of DRi. As expected, the source term is greater for higher values of diffusion coefficients DRi and higher saturation degree. The influence of this latter parameter appears slightly more pronounced than on Fig. 3 (up), which highlights its importance in the corrosion processes as modeled in this study.

99

2.2.4. Thickness of corrosion products Having specified the source term W O2 , we now detail the methodology for calculating the thickness of the corrosion products layers, which will be used to determine the induced mechanical loading. As explained above, the kinetics of formation of the stable layer directly relates to the kinetics of transformation of unstable rust into stable one, and depends on the oxygen concentration in the unstable layer. The corrosion rate can be expressed by means of the oxygen reduction rate. During a time interval Δt, the additional thickness of the stable layer ΔeR2 is given by the following relation [20,36]: ΔeR2 ¼

4 M Fe R1 W Δt 3 ρ Fe O2

ð16Þ

in which MFe and ρFe are respectively the molar mass of the iron and steel density (55.8 g/mol and 7800 kg/m3). The thickness of the unstable products layer evolves because of three competing mechanisms: it increases due to the oxidation of virgin steel, but at the same time decreases due to its transformation into stable products and by precipitation in the TM. In our calculations, this precipitation is reformulated by defining a loss of corrosion products representing the quantity of iron oxides that dissolve and migrate in the TM. The variation of the unstable layer thickness in a time interval Δt is expressed as follows, with αMT defining the proportion of unstable products migrating in the TM: ΔeR1 ¼

4 M Fe AC 4 M Fe R1 W ð1−α MT ÞΔt− W Δt: 3 ρ Fe O2 3 ρ Fe O2

ð17Þ

It could be noted that in this model, the corrosion rate is expressed through the macroscopic thickness of the rust layers, which directly takes into account the expansion at the interface with concrete. The influence of the parameters introduced in Eqs. (16)–(17) on the respective values of ΔeR1 and ΔeR2 will be analyzed in Section 3.3. 2.3. Coupled creep-damage model

0.25

/

0.5 1.0

( )

mol m

s

1.0

Corrosion due to carbonation of concrete subjected to the atmosphere occurs after several decades of service. In this context, it is important to take into account the time-dependency of the mechanical behavior of the material, which manifests itself by a time evolution of strains and stresses at fixed loadings. When subjected to the action of a moderate long-term loading, concrete behaves like a viscous material. This creep phenomenon is accounting for in the model by adopting a standard linear viscoelastic formulation for the matrix, while the aggregate phase remains elastic. On the other hand, cracking appearing in the matrix is reproduced via a classical isotropic damage model, as is explained below. The viscodamage approach developed in [22] is adopted in this study. It consists in a creep formulation based on two generalized Maxwell models (characterizing the time evolutions of hydrostatic and deviatoric strains, respectively) with N + 1 elements, coupled to the

0.5 0.25

(m m s

)

Fig. 3. Profiles of oxygen concentration within the two corrosion layers (up), and values of the source term W O2 as a function of the oxygen diffusion coefficient (down) for 3 different saturation degrees (0.25, 0.5 and 1).

Fig. 4. Generalized Maxwell model.

100

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

damage model due to Mazars [21] via the concept of pseudo-strains introduced in [37,38]. The generalized Maxwell models are composed of N parallel chains, consisting in a combination of one spring and one dashpot in series (Fig. 4); one additional chain in parallel comprises only one spring, which ensures a finite limit to the strain. Each branch of this model is characterized by two parameters: an elastic modulus Ei and a relaxation time τi linked by: τi ¼

ηi Ei

ð18Þ

with ηi the viscosity of the dashpot i. The time-dependent stress–strain relation defining the behavior of the matrix is classically expressed by means of the following nonaging linear viscoelastic formulation: Zt

m

k ðt−τ Þ

σðt Þ ¼ 3

dϵ m dτ1 þ 2 dτ

0

Zt μ m ðt−τÞ

dem dτ: dτ

ð19Þ

0

with σ(t), ϵm and em the macroscopic stress tensor, bulk strain and deviatoric strain tensor, respectively; 1 is the second order identity tensor. The matrix moduli km(t) and μm(t) take the following form according to the generalized Maxwell model: m

m

k ðt Þ ¼ k∞ þ

N X

ki e−t=τi ; μ m ðt Þ ¼ μ m ∞ þ m

i¼1

m

N X

−t=τ i μm i e

m

ð20Þ

i¼1

m m where km i , μi and τi are the bulk and shear moduli and relaxation time m th of the i chain of the generalized Maxwell models; km ∞ and μ ∞ are the moduli of the chains without dashpot. As mentioned above, the behavior law (20) is coupled to the isotropic damage Mazars model through the simplified pseudo-strain approach due to [37,38]. We simply recall hereafter the principles of this method in a 1D formulation; the readers are invited to refer to [37,38] and [22,39] for further details. In the onedimensional case, Eq. (19) simplifies into:

Zt σ¼

dε Eðt−τÞ dτ dτ

ð21Þ

0

where σ and ε are the stress and strain, E(t) is the relaxation function. This viscoelastic problem could then be converted into an elastic

equivalent formulation by using the pseudo-strains εR introduced in [37,38,40] and defined by:

εR ¼

σ 1 ¼ E R ER

Zt Eðt−τÞ

dε dτ dτ

ð22Þ

0

with ER a model parameter called reference modulus. These pseudostrains encompass the entire history of loading. Like in [40], the key idea of [22,39] is to introduce the effects of the damage variable in the equivalent elastic formulation involving the pseudo-strains; we then obtain: σ ¼ ð1−DÞER εR

ð23Þ

with D the damage variable. If we now go back to the 3D case, the damage evolutions are assumed to be governed by local extensions through an equivalent strain εeq. Based on [21], this strain is calculated from the positive eigenvalues of the pseudo-strain tensor (instead of the total strain tensor as in the original model) in the form [22,39]: εeq ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi hε R iþ : hεR iþ

ð24Þ

with the operator 〈 ⋅ 〉+ defining the Macauley brackets as  hε i iþ ¼ εi if ε i N0 . hε i iþ ¼ 0 if ε i ≤0 The asymmetry in the tensile and compressive stress–strain response of the material behavior is reproduced via two variables Dt (tension) and Dc (compression), whose evolutions are expressed as [21]: Dt ¼ 1−



ð1−At Þε d0 −At exp −Bt εeq −ε d0 εeq

ð25Þ

Dc ¼ 1−



ð1−Ac Þεd0 −Ac exp −Bc εeq −ε d0 εeq

ð26Þ

where At, Bt, Ac and Bc are model parameters and εd0 is the damage threshold. Note that in this first attempt, the effects of the drying shrinkage are not yet taken into account in the model. Moreover, the influence of the saturation degree on the creep characteristics, which are known to be of importance (see e.g. [41]), are not reproduced. These aspects, together

Exposed surfaces

Fig. 5. 3D mesostructure: steel rebar (red) and inclusions (gray), left; mesh of the mesostructure, right.

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107 Table 1 Material parameters and initial conditions for mass conservation Eqs. (25)–(28). Parameter

P0 (GPa)

m

K0 (m2)

ϕ0

ACO2 (m2 s−1)

nCO2

Values Parameter Values

45 AO2 (m2 s−1) 1.2 × 10−8

0.6 nO2 0.23

2.1 × 10−22 hr0 0.92

0.14 Cg0(m2 s−1) 0

0.8 × 10−8 P c 0 (Pa) 0

0.23

with the hysteretic nature of the sorption–desorption isotherm curve, will be the focus of a future work. 3. Numerical applications 3.1. Description of the reinforced concrete mesostructure 3.1.1. Generation of the 3D mesostructure The coupled carbonation–corrosion–damage model presented above is applied to a tridimensional mesostructure defining a representative concrete sample exposed to atmospheric conditions. This sample is a rectangular parallelepiped, which represents a quarter of a slice of a column-like concrete structure reinforced with one rebar. It is composed of four phases: the mortar matrix, the aggregates, the layers of corrosion products and the steel rebar. Taking into account the ability of computers and the size of the studied structure, sand aggregates are not modeled. Note that the explicit introduction of aggregates in the numerical models allows reproducing to a certain extent a part of the heterogeneous nature of the concrete, and differs from most of the other studies which consider a homogeneous concrete (see e.g. [42–44]). This is believed to be important in particular for estimating the initiation and propagation paths of the cracks, as the behavior of both mortar matrix and aggregates are dissimilar. The geometry and mesh of the mesostructure are generated via a script developed in python language, which makes use of the geometry module and meshing softwares plugged to the integration platform for numerical simulations Salome (www.salome-platform.org) [26,45,46]. Aggregates are modeled by Voronoi particles, whose shape is believed to be more realistic than simpler ones like spheres (see also e.g. [47]). Note that the influence of the shape of the aggregates on the carbonation was studied in [26] in terms of degradation kinetics and evolution of mineral profiles. From all the criteria prescribed by the user as the number or volume fraction of inclusive phase and steel, minimum distance between particles, maximum size of inclusions, the script creates a geometry by randomly dispersing the particles in a box of given size. In our study, only five inclusion characteristic dimensions (i.e., equal to the diameter of the sphere of equivalent volume) are selected: 8, 10, 12, 14 and 16 mm. The steel rebar is located at 2 cm from the exposed faces. A non-

1 year

101

overlapping condition between granular particles and exposed surfaces of the box is imposed, whereas the other surfaces are ‘cut’ at the required dimensions to reproduce symmetry conditions. Fig. 4 illustrates the various phases and the mesh of the mesostructure whose dimensions are chosen as 10 cm × 10 cm × 8 cm, the last dimension corresponding to the thickness. The volume fraction of the 440 aggregates in the concrete is 31%. The mesh contains about 300,000 and 105,000 tetrahedral elements in the matrix and aggregate phase, respectively. Although the dimensions and number of aggregates of the sample are necessarily limited to be numerically tractable, we consider that it will nevertheless permit an insightful analysis at the mesoscale of the damage induced by corrosion in carbonated concrete, and of the corresponding effects of the aggregates. As a perspective, a more accurate analysis of the effects of the aggregate packing around the rebar by e.g. specifically describing this zone with finer aggregates would deserve to be carried out. Moreover, the impact of the mesh size has not been investigated here as the computations are yet quite long with the moderate size mesh (about 3–5 days), and to limit the length of the paper. We have assumed that the mesh used for the simulations is sufficiently fine to provide a correct description of the different variables of the problem and of the crack pattern. 3.1.2. Numerical representation of the corrosion layers As explained in Section 2.2, the growth of the different layers of corrosion products is estimated directly via their macroscopic thicknesses, without introducing any expansion coefficient. From the mechanical viewpoint, this progressive growth generates strains and induced stresses at the mortar–steel interface. At the mesoscale as considered in this study, the thickness of these layers appears however very small, with sizes of about 2 orders of magnitude lower than those of aggregates and steel rebar. To avoid explicitly meshing in the 3D numerical model the expansive corrosion layers, the approach proposed in [7, 8] is retained. It consists in using specific interface finite elements called ‘joint’ elements, composed of two superposed surfaces, which are placed between the mortar and steel (see also e.g. [43]). These elements are characterized by a normal kn and a tangential (shear) ks stiffness directly depending on their (virtual) thickness as [8]: kn ¼

Er Er ; ks ¼ er 2ð1 þ υr Þer

ð27Þ

in which Er , υr and er are respectively the elastic modulus, Poisson ratio and thickness of the interface element. We now specialize these simplified relations to the case of the layered corrosion products. According to the scenario of the existence of two layers of corrosion having different elastic modulus E1 (unstable layer) and E2 (stable layer), the equivalent

20 years Fig. 6. Portlandite volume fraction in the mesostructure at different times.

40 years

102

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

caused by both neighboring mortar and steel rebar, inducing radial compressive states of stresses. Obviously, these local stresses are superposed with the ones resulting from potential external mechanical loadings, or from particular boundary conditions. It should be noted that an initial transition zone similar to the ITZ may be introduced if needed in the above model to replicate a region of mortar surrounding the rebar with reduced mechanical properties and specific characteristics relating to corrosion products growth. Such zone is not accounted for in this study. Finally, steel is assumed to behave as an elastic linear material with young modulus Es = 210 × 106 Pa and Poisson ratio υs = 0.3, and the reduction of rebar diameter due to corrosion is neglected.

Fig. 7. Evolution of carbonation depth.

elastic modulus Er of these two layers may be expressed by the following simple series model: er Er ¼ e e2 : 1 þ E1 E2

ð28Þ

The total thickness er is given by er = e1 + e2. Note that since the Poisson ratio υi is assumed equal to 0.37 for both layers, Eqs. (27)–(28) hold for the composite rusts with υr = 0.37. Taking into account the experimental results on the elastic modulus of the corrosion products mentioned in the Introduction, the very small value E1 = 1 GPa is associated to the unstable rust layer, while a relatively high value of E2 = 150 GPa is adopted for the stable one. The impact of the respective values of these elastic moduli on the mechanical response and the crack initiation would deserve to be investigated in details; this aspect is not treated in this paper and is left to a further analysis. Following the numerical method developed in [7,8], the progressive expansion of the layers of corrosion products is reproduced by making use of a thermal analogy and applying a thermal loading to the layers such that: ð29Þ

αΔT ¼ er

where α and ΔT are the coefficient of thermal expansion and the increment of temperature, respectively. Note that in practice only the product αΔT has an importance; in our applications we have taken α = 1 and then imposed ΔT = er. This approach allows introducing an effective expansion of the corrosion products (i.e. in zero-stress conditions), keeping in mind that the real growth is likely to be smaller due to the weak mechanical properties of the unstable layer and to the confinement

1 month

1 year

3.1.3. Loading and external conditions The mesostructure is exposed to ambient atmospheric conditions on two surfaces: east and south, see Fig. 5. On the other surfaces, zero-flux (symmetry) conditions are prescribed for the variables of the mass conservation equations, i.e. oxygen, carbon dioxide and liquid water. It is worth noting that inclusions and steel are considered non-porous, and therefore the phenomena of drying, carbonation and oxygen diffusion occur only in the mortar matrix. Regarding the boundary conditions for the mechanical model, the exposed surfaces are free of stresses, while the perpendicular displacements and tangential stresses of the opposite surfaces are imposed at zero. A compressive stress of 10 MPa is applied to the faces perpendicular to the rebar, to reproduce the loading undergone by a column. This choice seems coherent with the aforementioned focus on the effects of corrosion products expansion on the initiation and propagation of cracks. Different results would probably be obtained if external loadings resulting from e.g. bending states of stresses were considered. These aspects are beyond the scope of this paper and will be investigated in a further contribution. 3.2. Simulation results of carbonation, drying and oxygen transport in the 3D mesostructure The coupled carbonation–drying model has been implemented in the finite element code Cast3M (www-cast3m.cea.fr). The mass conservation equations are solved successively and iteratively with a fully implicit numerical scheme for the time space. The material properties as well as the model parameters and initial conditions (in terms of initial relative humidity hr0, oxygen concentration Cg0 and partial pressure of CO2 P c0 ) are chosen for CEM I concrete with a water/cement ratio equal to 0.42 (for further information, see [25–28]). They are listed in Table 1. For atmospheric carbonation, we chose the following boundary conditions (i.e., conditions in the atmosphere): • Concentration of oxygen: 8.6 mol m−3. • Relative humidity in the atmosphere: 0.60. • Partial pressure of carbon dioxide: 40 Pa.

2 years

Fig. 8. Saturation degree (no unit) profiles at different times.

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

Fig. 9. Influence of the kinetic coefficient k on the corrosion product thicknesses.

103

We now present some results illustrating the capacities of the carbonation-drying model. Note that for conciseness these results, which are by nature spatially inhomogeneous and time-dependent, are not detailed nor analyzed here, as a complete similar analysis may be found in [26] in terms of local portlandite, calcite, porosity… evolutions as a function of time. The interested readers are then invited to refer to this study for more precisions. Fig. 6 shows the evolution of the portlandite volume fraction (no unit) in the mesostructure at different times. The decrease of this volume fraction versus time and distance from the exposed surfaces appears clearly. The effects of the presence of aggregates are also visible locally, i.e. in their immediate vicinity. It has been shown in [26] that these effects are almost negligible at the macroscale in terms of average carbonation front or average calcite and portlandite concentration, for mesostructures with nearly isotropic aggregates as considered in this study. Fig. 7 presents the evolutions of the mean carbonation depth versus time in the considered mesostructure. This depth corresponds to the average distance perpendicular to the exposed surfaces, and has been calculated from the overall results of portlandite volume fraction obtained from the 3D simulations, as defined in [26]. In the model, the change in pH in the concrete is not directly calculated, and the carbonation depth is estimated from the profile of portlandite as proposed in [26]. We consider that the concrete is degraded by carbonation (corresponding to a pH less than 9, according to phenolphthalein test) if the local concentration of portlandite is lower than 10% of its initial concentration. With this definition, Fig. 7 shows that the degradation of the concrete begins after about 36 years, though it is clear that calcite has progressively precipitated into the material from the beginning of the exposure. The carbonation depth reaches the surface of the steel bar in about 50 years. At this time, we consider that the mechanisms of corrosion become active. The necessary condition of a local carbonation (i.e. pH less than 9) for corrosion to be active makes it possible to reproduce a delay effect corresponding to the incubation period. The specimen is completely carbonated after about 90 years. The calculation results indicate that oxygen diffusion within the mesostructure is much more rapid than the carbonation processes. In fact, a quasi-homogeneous concentration within the mesostructure is obtained after about 1 month of simulation. Moreover, this quasihomogeneous distribution is negligibly affected when corrosion becomes active, i.e. when the sink term W O2 introduced in Eq. (4) is non-zero and governs the corrosion products development. This means that a (quasi) steady state is rapidly reached in the sample. However, this observation should be regarded cautiously since in our case no cyclic variations of ambient relative humidity are considered, which would probably lead to some continuous seasonal time variations of oxygen concentration. According to the model, oxygen transport from the exposed surfaces to the steel rebar does not appear to control the corrosion processes, since corrosion kinetics is much slower. In this study, as only drying takes place due to the assumed constant relative humidity in the atmosphere, the water transport is irreversible. A careful analysis of the inhomogeneous 3D numerical results (not detailed here for the sake of brevity) show that drying also occurs relatively

Fig. 10. Simulated total thickness of corrosion products around the rebar at different times.

104

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

quickly in the matrix when compared to the carbonation kinetics. After 1 month, the degree of saturation is significantly reduced in the region near the exposed surfaces (Fig. 8). More importantly, when the corrosion becomes active (after about 36 years), the simulations indicate that the saturation degree is nearly homogeneous within the mesostructure. It is worth noting that the saturation degree plays a very important role in the kinetics of corrosion, and that relative humidity changes in the atmosphere should be carefully considered to accurately reproduce corrosion processes (see e.g. [48]). This aspect will be investigated thoroughly in another study taking into account the drying–wetting cycles due to evolving climatic conditions, and the hysteresis in the absorption/desorption isotherm curve.

3.3. Evolutions of corrosion layer thickness: parametric study In the proposed corrosion model, the kinetics of formation of corrosion products layers is proportional to the kinetics of oxygen consumption in unstable rust and at the surface of the steel rebar. These mechanisms depend on all variables and parameters involved in the mass balance equation of the oxygen, in particular the saturation degree Sr, the oxygen concentration Cg, the diffusion coefficient of oxygen in the different layers DR1 and DR2, and the kinetics coefficients of oxygen reduction associated with interface and through-layer processes k and kR1. In practice, these parameters are very difficult to estimate or identify, and related experimental studies report variations over several orders of magnitude (e.g., [36]). Thus, it appears that the available data is not sufficiently reliable to determine all the parameters introduced in our deterministic model. In this context, we propose to carry out a numerical study focusing on the effects of some of these parameters on the thicknesses of the corrosion products' layers. In particular, we present the simulation results obtained with different values of the kinetic coefficient of oxygen reduction at the steel surface k, assuming that all other parameters are constant. Fig. 9 shows the thickness evolution versus time of the two layers for 4 values of the coefficient k ranging from 1 × 10−8 to 1 × 10−5 s−1. The parameters kR1, DR1 and DR2 are set to 5 × 10−10 m s−1, 1 × 10−10 m2 s−1 and 1 × 10−10 m2 s−1, respectively, after [20]. The ratio of unstable rust dissolved in the TM is assumed here at 20% of its thickness. It is worth noting that the quantity depicted corresponds to the calculated average thickness because corrosion around the rebar occurs non-homogeneously, as shown in Fig. 10 presenting the distribution of the total thickness at different instants. Indeed, the regions closer to the exposed surfaces are, not surprisingly, more corroded. For instance, at 70 years the differences between the most and least corroded points are as high as 40 μm. Moreover, small effects due to the presence of aggregates between the rebar and the surfaces, which modify locally the oxygen, carbon dioxide and water transport, are observable. As such, the model is able to reproduce a certain level of non-uniformity in the corrosion process, due to the non-uniformity of carbonation around the bar. It is interesting to notice that this inhomogeneous distribution of corrosion products as obtained with the proposed model is in general not taken into account in other simplified approaches, see e.g. [9,49–51]. However, a further degree of non-uniformity could be obtained by e.g. introducing in the model a statistical variability in some of its parameters (e.g. in the diffusion coefficients or initial portlandite volume fraction), reflecting local defects or inhomogeneities in the material. The simulations indicate that corrosion is activated at 47 years. We note that the higher the kinetic coefficient k, the greater the thickness of both stable and unstable layers of corrosion products. This is completely consistent because k represents the rate of oxidation of sound steel into unstable oxide-hydroxides. The reaction at the steel surface controls the quantity of unstable rust formed, which becomes the “feedstock” for the process of re-oxidation into stable rust. As expected, the thickness of the stable layer always increases with time, at a ratio depending on k. However, the same tendency is not observed in the unstable rusts. Indeed with a value of k greater than 1 × 10−7 m s−1, we remark that the thickness of the unstable layer increases up to about 60 years, and then decreases. This is explained

Table 2 Material mechanical parameters [52].

Fig. 11. Influence of the diffusion coefficient DR2 on corrosion products thicknesses.

Parameter

Parameter

At Ac εd0 Bt Bc

km ∞ μm ∞ km 1 μm 1 τm 1 km 2

0.8 1.4 1 × 10−4 17,000 1900

(Pa) (Pa) (Pa) (Pa) (days) (Pa)

Parameter 9

6.27 × 10 3.41 × 109 2.93 × 109 7.77 × 109 2 4.21 × 109

μm 2 (Pa) τm 2 (days) km 3 (Pa) m μ3 (Pa) τm 3 (days)

3.54 × 109 20 6.93 × 109 3.32 × 109 150

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

After 50 years

After 51 years

105

After 53 years

Fig. 12. Damage pattern in the matrix obtained with the viscodamage model at different times.

by the fact that the diffusion barrier effect of the layers becomes more and more significant when they progressively grow. Oxygen has to diffuse through a longer distance to reach the steel surface; correspondingly the formation kinetics of unstable rust decelerates. At 60 years, the rate of formation of unstable rust from virgin steel is lower than the transformation rate of unstable rust into stable one. This leads to the reduction of the thickness of unstable rusts. The same tendency has been noticed for the kinetic constant of oxygen reduction in the unstable layer kR1 (the results are not reported here for conciseness). We compared the total thickness of the corrosion layers (Fig. 9 down) with those observed on some existing concrete structures in France (age between 50 and 80 years). The points labeled BdT, CEA, and MdB are respectively the measured thickness of corrosion products of carbonated samples from Bourse du Travail, water house of CEA and Maison du Brésil structures [4]. We notice that a coefficient k ranging between 1 × 10−5 and 1 × 10−6 seems to give more accurate thickness values. When we fix the diffusion coefficient DRi of oxygen in the layer i of rust and vary that of the other layer, we observe a similar tendency (see Fig. 11 in the case of varying DR2 with k = 1 × 10−5 m s−1). Not surprisingly, a higher diffusion coefficient allows for better diffusion of oxygen in the layer, which promotes the cathodic reaction and thus the

formation of corrosion products. In Fig. 11, we note that the time evolution of the stable layer is almost linear. The thickness of the unstable layer increases, then decreases or becomes nearly constant for higher values of DR2 after around 70 years, as already observed in the above analysis. We also notice that for DR2 values greater than about 10−9 m2 s−1, the influence of a change of DR2 on the formation kinetics of corrosion products seems much lower than for lower values (see Fig. 11 down). Finally, a diffusion coefficient DR2 of about 10−9 m2 s−1 gives results of total thickness in relatively good agreement with experimental ones. It's noted that this value is higher than the one recommended in [36]. The above parametric study shows that the evolution of the thickness of corrosion products is very sensitive to the values of parameters introduced in the corrosion model. A further comparison with experimental results obtained for different configurations and environmental conditions would be necessary to precisely determine their order of magnitude. However, it appears that the measured thicknesses of corrosion are in general relatively dispersed: they usually range from several μm to hundreds μm in different existing structures at the same age, depending on the ambient conditions and the material (composition, microstructure etc). The thicknesses reported in different zones of the

Fig. 13. Damage pattern obtained after 51 years in the case with (left) and without (right) creep.

106

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107

same structure can also vary significantly due to the heterogeneous nature of concrete. Consequently, it is of great importance to have a good knowledge of the environmental conditions and of the composition of concrete for a meaningful comparison. 3.4. Mechanical response and crack pattern As mentioned in Section 3.1, the considered structural element models a slice of a quarter of column having axis parallel to the steel bar. A constant compressive load of 10 MPa is assigned on the surfaces perpendicular to steel axis as permanent load. The mortar matrix is modeled as explained in Section 2.3 by two generalized Maxwell viscoelastic models with 4 chains (including one chain purely elastic) coupled with the damage model of Mazars. The material parameters are given in Table 2 in terms of both bulk and shear moduli m m th km i and μi and relaxation times τi of the i chain of the generalized Maxwell models [52]. The parameters related to damage are also listed in Table 2. From the simulation of the coupled phenomena of drying, carbonation and oxygen transport as presented in the previous section, an evaluation of the corrosion products thickness as a function of time is inferred. The following set of parameters regarding the corrosion model has been selected: kR1 = 1 × 10−8 m s−1, k = 1 × 10−5 m s−1, DR1 = DR2 = 1 × 10−9 m2 s−1. As explained, the expansion of the steel-concrete interface generates a mechanical loading, which may damage the mortar. In our calculations, to simplify and because we are mainly interested in the description of the onset of cracking and on the main crack paths, no method has been implemented to regularize the localization of the damage. Fig. 12 shows the damage obtained in the matrix at 3 different instants (up) and the corresponding facies of cracking (down). As may be inferred from e.g. Fig. 10, it appears that damage does not develop uniformly around and along the reinforcement. To make a clearer visualization, we have cut a slice from the middle of the rebar during the image post-treatment to observe the cracking path on half of the structure. Damage initiates in the lower part and at the middle of the rebar, and then propagates along the rebar and to the exposed surfaces. This could be explained first by the non-uniform development of corrosion around the steel bar. As atmospheric boundary conditions are imposed on the south side and the east side of the structure, rusts form at first on southeast side. Further, crack paths are inevitably influenced by the presence of aggregates in the matrix, since they constitute nondamageable obstacles that have to be circumvented. Finally, the absence of a regularization method for the energy dissipated during the damage process may also affect the crack propagation paths. It is worth noting that this modeling corresponds to a perfect scenario in which mortar matrix has any defect. In practice, the network of micro-imperfections in the concrete promotes a random development of rust around the steel rebar. Furthermore, damage seems to reach more quickly the exposed surfaces (only 3 years after the appearance of the first rust) than what is usually observed in natural corrosion. This raises the question about the existence of an initial transitional zone with a large porosity (similar to the ITZ) between the steel and the concrete, which could be progressively filled at first by the growing corrosion products. This phenomenon would indeed lead to a further delay effect on the disruption of concrete by tensile stresses. Finally, to highlight the influence of the viscoelastic behavior of concrete on the initiation and propagation of the damage, an additional simulation with the same loading and boundary conditions has been carried out by considering the Mazars damage model without coupling with the generalized Maxwell models for creep (i.e. with an elasticdamage formulation). Fig. 13 presents the damage patterns obtained in the 2 cases, i.e. with and without creep, at the same instant (51 years). We observe that damage is less developed in the case of the coupled creep-damage model. This could be explained by the progressive relaxation of stresses over time consecutive to the viscous

properties of the material. Such relaxation leads to a reduction of stresses, in particular the tensile ones caused by the growth of corrosion products. The associated total strains do not necessarily diminish in magnitude, however the proposed formulation of the viscodamage model is such that the equivalent strain defined in Eq. (24) that controls the damage evolutions takes into account intrinsically the elastic strains (which indeed decrease) plus only a fraction of the viscous strains (see [22,52] for further details). This result illustrates the importance of the description of the material mechanical behavior, especially at longterm, for an accurate estimation of the cracking initiation and propagation. As already mentioned, a future improvement of the model will consist in including the shrinkage due to drying, and the corresponding effects of the saturation degree variations on the creep formulation. 4. Conclusions This paper proposes a simplified multi-physical approach for modeling the mechanical impact of steel rebar corrosion on concrete cover exposed to the atmosphere. This approach is based on existing studies that have investigated independently the phenomena of atmospheric carbonation, corrosion and creep-damage mechanical behavior. As such, relatively few new parameters are introduced in the resulting coupled model. Due to the carbonation of concrete in such ambient environmental conditions, the rate of rust formation is estimated by taking into account the following coupled mechanisms: carbonation, drying and oxygen diffusion. The assumption of a bi-layer composition of corrosion products as a stable and unstable rusts having different mechanical and diffusive properties is adopted. A parametric analysis regarding the kinetics of formation of these two layers shows a high sensibility of both kinetic and diffusion parameters introduced in the model. Their accurate identification would require a very large experimental database, with results collected from a great range of configurations, material composition and external conditions. The corrosion model is coupled to a viscodamage approach describing the creep behavior and estimating the appearance of cracking due to the progressive expansion of rust at the steel rebar–mortar matrix interface as a function of time. The model is applied to the simulation of a 3D heterogeneous mesostructure with polyhedral inclusions representing the coarse elastic aggregates, a steel rebar and the mortar matrix. The numerical results obtained in terms of damage initiation time and propagation paths consecutive to the growth of corrosion products highlight the influence of the aggregates on the crack pattern. The creep behavior of the material plays also an important role regarding damage propagation. In the future, it would be advisable to perform a detailed validation of the model in terms of cracking appearance on appropriate experimental data. Moreover, several improvements of the proposed model are planned, with some of them currently in progress. First, the drying shrinkage of the mortar will be included in the description of the long-term behavior, together with the effects of varying saturation degrees on the creep model. Second, as the external relative humidity exhibits significant cyclic evolutions within a day (between night and day) as well as between different seasons, a model taking into account the hysteretic behavior of the isotherm curve will be introduced in our approach. This will allow a better estimation of the saturation degree states and their evolutions in the material, which are known to be very influential on the transport mechanisms and the mechanical behavior. Acknowledgments ANDRA is gratefully acknowledged for its financial support. References [1] C.L. Page, K.W.J. Treadaway, Aspects of the electrochemistry of steel in concrete, Nature 297 (5862) (May 1982) 109–115.

T.T.H. Nguyen et al. / Cement and Concrete Research 74 (2015) 95–107 [2] K. Tuutti, Corrosion of steel in concrete, CBI Forsk/Research41982. [3] K. Bhargava, A.K. Ghosh, Y. Mori, S. Ramanujam, Model for cover cracking due to rebar corrosion in RC structures, Eng. Struct. 28 (8) (2006) 1093–1109 (Juillet). [4] V. L'Hostis, D. Neff, L. Bellot-Gurlet, P. Dillmann, Characterization of long-term corrosion of rebars embedded in concretes sampled on French historical buildings aged from 50 to 80 years, Mater. Corros. 60 (2) (Feb. 2009) 93–98. [5] V. Baroghel-Bouny, T.Q. Nguyen, P. Dangla, Assessment and prediction of RC structure service life by means of durability indicators and physical/chemical models, Cem. Concr. Compos. 31 (8) (Sep. 2009) 522–534. [6] F.J. Molina, C. Alonso, C. Andrade, Cover cracking as a function of rebar corrosion: part 2—numerical model, Mater. Struct. 26 (9) (1993) 532–548. [7] Q.T. Nguyen, S. Caré, A. Millard, Y. Berthaud, Analyse de la fissuration du béton armé en corrosion accélérée, C.R. Mec. 335 (2) (Feb. 2007) 99–104. [8] Q.T. Nguyen, S. Caré, Y. Berthaud, A. Millard, F. Ragueneau, Experimental and numerical behaviour of reinforced mortar plates subjected to accelerated corrosion, Int. J. Numer. Anal. Methods Geomech. 35 (10) (Jul. 2011) 1141–1159. [9] K. Suda, S. Misra, K. Motohashi, Corrosion products of reinforcing bars embedded in concrete, Corros. Sci. 35 (5–8) (1993) 1543–1549. [10] Y. Zhao, H. Dai, H. Ren, W. Jin, Experimental study of the modulus of steel corrosion in a concrete port, Corros. Sci. 56 (Mar. 2012) 17–25. [11] A. Dehoux, F. Bouchelaghem, Y. Berthaud, D. Neff, V. L'Hostis, Micromechanical study of corrosion products layers. Part I: experimental characterization, Corros. Sci. 54 (Jan. 2012) 52–59. [12] G.V. Samsonov, The Oxide Handbook, 7, Plenum, New York, 1973. [13] A. Ouglova, Y. Berthaud, M. Francois, F. Foct, Mechanical properties of an iron oxide formed by corrosion in reinforced concrete structures, Corros. Sci. 48 (12) (Dec. 2006) 3988–4000. [14] S. Caré, A. Raharinaivo, Influence of impressed current on the initiation of damage in reinforced mortar due to corrosion of embedded steel, Cem. Concr. Res. 37 (12) (Dec. 2007) 1598–1612. [15] T. Misawa, K. Hashimoto, S. Shimodaira, The mechanism of formation of iron oxide and oxyhydroxides in aqueous solutions at room temperature, Corros. Sci. 14 (2) (1974) 131–149. [16] J.M. Blengino, M. Keddam, J.P. Labbe, L. Robbiola, Physico-chemical characterization of corrosion layers formed on iron in a sodium carbonate-bicarbonate containing environment, Corros. Sci. 37 (4) (Apr. 1995) 621–643. [17] P. Refait, S.H. Drissi, J. Pytkiewicz, J.-M.R. Génin, The anionic species competition in iron aqueous corrosion: role of various green rust compounds, Corros. Sci. 39 (9) (Sep. 1997) 1699–1710. [18] G. Duffó, W. Morris, I. Raspini, C. Saragovi, A study of steel rebars embedded in concrete during 65 years, Corros. Sci. 46 (9) (2004) 2143–2157. [19] P. Dillmann, F. Mazaudier, S. Hœrlé, Advances in understanding atmospheric corrosion of iron. I. Rust characterisation of ancient ferrous artefacts exposed to indoor atmospheric corrosion, Corros. Sci. 46 (6) (Jun. 2004) 1401–1429. [20] B. Huet, V. L'hostis, G. Santarini, D. Feron, H. Idrissi, Steel corrosion in concrete: determinist modeling of cathodic reaction as a function of water saturation degree, Corros. Sci. 49 (4) (2007) 1918–1932. [21] J. Mazars, A description of micro-and macroscale damage of concrete structures, Eng. Fract. Mech. 25 (5–6) (1986) 729–737. [22] M.-Q. Thai, B. Bary, Q.-C. He, A homogenization-enriched viscodamage model for cement-based material creep, Eng. Fract. Mech. 126 (Aug. 2014) 54–72. [23] Cast3M, www-cast3m.cea.fr.. [24] B. Bary, A. Sellier, Coupled moisture—carbon dioxide–calcium transfer model for carbonation of concrete, Cem. Concr. Res. 34 (10) (2004) 1859–1872. [25] B. Bary, C. Mügler, Simplified modelling and numerical simulations of concrete carbonation in unsaturated conditions, Rev. euro. de génie civil 10 (9) (2006) 1049–1072. [26] T. de Larrard, B. Bary, E. Adam, F. Kloss, Influence of aggregate shapes on drying and carbonation phenomena in 3D concrete numerical samples, Comput. Mater. Sci. 72 (2013) 1–14. [27] V.G. Papadakis, C.G. Vayenas, M.N. Fardis, Physical and chemical characteristics affecting the durability of concrete, ACI Mater. J. 88 (2) (1991). [28] Y.F. Houst, Diffusion de gaz, carbonatation et retrait de la pâte de ciment durcie, Université de Lausanne, 1993. [29] M.T. van Genuchten, Closed-form equation for predicting the hydraulic conductivity of unsaturated soils, Soil Sci. Soc. Am. J. 44 (5) (1980) 892–898.

107

[30] P. Carman, Fluid flow through granular beds, Trans. Inst. Chem. Eng. 15 (1937) 150–166. [31] M. Thiery, G. Villain, P. Dangla, G. Platret, Investigation of the carbonation front shape on cementitious materials: effects of the chemical kinetics, Cem. Concr. Res. 37 (7) (2007) 1047–1058. [32] M.L. Schlegel, C. Bataillon, K. Benhamida, C. Blanc, D. Menut, J.L. Lacour, Metal corrosion and argillite transformation at the water-saturated, high-temperature iron– clay interface: a microscopic-scale study, Appl. Geochem. 23 (9) (2008) 2619–2633. [33] E. Burger, L. Legrand, D. Neff, H. Faiz, S. Perrin, V. L'Hostis, P. Dillmann, In situ structural characterisation of nonstable phases involved in atmospheric corrosion of ferrous heritage artefacts, Corros. Eng. Sci. Technol. 45 (5) (2010) 395–399. [34] A. Demoulin, C. Trigance, D. Neff, E. Foy, P. Dillmann, V. L'Hostis, The evolution of the corrosion of iron in hydraulic binders analysed from 46- and 260-year-old buildings, Corros. Sci. 52 (10) (Oct. 2010) 3168–3179. [35] W.J. Chitty, P. Dillmann, V. L'Hostis, C. Lombard, Long-term corrosion resistance of metallic reinforcements in concrete—a study of corrosion mechanisms based on archaeological artefacts, Corros. Sci. 47 (6) (2005) 1555–1581. [36] W.J. Chitty, P. Dillmann, V. L'Hostis, A. Millard, Long-term corrosion of rebars embedded in aerial and hydraulic binders—parametric study and first step of modelling, Corros. Sci. 50 (11) (2008) 3047–3055. [37] R.A. Schapery, Correspondence principles and a generalized J integral for large deformation and fracture analysis of viscoelastic media, Int. J. Fract. 25 (3) (1984) 195–223. [38] R.A. Schapery, A micromechanical model for non-linear viscoelastic behavior of particle-reinforced rubber with distributed damage, Eng. Fract. Mech. 25 (5–6) (1986) 845–867. [39] M.-Q. Thai, Modélisation micromécanique et simulation numérique du fluage des bétons avec prise en compte de la microfissuration et des effets thermohydriques(PhD thesis (in French)) Université Paris-Est, 2012. [40] S.W. Park, R.Y. Kim, R.A. Schapery, A viscoelastic continuum damage model and its application to uniaxial behavior of asphalt concrete, Mech. Mater. 24 (4) (1996) 241–255. [41] Z.P. Bažant, Y. Xi, Drying creep of concrete: constitutive model and new experiments separating its mechanisms, Mater. Struct. 27 (1) (1994) 3–14. [42] J. Ožbolt, G. Balabanić, M. Kušter, 3D numerical modelling of steel corrosion in concrete structures, Corros. Sci. 53 (12) (2011) 4166–4177 (Décembre). [43] B. Richard, F. Ragueneau, C. Cremona, L. Adelaide, J.L. Tailhan, A three-dimensional steel/concrete interface model including corrosion effects, Eng. Fract. Mech. 77 (6) (2010) 951–973 (Avril). [44] P.J. Sánchez, A.E. Huespe, J. Oliver, S. Toro, Mesoscopic model to simulate the mechanical behavior of reinforced concrete members affected by corrosion, Int. J. Solids Struct. 47 (5) (Mar. 2010) 559–570. [45] C. Bourcier, W. Dridi, L. Chomat, E. Laucoin, B. Bary, E. Adam, Combs: open source python library for RVE generation. Application to microscale diffusion simulations in cementitious materialspresented at the Joint International Conference on Supercomputing in Nuclear Applications and Monte Carlo 2013 (SNA + MC 2013), October 27–31, Paris, 2013. [46] B. Bary, M. Ben Haha, E. Adam, P. Montarnal, Numerical and analytical effective elastic properties of degraded cement pastes, Cem. Concr. Res. 39 (10) (Oct. 2009) 902–912. [47] I. Comby-Peyrot, F. Bernard, P.O. Bouchard, F. Bay, E. Garcia-Diaz, Development and validation of a 3D computational tool to describe concrete behaviour at mesoscale. Application to the alkali–silica reaction, Comput. Mater. Sci. 46 (4) (2009) 1163–1177. [48] C. Andrade, C. Alonso, J. Sarría, Corrosion rate evolution in concrete structures exposed to the atmosphere, Cem. Concr. Compos. 24 (1) (2002) 55–64 (Février). [49] T. El Maaddawy, K. Soudki, A model for prediction of time from corrosion initiation to corrosion cracking, Cem. Concr. Compos. 29 (3) (Mar. 2007) 168–175. [50] H.M. Shodja, K. Kiani, A. Hashemian, A model for the evolution of concrete deterioration due to reinforcement corrosion, Math. Comput. Model. 52 (9–10) (Nov. 2010) 1403–1422. [51] Y. Zhao, J. Yu, W. Jin, Damage analysis and cracking model of reinforced concrete structures with rebar corrosion, Corros. Sci. 53 (10) (2011) 3388–3397. [52] B. Bary, Q.C. He, M.Q. Thai, Coupled damage and multiscale creep model applied to cementitious materialspresented at the CONCREEP 09, Sept. 22–24, Boston, USA 2013, pp. 219–226.