Hysteretic moisture behavior of concrete: Modeling and analysis

Cement and Concrete Research 42 (2012) 1379–1388

Contents lists available at SciVerse ScienceDirect

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

Hysteretic moisture behavior of concrete: Modeling and analysis Hannelore Derluyn a, b,⁎, Dominique Derome b, Jan Carmeliet a, b, Eric Stora c, Rémi Barbarulo c a b c

ETH Zürich, Chair of Building Physics, Wolfgang-Pauli-Strasse 15, 8093 Zürich, Switzerland EMPA, Swiss Federal Laboratories for Materials Science and Technology, Laboratory for Building Science and Technology, Überlandstrasse 129, 8600 Dübendorf, Switzerland Lafarge Research Center, 95, rue de Montmurier, BP 15, 38291 St Quentin Fallavier Cedex, France

a r t i c l e

i n f o

Article history: Received 11 February 2012 Accepted 26 June 2012 Keywords: C. Durability C. Finite element analysis C. Transport properties E. Modeling Moisture hysteresis

a b s t r a c t Durability of reinforced concrete is crucial for civil and building infrastructures. Deterioration related to corrosion of the reinforcement can be caused by chloride diffusion or carbonation due to atmospheric CO2 ingress. Chlorides and CO2 transport coefficients are dependent on the saturation degree of the porous material. Thus, the durability of concrete elements submitted to cyclic atmospheric conditions could depend on the hysteretic moisture behavior of concrete. To assess this influence, an analytical model has been developed describing hysteretic moisture behavior based on the independent domain theory. This model is incorporated in a finite element code for coupled heat and mass transfer, allowing to calculate the concrete response to changing temperature and relative humidity conditions. A comparison is made between simulations using the hysteresis model and using either the main ad- or desorption curve. The results show that durability risks may be underestimated when omitting moisture hysteresis in the outer concrete layer. © 2012 Elsevier Ltd. All rights reserved.

1. Introduction The durability of reinforced concrete is crucial for civil and building infrastructures. The high alkalinity of the cement matrix promotes the formation of a passive layer on steel, protecting it against corrosion. Destabilization of this layer can happen because of the presence of chlorides, or because of a pH drop due to atmospheric carbonation [1,2]. Steel corrosion may then initiate and lead to concrete degradation. Chloride ingress occurs by diffusion through the liquid phase, or by advective transport induced by water movement in the porous medium. CO2 ingress happens by gas diffusion in the non saturated part of the porous medium. In both cases, the transport of chlorides or CO2 depends on the saturation degree of the porous material. Concrete exhibits a moisture hysteretic behavior that could thus affect the ingress of chlorides and CO2, with consequences for the durability of reinforced concrete. Moisture hysteresis refers to the phenomenon that, at the same relative humidity, the material experiences a different degree of moisture saturation or level of moisture content depending on its loading history. The difference in moisture content at the same relative humidity for cement or concrete during repeated cycles of ad- and desorption under isothermal conditions is illustrated in amongst others De Belie et al. [3], Anderberg and Wadsö [4], Baroghel-Bouny [5] and Maruyama and ⁎ Corresponding author at: ETH Zürich, Chair of Building Physics, Wolfgang-Pauli-Strasse 15, 8093 Zürich, Switzerland. Tel.: + 41 58 765 4722, + 41 44 633 3632; fax: + 41 44 633 1041. E-mail addresses: [email protected], [email protected] (H. Derluyn), [email protected] (D. Derome), [email protected], [email protected] (J. Carmeliet), [email protected] (E. Stora), [email protected] (R. Barbarulo). 0008-8846/$ – see front matter © 2012 Elsevier Ltd. All rights reserved. doi:10.1016/j.cemconres.2012.06.010

Igarashi [6]. These hysteresis measurements count a limited number of ad/desorption loops, as the sorption process of concrete is slow. However, when a concrete structure is exposed to environmental conditions, it undergoes a large number of sorption cycles. Thus the ability to predict the hysteretic response of a concrete structure would be valuable especially with respect to its durability. The assessment of the saturation degree in a concrete structure exposed to real climatic conditions requires a numerical solving of the transport equations. In this paper, we describe a methodology to incorporate the hysteretic moisture behavior in a finite element model for heat and mass transfer. In contrast to other moisture hysteresis models mentioned in literature based on pore network modeling and the ink-bottle effect [7–10] or multi-layer adsorption modeling [8], we describe hysteretic moisture behavior using the independent domain theory [11–15]. The proposed methodology could be useful for concrete practitioners to determine the impact of simulating only with a main adsorption or a main desorption curve or incorporating the hysteretic moisture behavior. Advanced durability models [2,16–18] are being developed to extrapolate accelerated lab tests to the long-term behavior of concrete in real life conditions. These approaches give encouraging results in well-defined controlled conditions but the application of these models to situations where many relative humidity cycles are involved is a challenge. An insight into whether or not taking into account moisture hysteresis results in the over- or underestimation of degradation risks would be interesting. As moisture hysteresis can be of interest to other applications than concrete, the model is developed in a general way so that it can also be applied for other porous materials. In the first part of the paper, the analytical moisture hysteresis model, based on the independent domain theory, is explained. In

1380

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

the second part, the implementation of this model in a finite element framework for heat and mass transfer is described. The third part summarizes the heat and moisture concrete properties used for the simulation examples, which are discussed in the fourth part. 2. Moisture hysteresis model 2.1. Independent domain theory and IPM model Following the independent domain approach [11–15], a porous material contains a number of domains or sorption sites, which behave independently. Fig. 1 describes the typical behavior of a sorption site. With increasing relative humidity, an initially dry sorption site becomes filled at a critical relative humidity ϕa resulting in a moisture content mf. When the relative humidity decreases again, the site unfills at a lower relative humidity ϕd, resulting in a jump back to the unfilled state. Since ϕd ≤ ϕa hysteresis occurs between adsorption and desorption, ϕa and ϕd represent respectively the relative humidity in the adsorption (a) or desorption (d) regime. Adsorption includes both multilayer adsorption of water molecules on the pore wall surfaces and capillary condensation. Current physical explanations for the existence of hysteresis for a single sorption site are, for example, capillary condensation hysteresis, contact angle hysteresis and the ink-bottle effect due to entry pores with small diameter. A porous material consists of a number of sorption sites. Thus changes in moisture content are the outcome of the behavior of an assemblage of independent sorption sites. An efficient way to represent sorption sites is to map a site characterized by (ϕa, ϕd), on the half ϕa − ϕd space, also called the PM space (Preisach–Mayergoyz [15,19]). The PM space is represented in Fig. 2. On the PM space, a frequency density distribution ρ(ϕa, ϕd), called the PM distribution, represents the number of sorption sites. The PM space is triangular and filling occurs in its lower surface, since ϕd ≤ϕa. Suppose the specimen is taken through a relative humidity (RH) history RH=0, RHA, RHB, RHC (=RHA), RHD ≈1, RHE (=RHA), RHF (=RHB), 0 (Fig. 3a). The series of Fig. 3b show the states of the sorption elements in the PM space for each RH (filled elements are gray). The resulting sorption curve is shown in Fig. 3c. Initially the porous material is at RH=0, all sorption sites are unfilled and the moisture content is 0. When RH is increased towards RHB, the sorption sites become progressively filled from left to right in the PM space, and the main adsorption curve is followed. As RH decreases from RHB to RHC the elements empty from top to bottom and a primary desorption curve is formed. We observe a different moisture content at RHA and RHC (=RHA) resulting from the difference in filling (A→B) and unfilling of elements (B→C). Comparing the gray areas representing the filling in the PM space for state A and C, we see that the

Fig. 2. The triangular PM space on which a density distribution function ρ(ϕa, ϕd) is defined. The gray triangle represents the area over which integration is needed to obtain the function value H(RHF, RHE) in the integrated PM space (see Fig. 4).

material is at a different filling state. When RH is further increased from RHC, a loop is formed, which closes at RHB. From then on, the main adsorption isotherm is again followed towards RHD. When decreasing RH, the main desorption isotherm is followed, emptying elements from top to bottom. While increasing and decreasing again RH, a second loop is formed (E→F→E), which is congruent with the first loop (B→C→B) formed during the adsorption. The region of filled sites is denoted Ω. After a certain RH loading history, Ω is typically bordered at the right side by a staircase line. Performing integration over the Ω domain gives the moisture content M: MðϕÞ ¼ ∬ ρðϕa ; ϕd Þmðϕa ; ϕd Þdϕa dϕd :

ð1Þ

Ω

To determine the integrals, the PM density ρ(ϕa, ϕd) and the moisture function m(ϕa, ϕd) have to be known. Instead of trying to directly determine these functions, the integrated PM (IPM) approach is adopted. In the IPM approach, an integral function H is defined as: ϕa ϕa

H ðϕa ; ϕd Þ ¼ ∫ ∫ ρðx; yÞmðx; yÞd xd y

ð2Þ

ϕd y

where for example H(RHE, RHF) equals the integral over the triangle (RHE, RHE), (RHF, RHE), (RHF, RHF) (shown in Fig. 2). By the introduction of the IPM function H, the integral M in Eq. (1) can be calculated as: M ðϕÞ ¼ H ðϕa1 ; ϕd1 Þ−Hðϕa1 ; ϕd2 Þ þ H ðϕa2 ; ϕd2 Þ−H ðϕa2 ; ϕd3 Þ þ ::::: m n X X   ¼ H ðϕai ; ϕdi Þ− H ϕai ; ϕdiþ1 i¼1

i¼1

ð3Þ with m and n the number of vertical respectively horizontal boundary segments characterizing the staircase border. Thus the moisture content M can be determined by summing and subtracting values of H. The IPM function H covers the triangle of the PM space, attains its maximum at the bottom and right borders, decreases towards the diagonal and becomes zero on the diagonal (Fig. 4). 2.2. Describing the H-function in the IPM space

Fig. 1. Characteristic behavior of a hysteretic sorption site. During increasing relative humidity a jump in moisture content occurs to mf at a relative humidity ϕa. During decreasing relative humidity, unfilling occurs at ϕd, jumping back from mf to the dry state.

The IPM function H(ϕa, ϕd) is determined based on experimental data of the moisture content M (Eq. (3)). An analytical expression is proposed where the parameter values are found by a best fitting procedure on the experimental data. First we describe the main boundary curves in the

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

1381

Fig. 3. Typical snapshots during (a) a relative humidity loading history of (b) the PM space and (c) the moisture content M. ϕa is the relative humidity in adsorption, ϕd is the relative humidity in desorption and Ω is the region of domains completely filled with water.

IPM space. The main adsorption isotherm is given by Ma(ϕa)=H(ϕa,0), the boundary curve at the bottom side in the IPM space. The main desorption isotherm Md(ϕd)is located at ϕa =1 and is given by: Md ðϕd Þ ¼ H ð1; 0Þ−H ð1; ϕd Þ:

ð4Þ

For ϕd = 1, the desorption curve attains its maximum or Md(1) = Mmax. Applying Eq. (4) and realizing that H(1,1) = 0, we find that H(1,0) = Mmax. Eq. (4) can then be rewritten as: H ð1; ϕd Þ ¼ Mmax −Md ðϕd Þ

ð5Þ

with H(1,ϕd) the boundary curve at the right side in the IPM space, which is defined by the main desorption isotherm. In order to describe the S-shaped sorption curves typically observed for materials like concrete or stone, the main desorption isotherm Md(ϕd)is described by the summation of a concave and a convex curve, or Md = M1d + M2d. We propose for M1d and M2d respectively: ð1 þ ad Þ⋅ϕd M1d ðϕd Þ ¼ M 1max ⋅ 1 þ ad ⋅ϕd

ð6Þ

M2d ðϕd Þ ¼ M 2max ⋅

1− expðbd ⋅ϕd Þ 1− expðbd Þ

ð7Þ

where ad and bd are both positive values and M1max and M2max are defined as: M1max ¼ l⋅Mmax M2max ¼ ð1−lÞ⋅Mmax

ð8Þ

with 0 ≤ l ≤ 1 so that M1max + M2max = Mmax, the maximal moisture content. According to Eq. (5), H(1,ϕd) can then also be expressed as a sum of two functions: H ð1; ϕd Þ ¼ H 1 ð1; ϕd Þ þ H 2 ð1; ϕd Þ

ð9Þ

with: H 1 ð1; ϕd Þ ¼ M1max −M1d ðϕd Þ

ð10Þ

H 2 ð1; ϕd Þ ¼ M2max −M2d ðϕd Þ

ð11Þ

Eq. (9) includes three parameters: l, ad and bd. To describe the decay from the boundary curve H(1,ϕd) to the diagonal in a direction parallel to the adsorption (ϕa) axis and so to describe the full IPM space, we multiply H1(1,ϕd) and H2(1,ϕd) with the functions f1(y) and f2(y) respectively: H 1 ðϕa ; ϕd Þ ¼ ðM1max −M1d ðϕd ÞÞ⋅f 1 ðyÞ

ð12Þ

H 2 ðϕa ; ϕd Þ ¼ ðM2max −M2d ðϕd ÞÞ⋅f 2 ðyÞ

ð13Þ

The variable y is defined as: yðϕa ; ϕd Þ ¼

Fig. 4. Representation of the function H(ϕa, ϕd) in the IPM space. The IPM function is zero on the diagonal and attains maximal values at the boundaries.

ϕa −ϕd 1−ϕd

ð14Þ

and gives the normalized distance from the diagonal in horizontal direction. Special cases in the IPM space are: y = 0 describing the location of the diagonal (ϕa = ϕd), y = ϕa representing the location of main adsorption (ϕd = 0), and y = 1 describing the location of main desorption (ϕa = 1). The function H(ϕa, ϕd) is by definition zero on the diagonal, or fi(0) = 0. To recover the boundary curve H(1,ϕd) we have as supplementary condition fi(1) = 1. To satisfy the conditions

1382

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

fi(0) = 0 and fi(1) = 1 and to keep the combination of a concave and a convex curve, we propose: f 1 ðyÞ ¼

ð1 þ k1 Þ⋅y 1 þ k1 ⋅y

ð15Þ

f 2 ðyÞ ¼

1−expðk2 ⋅yÞ 1−expðk2 Þ

ð16Þ

with k1 and k2 positive parameters. To add flexibility to H(ϕa, ϕd), we assume ki to be linearly dependent on ϕd: k1 ðϕd Þ ¼ aa þ ch ⋅ϕd

ð17Þ

k2 ðϕd Þ ¼ ba þ dh ⋅ϕd

ð18Þ

Concluding, the H-function is described in the integrated PM space by: H ðϕa ; ϕd Þ ¼ H 1 ðϕa ; ϕd Þ þ H 2 ðϕa ; ϕd Þ

ð19Þ

with H1 and H2: H 1 ðϕa ; ϕd Þ ¼ ðM 1max −M1d ðϕd ÞÞ

ð1 þ k1 ðϕd ÞÞ⋅yðϕa ; ϕd Þ 1 þ k1 ðϕd Þ⋅yðϕa ; ϕd Þ

1− expðk2 ðϕd Þ⋅yðϕa ; ϕd ÞÞ H 2 ðϕa ; ϕd Þ ¼ ðM 2max −M2d ðϕd ÞÞ 1− expðk2 ðϕd ÞÞ

ρd

∂u ∂pc ∂p ¼ C c ¼ −∇g v ∂t ∂pc ∂t

g v ¼ −δv ∇pv

ð23Þ ð24Þ

with ρd the dry density of the material, u the moisture content (moisture mass versus dry material mass), gv the vapor flow, pv the water vapor pressure and δv the water vapor permeability, which depends on the vapor pressure. The vapor pressure is expressed by Kelvin's equation as:   pc pv ¼ pvsat ⋅exp ρw Rv T

ð25Þ

with pvsat the saturated vapor pressure, ρw the density of water, Rv the gas constant for water vapor and T the temperature. C, the moisture capacity, is a measure for the ability of the material to retain water. When describing hysteretic moisture behavior, the moisture content u is expressed by the function derived in Section 2:

ð20Þ

   pc uðpc Þ ¼ M ϕ ¼ exp ρw Rv T IPM

ð21Þ

with TIPM the temperature at which the moisture isotherms have been measured to determine the parameters describing the IPM function H. The capacity can then be calculated as:

M1d and M2d are given by Eqs. (6) and (7); k1 and k2 by Eqs. (17) and (18). The main adsorption isotherm Ma(ϕa) = H(ϕa,0) is given by: Ma ðϕa Þ ¼ H ðϕa ; 0Þ ¼ M1max ⋅

The vapor transport is represented through the following differential equations [22]:

ð1 þ aa Þ⋅ϕa 1−expðba ⋅ϕa Þ ð22Þ þ M2max ⋅ 1−expðba Þ 1 þ aa ⋅ϕa

and depends on three parameters l, aa and ba, with aa and ba positive values. The H-function itself needs 7 parameters: l, aa, ba, ad, bd, ch and dh, which have to be experimentally determined. The values have to be chosen so that the H-function is monotonically increasing with increasing ϕa and decreasing ϕd (as the moisture capacity or derivative of the moisture content to the relative humidity cannot be negative). As the parameters are based on real ad- and desorption data, phenomena such as film forming and capillary condensation are implicitly incorporated. An extended domain theory with explicit description of film adsorption can be found in [20]. 3. Numerical implementation Concrete after curing shows a high initial moisture content. In this paper, we assume that the concrete is protected from additional liquid ingress (such as rain or ground water) and is only subjected to changes in relative humidity and temperature. Moisture transport can then be described as vapor transport using a vapor permeability accounting for the enhancement of macroscopic vapor transport by microscopic liquid flow in pore regions showing capillary condensation [21]. Thus the vapor permeability includes transport in the gaseous and liquid phase, when capillary condensation occurs. Macroscopic liquid flow is not included, as we consider the concrete to be protected from liquid ingress. We further adopt the moisture and heat transfer equations for porous building materials as given in [22], with capillary pressure pc and temperature T as driving potential. The capillary pressure pc is defined as the difference in pressure between the liquid and the gaseous phase in the pore space. The derivation assumes that moisture storage depends only on capillary pressure.

C ¼ ρd

∂u ∂M ∂ϕ ¼ ρd ∂pc ∂ϕ ∂pc

ð26Þ

ð27Þ

and is numerically determined using the finite difference method applying a small perturbation on the capillary pressure. Eq. (25) allows rewriting Eq. (24) as [22,23]: g v ¼ −δv

  pv Lv lnðpv =pvsat Þ ð Þ Tγ−1 ∇T ∇pc −δv pv þ T ρw Rv T Rv T 2

ð28Þ

with Lv the latent heat of vaporization and γ the normalized derivative of the surface tension to the temperature. Heat transport is expressed by the conservation of energy equation [22]:   ∂h ∂T ∂u ¼ ρd ðcd þ cw uÞ þ cw ðT−T 0 Þρd ¼ −∇ g h;c þ g h;a ∂t ∂t ∂t

ð29Þ

with h the enthalpy, cd the specific heat of the dry material and cw the specific heat of water. The variation of the vapor and dry air enthalpy is considered negligible over time. T0 is the reference temperature for the enthalpy, being 0 °C (273.15 K). The heat flow is a combination of a conductive part gh,c = − λ ∇ T, with λ the thermal conductivity of the porous material, and an advective part gh,a = (cv(T − T0) + Lv) ⋅ gv, where cv is the specific heat of vapor. The boundary condition for vapor transport is given by: qv ¼ βðpve −pvs Þ

ð30Þ

with β the convective water vapor transfer coefficient, pve the water vapor pressure of the environment and pvs the water vapor pressure at the surface of the material. For heat transport, the boundary condition incorporates convective and radiative heat exchange by using a heat transfer coefficient α, and latent and sensible heat transport due to vapor exchange: qh ¼ α ðT e −T s Þ þ ðcv ðT s −T 0 Þ þ Lv Þ⋅qv

ð31Þ

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

with Te the environmental temperature and Ts the material surface temperature. The water vapor permeability δv depends on the vapor pressure and can be measured for different vapor pressure or relative humidity ranges under isothermal conditions using cup tests [24]. The nonlinear vapor permeability is described as function of RH by: δv ðϕÞ ¼ δair ⋅ða þ b·expðc ϕÞÞ

ð32Þ

with a, b and c parameters. The vapor permeability of air δair is given by Schirmer's equation [25,26]. The increase of the water vapor permeability with relative humidity is commonly attributed to the enhanced microscopic liquid water transport in water filled pores due to capillary condensation. Hence we can assume that the water vapor transport coefficient for hysteretic materials is a unique function of the amount of water present in the pore space, thus a unique function of the saturation degree or moisture content u [27,28]. Using the governing sorption isotherm during the vapor permeability measurements M(ϕ) gives:   −1 δv ðuÞ ¼ δv M ðuÞ

ð33Þ

The vapor resistance factor μ, commonly used in building physics, is then defined as the ratio between the vapor permeability of air and the vapor permeability of the material: μ ðuÞ ¼

δair δv ðuÞ

ð34Þ

As the vapor permeability is assumed to be a unique function of the moisture content, the vapor resistance factor is the same for all sorption curves when expressed in function of the moisture content. The transport Eqs. (23) and (29) are solved using the finite element method with capillary pressure and temperature as primary variables. In order to obtain a mass and energy conservative system of equations, a mixed form of the capacitive terms is used, as described in Janssen et al. [22]. A Newton–Raphson iterative scheme is employed and an adaptive time stepping algorithm is embedded to optimize the efficiency of the numerical solver. In order to calculate the hysteretic moisture content and capacity as per Eqs. (3) and (27) respectively, we need to keep track of the coordinates (ϕai, ϕdj) of the staircase describing the border between filled and unfilled domains for every integration point and nodal point of an element in the finite element mesh. We keep track of the real history at the last converged time step and the temporary history during the iterative procedure of the next time step building up a real and temporary history matrix. The moisture content and capacity are calculated for every iteration step based on the updated temporary history matrix. When convergence is reached, the real history is updated, i.e. it is replaced by the temporary history that led to the converged solution. A more extensive explanation on how to keep track of the history is given in Appendix A.

1383

in Table 1. The agreement with the measurement data is depicted in Fig. 5. The vapor permeability is taken from [30] for a concrete with a w/c ratio close to the concrete selected from [5] for the sorption data. The data were determined by cup tests [24]. The samples were first cured for more than 28 days. Three couples of RH were used for the cup tests: (90%, 75%), (75%, 54%) and (54%, 33%). For the couple (90%, 75%), the samples were equilibrated at 90% RH until equilibrium before the cup test. For the couple (75%, 54%), the samples were equilibrated at 90% RH for 2 months and then at 75% RH until equilibrium. For the couple (54%, 33%), the samples were equilibrated at 90% RH for 2 months, consecutively at 75% RH for 2 months and finally at 54% RH until equilibrium. Thus the vapor permeability was measured in the desorption regime. The parameters for Eq. (32) are given in Table 1 and the corresponding vapor resistance factor in function of the moisture content, using the main desorption isotherm (cf. Eqs. (33) and (34)), is presented in Fig. 6. This is the unique curve for the vapor resistance factor in function of the moisture content, applicable for all sorption curves. The liquid permeability is not considered, since we are studying wetting and drying processes restricted to the hygroscopic regime. 5. Simulation results and discussion First, we present simulation results to exemplify the proposed model assuming concrete exposed to an environment of constant temperature, 23 °C, and sinusoidally varying relative humidity. The environment determines the vapor pressure pve and the temperature Te in Eqs. (30) and (31). The influence of a sinusoidal variation with a period of 1 day, 6 days and 1 year is studied. The relative humidity of the environment is chosen to vary between 30 and 90% as this is the range in which the largest difference between main adsorption and desorption occurs for concrete (see also Fig. 5). As initial condition, two extremes are considered: (1) starting from low moisture content (‘dry state’), situated on the main adsorption curve at 30% relative humidity or (2) from high moisture content (‘wet state’), situated on the main desorption curve at 90% relative humidity. The latter case is analogous to the case of a concrete cured in sealed conditions and then exposed to variations of relative humidity. For initial condition ‘dry state’ we compare the reference material response using the hysteresis model with a simplified simulation only using the main adsorption curve; while for initial condition ‘wet state’, we compare the reference response using the hysteresis model with a simulation only using the main desorption curve. In a second example, the

4. Thermal and hygric properties of concrete We compose a material property dataset for a common concrete based on material properties from different literature sources. The concrete in this study has a density of 2270 kg/m3 and a thermal capacity of 940 J/kgK [29]. The thermal conductivity is moisture dependent and calculated as λ=1.5+0.011· w W/mK [29], with w the moisture content in kg/m3. For the hysteretic moisture behavior, sorption data for an average structural concrete are taken (i.e. BO-concrete in [5]). Main ad- and desorption curves as well as some intermediate scanning curves were measured at 23 °C as shown in Fig. 5. The parameters for the H-function (Eq. (19)) are determined based on these data and given

Fig. 5. Hysteretic moisture data (points) of concrete [5] and the IPM model fit (lines).

1384

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

Table 1 Parameters for the hysteretic moisture retention and the vapor resistance factor. Parameters IPM model (Eq. (19))

Parameters vapor resistance model (Eq. (32))

Mmax l ad bd aa ba ch dh

a b c

3.86 0.347 2 1.57 1.1 10.6 −1.05 −10.55

0.0048 7.895e-6 9.247

behavior of concrete under real climatic conditions is calculated, using hourly data of relative humidity and temperature in Essen, Germany. As concrete is initially in a wet state if sufficiently well cured on site, an initial ‘wet state’, equal to the initial outdoor conditions in Essen (91% RH, 0.5 °C), is considered. The simulations are performed with a one-dimensional mesh. The first element has a length of 20 μm, the element size is gradually increasing with a factor 1.048. The total length amounts to 5 cm (100 elements) when applying sinusoidal relative humidity fluctuations with a 1-day or 6-day period and to 30 cm (140 elements) when imposing relative humidity fluctuations with a 1-year period and climatic boundary conditions. A heat transfer coefficient α (Eq. (31)) of 7.7 W/m2K is typically used at indoor surfaces (i.e. for the isothermal simulations at 23 °C) and a value of 25 W/m2K outdoor (i.e. for the simulations with climatic boundary conditions) [31]. The convective part of these heat transfer coefficients amounts 2.5 W/m2K and 20 W/m2K for indoor and outdoor conditions respectively [32]. Using the Chilton–Colburn analogy [33], a value of 1.8∙10−8 s/m (for the isothermal simulations) or 14∙10−8 s/m (for the simulations with the climatic boundary conditions) is used for the mass transfer coefficient β (Eq. (30)). The environmental conditions are applied on one side of the mesh; on the other boundary, no flow conditions (i.e. zero flow) are imposed. The simulations are run until a quasi-steady state is reached inside the concrete material. These quasi-steady state conditions are used to analyze the different results using the hysteresis model and main adsorption or main desorption curves. For the sinusoidal relative humidity fluctuations with a period of 1 day and 6 days, the simulations are done for 600 days, for the sinusoidal relative humidity fluctuations with a period of 1 year and the Essen climate, the simulations comprise 40 years. The calculated response at different positions in the 5 cm thick concrete after being subjected for 600 days to the sinusoidal relative

Fig. 6. Vapor resistance factor versus moisture content.

humidity fluctuations with a 1-day period, is given in Fig. 7, for the two different initial conditions, starting from the ‘dry state’ (30% RH, main adsorption curve) and from the ‘wet state’ (90% RH, main desorption curve). The figures are based on the results of the last (600th) day. The moisture content u, the corresponding vapor resistance factor μ and moisture capacity C are shown as function of the relative humidity. Comparing the moisture content variations as determined by the hysteresis model and the main adsorption curve (Fig. 7a, initial ‘dry’ material), we observe an important difference for the concrete surface; deeper inside the concrete, the material approaches the main adsorption curve. The difference is larger when comparing the moisture content variations between hysteresis model and main desorption curve (Fig. 7d, initial ‘wet’ material), where the material response of the first 5 mm clearly deviates from the main desorption curve. The differences can be explained by the fact that the moisture capacity differs significantly from the main ad/desorption curves (Fig. 7c and f) and especially from the main desorption curve. Consequently the vapor resistance factor, thus moisture and gas permeability, which highly depend on the moisture saturation, also differ from the main ad/desorption values (Fig. 7b and e). The vapor resistance factor deviates less from the main adsorption curve than from the main desorption curve (Fig. 7b and e) because of the broader range of μ-values covered by the main desorption curve in the range of 30–90% relative humidity. Fig. 8 shows the material response in function of depth after exposure to sinusoidal relative humidity fluctuations with a 6-day period, starting from the ‘wet state’, for the last 6 days of the simulation. The results using the hysteresis model are given in Figs. 8a–b–c, the responses using the main desorption curve are shown in Figs. 8d–e–f. The profiles show that the fluctuations in moisture properties decrease with depth, thus that the relative humidity fluctuations decrease deeper inside the material. This allows to estimate the penetration depth, i.e. an indication of the depth penetrated by moisture into the concrete structure under transient and cyclic conditions [34]. We observe that the penetration depth is about 15 mm for a variation with a period of 6 days and the depth is independent of simulating with the hysteresis model or the main sorption curve. The moisture capacity values obtained when using the hysteresis model, given in Fig. 8c, show that beyond a depth of 33 mm the material still follows the main desorption curve; thus beyond this depth the material did not experience any adsorption step during the 600 days of simulation. The same curves can be drawn for periods of 1 day and 1 year. The penetration depths are respectively 6 mm and 150 mm. Thus relative humidity cycles with longer periods lead to variations penetrating deeper into the material. Figs. 9a–b–c show the same curves as Figs. 7d–e–f, but for the relative humidity fluctuations with a 1-year period, starting from the ‘wet state’. The loops represent the response during the final calculated year (40th year). The same observations as described for Fig. 7 hold, but the yearly conditions penetrate much deeper into the structure. Figs. 9d–e–f give the same graphs, when the structure is exposed to the real climatic conditions of Essen. The climatic variations penetrate less deep than the isothermal relative humidity fluctuations with a 1-year period as a climate comprises a superposition of daily and yearly variations. The penetration depth is around 100 mm. Again a significant difference is observed between simulating with the main sorption curves or the hysteresis model. The moisture behavior is characterized by a cloud of values and do not form a single loop as in the isothermal simulations. The deeper inside the concrete material, the thinner the cloud gets. Remark that the simulation results using the main ad- or desorption curve also form a thin cloud and not a single line: due to the temperature fluctuations, different moisture contents can be found at the same relative humidity. Depending on the climate, hysteresis will be important over a depth of several centimeters in a concrete structure. Especially for locations within the penetration depth, the difference in moisture behavior can be important when using only the main sorption curve compared to applying the hysteresis model. It is also impossible to

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

1385

Fig. 7. Moisture content, vapor resistance factor and moisture capacity in function of relative humidity in a 5 cm thick concrete sample during the last day of 600 days exposure to relative humidity fluctuations with a 1-day period: (a)–(b)–(c) starting from ‘dry state’ (30% RH, main adsorption) with fluctuations between 30 and 90% RH; (d)–(e)–(f) starting from ‘wet state’ (90% RH, main desorption) with fluctuations between 90 and 30% RH.

determine a single intermediate sorption curve that would cover the hysteretic behavior for all positions uniquely: although the capacities lie in the same range and the hysteretic loops are not very broad, the moisture content values at different depths do not approach a single curve. As a consequence, following Eq. (34), the vapor resistance factors differ. Thus a hysteresis model is indispensable for predicting accurately the material behavior when being exposed to environmental fluctuations. Degradation of concrete occurs in the outer concrete layer. It is in this same outer layer that hysteresis plays its major role and influences the permeability of CO2 and chlorides. The moisture content, or saturation

degree, is a good parameter to evaluate the risk for carbonation or chloride ingress, as the water permeability is directly related to the moisture content through Eq. (33). With increasing saturation degree, the water permeability will increase and consequently gas permeability will decrease since more pores are filled with water. Thus the diffusion of CO2 can be underestimated due to an overestimation of the saturation degree and consequently the risk for carbonation will not be predicted correctly. This means that if one would use the main desorption curve to simulate the behavior of an initial ‘wet’ concrete, the diffusion of CO2 is underestimated compared to using the hysteresis model. On the other hand, the liquid permeability will be underestimated when the saturation

1386

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

Fig. 8. Moisture content, vapor resistance factor and moisture capacity as function of depth for a 5 cm thick concrete sample during the last 6 days of 600 days exposure to relative humidity fluctuations between 90 and 30% RH with a 6-day period: (a)–(b)–(c) using the hysteresis model and (d)–(e)–(f) using the main desorption curve to model the hygroscopic response. The response after 594 and 600 days is the same as an equilibrium state is reached.

degree is underestimated and consequently the risk for ingress of chlorides is underestimated. Therefore if one would use the main adsorption curve instead of the hysteresis model to describe the behavior of an initial ‘dry’ concrete, the moisture permeability is underestimated and so is the corrosion risk due to chloride infiltration. 6. Conclusion An analytical model has been developed to describe hysteretic moisture behavior in porous materials based on the independent domain theory. This model has been integrated in a finite element code for coupled heat and mass transfer. It allows assessing the influence of moisture hysteresis in a material when exposed to changing temperature and

relative humidity conditions. The methodology has been successfully applied on examples of concrete subjected to different environmental conditions. It illustrates how carbonation or chloride ingress is difficult to estimate when omitting moisture hysteresis in the outer concrete layer. Future concrete research could focus on the acquisition of more experimental data on moisture hysteresis and the exact relationship between gas permeability and carbonation risks, chloride diffusion in non-saturated conditions and corrosion risks. The presented simulation results demonstrate that the developed model would then be able to quantify the behavior of different concrete mixtures under different climatic conditions. The proposed approach could provide additional information when selecting a concrete for specific environmental conditions.

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

1387

Fig. 9. Moisture content, vapor resistance factor and moisture capacity in function of relative humidity in a 30 cm thick concrete sample during the last year of 40 years exposure. The sample is initially in a ‘wet state’ and subjected to (a)–(b)–(c) relative humidity fluctuations between 90 and 30% RH with a 1-year period; (d)–(e)–(f) the outdoor climate of Essen, Germany.

Role of the funding source This work has been funded in part by the Lafarge Research Center within the project ‘An investigation of modeling moisture hysteresis in the hygroscopic regime for porous materials’. The authors are grateful for the financial support of the Swiss National Science Foundation (SNF) under Grant No. 125184. Appendix A For each integration point and each nodal point of an element in the finite element mesh, the real history at the last converged time step,

and the temporary history during the iterative procedure of the next time step is stored. To keep track of the coordinates of the staircase describing the point's history in the PM-space, a real history and temporary history matrix is defined for each point. Each matrix consists of two columns and an adaptable number of rows. The first column consists of the ϕa-values, the second column of the ϕd-values. The uneven rows represent an adsorption step, so the corresponding H-value needs to be added to the moisture content calculation; the even rows represent a desorption step, meaning that the corresponding H-value needs to be subtracted. The total moisture content is found by summing and subtracting all H-values calculated with the coordinates given in each row. Let's define ϕ0 as the relative humidity of the last converged time step

1388

H. Derluyn et al. / Cement and Concrete Research 42 (2012) 1379–1388

and ϕa,k,0 and ϕd,k,0 as the relative humidity values on the last row k of the real history matrix. During every iteration step it the temporary history matrix, containing the elements ϕa,m and ϕd,m (with m being the row number), is built up as follows:

If ϕit > ϕ0 l=k do while ϕa,l,0 b ϕit l = l − 1 (delete lth row) enddo l=l+1 for m = 1,l − 1 ϕa,m = ϕa,m,0 ϕd,m = ϕd,m,0 end ϕa,l = ϕit ϕd,l = ϕd,l−1,0 endif

If ϕit b ϕ0 l=k do while ϕd,l,0 >ϕit l = l − 1 (delete lth row) enddo l=l+1 for m = 1,l − 1 ϕa,m = ϕa,m,0 ϕd,m = ϕd,m,0 end ϕa,l = ϕa,l−1,0 ϕd,l = ϕit endif

When convergence is reached, the real history is updated, i.e. it is replaced by the temporary one. Depending on the relative humidity loading protocol, the history matrix can vary from having only 1 row to a large number of rows. By updating the temporary and real history matrix after respectively every iteration step and every time step, their sizes are limited to the minimum and an optimal memory usage is ensured. Other applications of the PM-space, namely in the field on non-linear mechanics, are described in [35] and [36]. References [1] F.P. Glasser, J. Marchand, E. Samson, Durability of concrete — degradation phenomena involving detrimental chemical reactions, Cem. Concr. Res. 38 (2008) 226–246. [2] T. Ishida, K. Maekawa, M. Soltani, Theoretically identified strong coupling of carbonation rate and thermodynamic moisture states in micropores of concrete, J. Adv. Concr. Technol. 2 (2) (2004) 213–222. [3] N. De Belie, J. Kratky, S. Van Vlierberghe, Influence of pozzolans and slag on the microstructure of partially carbonated cement paste by means of water vapour and nitrogen sorption experiments and BET calculations, Cem. Concr. Res. 40 (2010) 1723–1733. [4] A. Anderberg, L. Wadsö, Moisture properties of self-levelling flooring compounds. Part II. Sorption isotherms, Nord. Concr. Res. 32 (2) (2004) 16–30. [5] V. Baroghel-Bouny, Water vapour sorption experiments on hardened cementitious materials. Part I: essential tool for analysis of hygral behaviour and its relation to pore structure, Cem. Concr. Res. 37 (2007) 414–437. [6] I. Maruyama, G. Igarashi, Mechanism of moisture sorption hysteresis of hardened cement paste, Cem. Sci. Concr. Technol. (Trans. JCA) 64 (2011) 96–102. [7] H. Ranaivomanana, J. Verdier, A. Sellier, X. Bourbon, Toward a better comprehension and modeling of hysteresis cycles in the water sorption–desorption process for cement based materials, Cem. Concr. Res. 41 (2011) 817–827. [8] F. Brue, C.A. Davy, F. Skoczylas, N. Burlion, X. Bourbon, Effect of temperature on the water retention properties of two high performance concretes, Cem. Concr. Res. 42 (2012) 384–396. [9] T. Ishida, K. Maekawa, T. Kishi, Enhanced modeling of moisture equilibrium and transport in cementitious materials under arbitrary temperature and relative humidity history, Cem. Concr. Res. 37 (2007) 565–578. [10] R.M. Espinosa, L. Franke, Inkbottle pore-method: prediction of hygroscopic water content in hardened cement paste at variable climatic conditions, Cem. Concr. Res. 36 (2006) 1954–1968.

[11] D.H. Everett, W.I. Whitton, A general approach to hysteresis, Trans. Faraday Soc. 48 (1952) 749–757. [12] D.H. Everett, F.W. Smith, A general approach to hysteresis part 2: development of the domain theory, Trans. Faraday Soc. 50 (1954) 187–197. [13] D.H. Everett, A general approach to hysteresis part 3: a formal treatment of the independent domain model of hysteresis, Trans. Faraday Soc. 50 (1954) 1077–1096. [14] D.H. Everett, A general approach to hysteresis part 4: an alternative formulation of the domain model, Trans. Faraday Soc. 51 (1955) 1551–1557. [15] F. Preisach, Über die magnetische Nachwirkung, Z. Phys. 94 (1935) 277–302. [16] B. Bary, A. Sellier, Coupled moisture – carbon dioxide – calcium transfer model for carbonation of concrete, Cem. Concr. Res. 34 (2004) 1859–1872. [17] M. Thiery, Modélisation de la carbonatation atmosphérique des matériaux cimentaires: Prise en compte des effets cinétiques et des modifications microstructurales et hydriques, PhD Thesis, ENPC, France, 2005. [18] E. Stora, C. Conrardy, R. Barbarulo, J. Chen, F. Thoué, A reactive transport model based on thermodynamic computations to estimate the atmospheric carbonation of cementitious materials, In: V. L'Hostis, R. Gens, C. Gallé (Eds.), NUCPERF 2009 — Long Term Performance of Cementitious Barriers and Reinforced Concrete in Nuclear Power Plants and Waste Management, Cadarache, France, 2009, pp. 193–202. [19] J.D. Mayergoyz, Hysteresis models from the mathematical and control theory points of view, J. Appl. Phys. 57 (1985) 3803–3805. [20] B. Coasne, K.E. Gubbins, R.J.-M. Pellenq, Domain theory for capillary condensation hysteresis, Phys. Rev. B 72 (2005) 024304. [21] Z.P. Bazant, L.J. Najjar, Nonlinear water diffusion in nonsaturated concrete, Mater. Struct. 5 (1972) 3–20. [22] H. Janssen, B. Blocken, J. Carmeliet, Conservative modelling of the moisture and heat transfer in building components under atmospheric excitation, Int. J. Heat Mass Transfer 50 (2007) 1128–1140. [23] C.E. Hagentoft, In: HAMSTAD — WP2 Modeling, Version 4, Report R-02:9, Department of Building Physics, Chalmers University of Technology, Gothenburg, Sweden, 2001, p. 12. [24] EN ISO 12572, Hygrothermal Performance of Building Materials and Products — Determination of Water Vapour Transmission Properties, CEN, Brussels, 2001. [25] R. Schirmer, Die Diffusionszahl von Wasserdampf-Luftgemischen und die Verdampfungs-geschwindigkeit, ZVDI Beheift Verfahrenstech. 6 (1938) 170. [26] F. Ochs, H. Heidemann, H. Müller-Steinhagen, Effective thermal conductivity of moisture insulation materials as a function of temperature, Int. J. Heat Mass Transfer 51 (2008) 539–552. [27] S. Roels, H. Janssen, J. Carmeliet, J. Diepens, M. de Wit, Hygric buffering capacities of uncoated and coated gypsum board, In: Proceedings of the Third International Building Physics Conference, Montreal, Canada, 2006, pp. 27–32. [28] H. Derluyn, H. Janssen, J. Diepens, D. Derome, J. Carmeliet, Hygroscopic behavior of paper and books, J. Build. Phys. 31 (2007) 9–34. [29] M.K. Kumaran, Heat, air and moisture transfer in insulated envelope parts, In: Final Report IEA Annex 24, Vol. 3, Task 3: Material Properties, Laboratory of Building Physics, Civil Engineering Department, KULeuven, Belgium, 1996, pp. 14–15. [30] G. Arliguie, H. Hornain, et al., In: GranDuBé: Grandeurs associées à la Durabilité des Bétons, Presses de l'Ecole Nationale des Ponts et Chaussées, Paris, 2007, p. 102. [31] SN 520 180, In: Isolation thermique et protection contre l'humidité dans les bâtiments, Société suisse des ingénieurs et des architectes, Zürich, 1999, p. 29. [32] prEN15026, Hygrothermal performance of building components and building elements, In: Assessment of Moisture Transfer by Numerical Simulation, CEN, Brussels, 2004, p. 15. [33] T.H. Chilton, A.P. Colburn, Mass transfer (absorption) coefficients, Ind. Eng. Chem. 26 (11) (1934) 1183–1187. [34] M.J. Cunningham, Effective penetration depth and effective resistance in moisture transfer, Build. Environ. 27 (1992) 379–386. [35] K.E.A. Van Den Abeele, J. Carmeliet, P.A. Johnson, B. Zinszner, Influence of water saturation on the nonlinear elastic mesoscopic response in earth materials and the implications to the mechanism of nonlinearity, J. Geophys. Res. Solid Earth 107 (6) (2002) 2121. [36] V. Aleshin, W. Desadeleer, M. Wevers, K. Van Den Abeele, Characterization of hysteretic stress–strain behavior using the integrated Preisach density, Int. J. Nonlinear Mech. 43 (2008) 151–163.