Numerical investigation of the effects of freezing on micro-internal damage and macro-mechanical properties of cement pastes

Cold Regions Science and Technology 106–107 (2014) 141–152

Contents lists available at ScienceDirect

Cold Regions Science and Technology journal homepage: www.elsevier.com/locate/coldregions

Numerical investigation of the effects of freezing on micro-internal damage and macro-mechanical properties of cement pastes Lin Liu a,⁎, Shengxing Wu a, Huisu Chen b, Zhao Haitao a a b

College of Civil and Transportation Engineering, Hohai University, Nanjing 210098, China Jiangsu Key Laboratory of Construction Materials, School of Material Science and Engineering, Southeast University, Nanjing 211189, China

a r t i c l e

i n f o

Article history: Received 13 August 2013 Accepted 8 July 2014 Available online 15 July 2014 Keywords: Freezing effects Hydraulic pressure Crystallization pressure Micro-internal damage Cement paste

a b s t r a c t It is of significance to investigate the deterioration of cement-based systems subjected to freezing temperatures including intrinsic microstructure change and reduction in mechanical resistance. From the microstructure of materials in question, few computer models take the coupled effects of water expulsion, in-pore crystallization and cryosuction into account to assess the realistic internal damage due to frost action. By considering the coupled effects of thermal action, hydraulic pressure and crystallization pressure, this study has established a microstructure-based model to investigate the effects of freezing on micro-internal damage and the macromechanical property changes of cement-based systems. Two parts of work are carried out: (1) by coupling a microstructure model HYMOSTRUC3D with a 3D lattice fracture model, the freezing of water in pores and the response of cement paste to the water freezing are captured. Changes in the microstructure of cement paste, represented by the creation of microcracks, are illustrated. (2) Changes in the mechanical properties of cement paste after freezing/thawing are estimated, and the estimated Young's modulus and volume deformation during freezing are compared with those by experiment. In the discussion chapter, from the aspects of thermodynamic considerations, damage sources, crack pattern and deformation, the hydraulic pressure model, ice crystallization pressure model and the newly established model in this study are discussed. © 2014 Elsevier B.V. All rights reserved.

1. Introduction Frost action (freeze–thaw cycles) is one of the major causes of the deterioration of concrete in cold wet climates, especially for hydraulic concrete. Concrete in some service environments can be subjected to freezing/thawing cycles coupled with exposure to deicer (Shi et al., 2010) or sulfate attack (Miao et al., 2002). To mitigate the premature deterioration of such concrete, enhanced understanding of frost action would also be very beneficial. Investigations on the deterioration of cement-based materials exposed to freezing/thawing have revealed a reduction in the mechanical resistance and an increase in the permeability, which are associated with the intrinsic changes of the microstructure (Diao et al., 2013; Liu and Wang, 2012; Liu et al., 2013; Shang and Song, 2013; Wardeh et al., 2010; Yang et al., 2013). The microstructure changes in hardened cement-based systems after freeze– thaw cycles have been studied by using a variety of experimental techniques, such as scanning electron microscope (SEM) (Skripkiūnas et al., 2013), X-ray diffraction (XRD) (Skripkiūnas et al., 2013), X-ray computed tomography (X-ray CT) (Liu, 2012), ultrasonic imaging (Molero et al., 2012) and so on. In addition, instead of obtaining the microstructure ⁎ Corresponding author. Tel./fax: +86 25 83786183. E-mail addresses: [email protected] (L. Liu), [email protected] (S. Wu), [email protected] (H. Chen), [email protected] (Z. Haitao).

http://dx.doi.org/10.1016/j.coldregions.2014.07.003 0165-232X/© 2014 Elsevier B.V. All rights reserved.

images of cement-based materials, Li and Weiss et al. (Li et al., 2012) have evaluated the internal damage of mortar during freezing–thawing processes by acoustic emission. Studies (Powers, 1949; Scherer, 1999; Scherer and Valenza, 2005; Zeng, 2012) on the mechanisms governing frost damage of hydrated cement-based systems for several decades have consistently found that harmful stresses could result from: a) hydraulic pressure due to a 9% expansion in volume accompanying ice formation; b) crystallization pressure, induced by the growth of crystals in pores and their interaction with pore walls; c) the mismatch of thermal effects between ice and solid phases, etc. According to the above frost damage mechanisms, numerous recent papers have attempted to simulate the freezing behavior and to predict the deterioration of cement-based materials (Bazant et al., 1988; Coussy and Monteiro, 2008; Dai et al., 2013; Duan et al., 2013; Hain and Wriggers, 2008; Liu et al., 2011, 2014; Wardeh and Perrin, 2008; Zeng et al., 2010, 2011; Zuber and Marchand, 2000). By linking thermodynamics to poromechanics, the macroscopic deformation during freezing (Zeng et al., 2010, 2011), the influences of pore size distribution and air-void spacing (Coussy and Monteiro, 2008) have been analyzed. By mathematical models, the deterioration of cement-based materials has been investigated in the literature (Bazant et al., 1988; Wardeh and Perrin, 2008; Zuber and Marchand, 2000). By coupling a microstructure of cement-based materials with a fracture model, the micro-internal damage of cement-based materials

142

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

due to frost action has been captured in the literature (Dai et al., 2013; Duan et al., 2013; Hain and Wriggers, 2008; Liu et al., 2011, 2014). Although these microstructure-based models, in which the microstructure of cement-based materials is taken into account, can provide the internal damage assessment, the damage sources are assumed too simplistic where either crystallization pressure (Dai et al., 2013; Liu et al., 2011) or hydraulic pressure (Liu et al., 2014) is playing a dominant role. Few models take the coupled effects of hydraulic pressure and crystallization pressure into account to assess the realistic internal damage due to frost action. By considering the coupled effects of thermal action, hydraulic pressure and crystallization pressure, this study is going to establish a microstructure-based model to assess the micro-internal frost damage, and to estimate the mechanical property changes of cement-based systems after freezing/thawing. In the following, based on thermodynamic and porous medium considerations, theoretical backgrounds are first presented, including a description of in-pore phase transition, a conceptual explanation for the freezing behavior of porous solids and a discussion about hydraulic pressure model and crystallization pressure model. Then, two parts of work are carried out: (1) by coupling a microstructure model HYMOSTRUC3D with a 3D lattice fracture model, the freezing of water in pores and the response of cement paste to the water freezing are captured. Changes in the microstructure of cement paste, represented by the creation of microcracks, are illustrated. (2) Freezing effects on the mechanical properties of cement paste are investigated, and the simulated Young's modulus and volume deformation is compared to experiments. Moreover, the micro-internal damage of cement paste under different frost mechanisms is discussed. 2. Theoretical background

Fig. 1. In-pore phase transition of ice and water. rp — the radius of a pore entry. δ — thickness of the liquid film between the crystal and the pore wall. κf — curvature of the crystal/liquid interface at the free end of the crystal. κn — curvature of the crystal/liquid interface at the non-free sides of the crystal. PD — disjoining pressure. Reproduced after Liu et al. (2011).

Owing to the presence of intermolecular forces, the disjoining pressure PD is introduced by a layer of liquid film between the crystal and the pore wall with a thickness of δ. Combining Eq. (2) with Eq. (3), PD can be obtained, P D ¼ γCL ðκ f −κ n Þ:

On the assumption of PL = Patm, the well-known Gibbs–Thomson equation is given by (Scherer, 1999), rc ¼ −

The thermodynamics of in-pore crystallization can provide us a fundamental understanding for the freezing behavior of porous solids. Therefore, a description of in-pore crystal–liquid phase transition is first presented in this chapter. A conceptual explanation for the freezing behavior of porous solids is followed by the presentation about hydraulic pressure model and crystallization pressure model.

ð4Þ

2γCL cosθ ΔSm ðT m −T Þ

ð5Þ

where rc is the radius of the crystal at its free end at temperature T. θ is the contact angle between the crystal and the pore wall. The Gibbs– Thomson equation provides a scientific basis for calculating the freezing order of pores as temperature decreases. Assuming that the contact angle is 180° and by considering the thickness of the liquid film between the crystal and the pore wall δ, Eq. (5) can be rewritten as:

2.1. In-pore crystal–liquid phase transition The in-pore crystal–liquid phase transition is totally different from that between bulk water and ice. At low temperatures, in order to maintain the equilibrium of chemical potential between ice crystals and liquid water, the Thomson equation (Coussy, 2010) is satisfied with, 0

0

P C −P atm −ðP L −P atm ÞV L =V C ¼ ðT m −T ÞΔSm

ð1Þ

where PL and PC are the pressure in the liquid and in the crystal, MPa. Patm is the atmospheric pressure, MPa. Tm is the freezing temperature where both liquid and crystal phases are at the atmospheric pressure, 0 K. ΔSm is the melting entropy per unit volume of crystal, MPa/K. V L 0 and V C are molar volumes of liquid water and ice. For a crystal of ice confined in a small pore (see Fig. 1), at its free end, it satisfies with the Laplace equation's (Liu et al., 2011; Scherer, 1999; Scherer and Valenza, 2005), P C ¼ P L þ γCL κ f

ð2Þ

where γCL is the energy of the crystal/liquid interface, J/m2. κf is the curvature of the crystal/liquid interface at the free end of the crystal, m−1. At the non-free sides of the crystal where its curvature κn is far smaller than that at the free end κf, the mechanical equilibrium of the crystal and liquid interface is maintained by involving a disjoining pressure PD, P C ¼ P L þ γCL κ n þ P D :

ð3Þ

rp ¼

2γCL þδ ðT m −T ÞΔSm

ð6Þ

where rp is the radius of a pore entry. Eq. (6) gives the penetration threshold of ice crystals at temperature T. The lower the temperature is, the smaller the pore will be which ice crystals can grow into, and a more significant surface energy cost is required with regard to the volume to be frozen. 2.2. Freezing behavior of porous solids For water-saturated porous solids, besides the thermal effect, their freezing behavior could be involved in the coupled effects of expulsion (i.e., the buildup of hydraulic pressure) and cryosuction during temperature decreasing. As illustrated in Fig. 2a, when ice crystal nucleates at the temperature Ta (Ta b Tm), a hydraulic pressure PL is built up (PL N 0) because a crystal of ice is formed in the biggest pore accompanying a 9% volume expansion of water-to-ice. The unfrozen water is expulsed from the freezing site to unfrozen sites. The pore wall which is in contact with the ice crystal is subjected to stress σ1 with, σ 1 ¼ PL þ PD :

ð7Þ

The pore wall remaining unfrozen is subjected to stress σ2 with, σ 2 ¼ PL :

ð8Þ

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

143

(expulsion effect). As a result, both the pores filled with ice 1 and ice 2 are enlarged. As to the pore remaining unfrozen, it could be enlarged or shrink depending on the liquid pressure PL. The quantification of PL is associated with the thermodynamics of in-pore phase transition and the mass balance of water. At any given temperatures, the mass balance of water can be expressed as, nw ¼ nC þ nL

ð9Þ

where nw, nC, and nL are the total molar of water, the molar of ice crystals and the molar of liquid water, respectively. Take the case in Fig. 2 as an example, from temperature Tb to Tc (Tb N Tc), the molar change of nC, and nL can be expressed as, nC;T c ¼ nC;T b þ nL→C;T b →T c ¼ SC;T c nw

ð10Þ

nL;T c ¼ nL;T b −nL→C;T b →T c ¼ SL;T c nw

ð11Þ

where nL → C,Tb → Tc is the molar change of water from liquid to crystals when the temperature is changed from Tb to Tc. nL → C,Tb → Tc is related to the phase transition in small pores (i.e., formation of ice 2 in pore 2) and the cryosuction in big pores (i.e., ice 1 sucks water). 2.3. Hydraulic pressure model and crystallization pressure model

Fig. 2. The freezing behavior of porous solids. (a) Ice nucleates at temperature Ta (Ta b Tm), hydraulic pressure is built up because of the expulsion of unfrozen water. (b) Temperature decreases to Tb (Tb b Ta), the pressure in the ice crystal PC increases and the volume of the ice crystal increases by sucking liquid water in order to maintain the equilibrium of chemical potential. As a result, the liquid pressure PL decreases because of the water flow to the ice crystal, that is, cryosuction effect. (c) Temperature subsequently decreases to Tc (Tc b Tb b Ta), the pressure in the ice crystal PC increases continuously. Ice 1 penetrates into its connected pore and ice 2 is formed. Ice 1 would suck water by the liquid film. Ice 2 may expulse liquid water because of its volume expansion at phase transition. The liquid pressure PL could be increased or decreased due to the coupled effects of cryosuction and expulsion.

When the temperature decreases to Tb where Tb b Ta, see Fig. 2b, the pressure in the ice crystal PC increases (according to Eq. (1)) and the volume of the ice crystal increases by sucking liquid water in order to maintain its chemical potential equilibrium. As a result, the liquid pressure PL decreases because of the water flow to the ice crystal, that is, cryosuction effect. Compared to the state at temperature Ta, the pore filled with ice is enlarged and the pores filled with unfrozen water may get smaller. As the temperature subsequently decreases to Tc where Tc b Tb b Ta, the pressure in the ice crystal PC increases continuously, see Fig. 2c. At this point, the original ice (named ice 1) penetrates into its connected pore and a new ice (named ice 2) is formed. The ice 1 will suck water to get larger by the liquid film (cryosuction effect), while ice 2 may expulse liquid water because of the volume expansion at phase transition

The hydraulic pressure model proposed by Powers (1949) indicates that the expulsion of excess water accompanying the water-to-ice transformation causes a pressure built up in pore, that is, hydraulic pressure. Hydraulic pressure is mainly attributed to the frost damage of sealed porous solids where the extra liquid cannot be driven out. In the case of nucleation at a low temperature, hydraulic pressure plays a dominant role. This is because at nucleation a large amount of ice crystals are immediately formed, and a high hydraulic pressure PL is suddenly produced. Ice nucleation in porous solids is classified as two categories: homogeneous nucleation and heterogeneous nucleation. Homogeneous nucleation is unlikely in porous solids because of a high energy barrier to overcome. Coussy et al. indicated that homogeneous nucleation theoretically occurs below −45 °C (Coussy, 2010). Heterogeneous nucleation occurs mostly in the freezing porous solids, which depends on the contact angle between the ice and the nucleating substrate. For mortars with and without nucleating agent, heterogeneous nucleation was observed at − 6 °C and − 10.6 °C, respectively (Sun and Scherer, 2010a,b). Defining the molar percentage of water transformed into ice as crystal saturation degree SC, and liquid saturation degree SL equals to 1 − SC, the hydraulic pressure PL and crystallization pressure PC are satisfied with, PL ¼ K L

PC ¼ K C

SL V 0 −V 0L V 0L

ð12Þ

  0 0 1 þ V C =V L SC V 0 −V 0C V 0C

ð13Þ

where KL and KC are the bulk modulus of water and ice, respectively. V0 is the total volume of pores before freezing; VL′ is the volume of pores remaining unfrozen after freezing; and VC′ is the volume of frozen pores after freezing. If the total volume of pores after freezing is equal to that before freezing, 0

0

VC þ VL ¼ V0

ð14Þ

144

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

combined with Eq. (1), the hydraulic pressure built up in a watersaturated porous solids can be quantified as (Coussy, 2010),

σ1 which is equal to the disjoining pressure PD. According to Eq. (4), we have (Liu et al., 2011),

  0 0 P L ¼ SC 1−V L =V C

σ 1 ¼ P D ¼ γCL κ f ð1−λÞ

KLKC : SL K C þ SC K L

ð15Þ

The hydraulic pressure PL could be overestimated in Eq. (15), because the pore volume after freezing could be enlarged, i.e. VC′ + VL′ N V0. Hydraulic pressure model is capable of providing an order of magnitude of stress and of the critical spacing factor for entrained air bubbles. However, it alone cannot explain for the subsequent work by Powers and co-works, that is, air-entrained pastes exhibit shrinkage at a given temperature below freezing (Coussy and Monteiro, 2008). In fact, the hydraulic pressure model only focuses on the expulsion effect and the cryosuction effect is not taken into account. If the cryosuction effect is ignored, the decrease of hydraulic pressure PL at low temperatures (for example, the case in Fig. 2b) cannot be captured. The crystallization pressure model (Scherer, 1999; Scherer and Valenza, 2005) indicates that the crystallization pressure in mesopores is mainly responsible for the frost damage of porous solids. Serious damage occurs at low temperatures, at which a high crystallization pressure PC is maintained. Evidences for the importance of crystallization pressure are the experimental results on freezing cement pastes filled with benzene (Beaudoin and MacInnis, 1974). Since the volume of benzene decreases at nucleation, no hydraulic pressure can be produced during the freezing process. The dilation of samples after nucleation can be attributed to crystallization pressure (Scherer and Valenza, 2005). For air-entrained porous solids or non-sealed porous solids where the excess water can be driven out, the hydraulic pressure can be negligible, that is, PL = Patm. Take Patm = 0 as a reference, the crystallization pressure PC is proportional to the temperature decrease (Tm − T) according to Eq. (1). In this case, only the frozen pore walls are subjected to stress

ð16Þ

where λ is a pore shape factor, which is in the range of 0–1. The pore shape factor reflects the curvature differences between the frozen pore wall and the free end of the crystal. For a large pore with small entries, its pore shape factor can approach to 0 at low temperatures when ice crystals can penetrate into the pore. Crystallization pressure model subjected to surface energy effect can consider the cryosuction effect, while the expulsion effect due to the volume expansion of water-to-ice has not been taken into account. Both hydraulic pressure model and crystallization pressure model are special cases for the frost damage of porous solids. It is of significance to study the freezing behavior of cement-based systems under the coupled effects of hydraulic pressure and crystallization pressure. This paper is going to simulate the internal damage of cement paste as well as its macro-mechanical property change due to the actions of hydraulic pressure and ice crystallization pressure by considering the coupled effects of expulsion and cryosuction. 3. Methodology For the purpose of simulating the frost damage of cement paste, mainly four aspects are considered as follows: (1) thermodynamics of in-pore phase transition; (2) pore structure of cement paste; (3) mass balance of water; and (4) mechanical equilibrium of frozen cement paste. From the thermodynamics of in-pore phase transition, coupled with the pore structure of cement paste, the freezing of water in cement paste pores can be configured. Based on the mass balance of water and the mechanical equilibrium of frozen cement paste, the response of

Fig. 3. Computational procedure for simulating the internal damage of cement paste subjected to frost action.

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152 Table 1 Microstructure parameters for cement paste generated by HYMOSTRUC3D. Cement type

Portland CEM I 42.5 N

Mineralogical composition of cement Fineness (Blaine surface area value) of cement Minimum diameter of cement particle Particle size distribution

64% C3S, 13% C2S, 8% C3A, 9% C4AF by mass 420 m2/kg 1 μm Rosin–Rammler distribution

Size of cement paste specimen w/c Curing temperature Hydration time Degree of hydration (α) Porosity

F ðxÞ ¼ 1−e−bx with n = 1.698 and b = 0.04408 75 × 75 × 75 μm3 0.4 20 °C 26 days 0.706 0.12 n

cement paste to the freezing of water in pores then can be obtained. In this study, a microstructure model HYMOSTRUC3D (van Breugel, 1991; Ye, 2003) is utilized to obtain a virtual microstructure of cement paste and a 3D lattice fracture model (Qian, 2012; Schlangen, 1993) is employed to analyze the mechanical response of cement paste. The freezing process from 0 to −18 °C is concerned since no phase transition occurs above 0. The computational procedure for simulating the internal damage of cement paste subjected to frost action is illustrated in Fig. 3. A microstructure of cement paste is first generated by HYMOSTRUC3D, as shown in Fig. 3a. Then, from the microstructure of cement paste, its pore structure is extracted and the freezing of water in this pore structure is simulated, see Fig. 3b. Next, a 3-D lattice structure corresponding to the freezing cement paste is constructed which will be used to analyze internal stress distributions, as illustrated in Fig. 3c. Finally, the internal damage of cement paste, in terms of crack pattern and microcrack propagation, is shown in Fig. 3d. 3.1. Freezing of water in cement paste pores The freezing of water in cement paste pores involves the thermodynamics of in-pore crystallization and the related pore structure. According to the freezing order of pores given by Eq. (6), the penetration of ice crystals in a virtual microstructure of cement paste at a given temperature can be obtained. For the cement paste generated by HYMOSTRUC3D, as shown in Fig. 3a, its microstructure parameters are listed in Table 1. The microstructure of cement paste is represented by vectors which give the x, y, and z values of each hydrated cement particle in the 3-D coordinate. For each hydrated cement particle, unhydrated cement and hydration products (i.e., inner and outer CSH) are represented by concentric spherical shells. To guarantee that all the capillary pores involved in the freezing process (in the temperature range of 0 to − 18 °C) can be captured, a resolution of 0.01 μm/voxel (estimated by Eq. (6)) is required if we digitalize the microstructure of cement paste. It is extremely memory

145

expansive to digitalize a microstructure with a dimension of 75 μm into a resolution of 0.01 μm/voxel (7500 × 7500 × 7500 voxels). For the purpose of saving computational memory and improving calculation efficiency, a multi-step digitalization algorithm is proposed (Liu et al., 2011). The pore structure of cement paste is digitalized as pore voxels of various sizes on different pore levels, see Fig. 4. For example, at a resolution of 1 μm/voxel, the microstructure is first digitalized as pore, solid and mixed voxels where the pore and solid phases coexist. Then, the mixed voxels are subdivided into new voxels at a resolution of 0.01 μm/voxel. This multi-step digitalization algorithm was also utilized in calculating the pore size distribution of cement paste and the reliability of pore structures which have been discussed in previous studies (Liu et al., 2012). After the pore structure of cement paste is characterized as multisized pore voxels, the freezing of water in pores is simulated according to the illustration in Fig. 5. The nucleation of ice occurs from a surface of the microstructure, named as freezing front. The penetration of ice crystals takes place through the pore network from the freezing front in. The penetration threshold of ice crystals at temperature T1 and T2 is calculated according to Eq. (6). Parameters associated with the thermodynamics of ice and water are listed in Table 2. Fig. 3b shows the frozen and unfrozen pores at −0.6 °C at a resolution of 1 μm/voxel. The curve of crystal saturation degree SC and frozen porosity versus temperature is plotted in Fig. 6. It should be noted that the freezing of gel pores is not considered while its contribution to the calculation of SC is taken into account. Gel porosity of 0.36 in low density C–S–H indicated in the literature (Jennings et al., 2007) is adopted in this study. By experimental techniques, the amount of ice can be investigated by using AC impedance spectroscopy (ACIS) techniques (Perron and Beaudoin, 2002), a dielectric capacitive apparatus (Fabbri and Fen-Chong, 2013) or a scanning calorimeter (Johannesson, 2010). A comparison of SC by simulation and by experiments is discussed in the previous study (Liu et al., 2014), indicating that the calculation of SC by simulation is reasonable. From the simulation of water freezing in pores, the spatial distribution of ice crystals in the microstructure of cement paste and the ice crystal saturation SC are obtained. The next is to get the response of cement paste to the freezing of water in pores. 3.2. Response of cement paste to the freezing of water in pores In this section, the 3-D lattice fracture model is employed to simulate the crack pattern and microcrack propagation that resulted from internal stresses on pore walls as well as the volume change of cement paste due to frost action. The 3-D lattice fracture model (Qian, 2012) can be used to simulate the stress–strain response, crack pattern and microcrack propagation, based on the microstructure or mesostructure of the material in question. Three stages are defined to make the modeling procedures clear: pre-processing, fracture processes simulation and post-processing, as shown in Fig. 7. In the pre-processing stage, a lattice network is first constructed. As illustrated in Fig. 8, each lattice node is generated in the center of a solid voxel, and all these nodes are connected by lattice

Fig. 4. Schematic diagram of multi-step digitalization algorithm. Reproduced after Liu et al. (2011).

146

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

area of each voxel (A = 1.0 μm2). On the basis of the microstructure of frozen cement paste, considering the relation between PC and PL (Eq. (1)), the mass balance of water (Eqs. (9)–(11)) and the compressibility of ice and liquid water (Eqs. (12)–(13)), the internal stresses on pore walls σ1 and σ2 can be estimated (Eqs. (7)–(8)). The fracture processes are performed by removing the critical element from the system one by one, representing the microcrack occurrence. The critical element is the element with the highest stress/strength ratio and whose tensile stress is larger than its tensile strength. Since microcracks play a role in the expansion of reservoir by providing new space for excess water to move in, the fracture of elements would affect the magnitude and distribution of internal loads in the lattice structure, and consequently, will influence the overall response of cement paste to freezing. The internal loads are modified by taking the volume of microcracks into account in the calculation of σ1 and σ2. Fig. 8b shows how to impose the modified internal loads in the lattice structure with microcracks. The frozen cement paste is a dynamic system containing unfrozen pores, frozen pores, microcraks and solid skeleton. Iterative fracture process simulation which was used to approach the real value of hydraulic pressure (represented by PL,E) in Liu et al. (2014) is employed in this study to approach the real σ1 and σ2. The fracture processes are finished when no critical element is found. In the post-processing stage, the crack pattern and microcrack propagation, the changes in the volume and the elasticity of specimens can be obtained, which will be presented in the following chapters.

a)

b)

4. Changes of microstructure and mechanical property

Unfrozen

c)

Fig. 5. Illustration of the freezing of water in multi-sized pores as temperature decreases.

elements. For cement paste with 75 × 75 × 75 voxel3, its 3D lattice structure is shown in Fig. 3c. The local mechanical properties (Young's modulus E, shear modulus G and tensile strength ft) assigned to every lattice element are determined as the demonstration in Fig. 9 (Liu et al., 2014). The mechanical properties of hydration products and ice crystals are listed in Table 3. Free boundary conditions are set when no external constrain is applied on the freezing sample. In the fracture process simulation, internal loads are imposed in the lattice structure, as shown in Fig. 8. The internal loads are calculated according to the internal stresses on pore walls (σ1 and σ2) and the

4.1. Microstructural changes during freezing For a CEM I 42.5 N cement paste with w/c of 0.4 at 0.71 degree of hydration, Table 4 list the parameters of the cement paste at different nucleation temperatures. At nucleation temperature, a hydraulic pressure is built up in unfrozen pores because of the expulsion of excess water. The maximum hydraulic pressure at nucleation by Eq. (15), given by PL,max, increases with decreasing temperature and increasing SC. Coupled with the crystallization pressure in frozen pores PC, internal damage is caused in cement paste, represented by microcracks. The newly created microcracks provide new space for excess water and release hydraulic pressure, and hydraulic pressure at equilibrium, PL,E, where no more microcracks is created, is obtained by iterative fracture process simulation. At − 18 °C, the maximum hydraulic pressure PL,max equals 96.57 MPa, while the hydraulic pressure at equilibrium PL,E approaches to 2.6 MPa after 12,083 microcracks are created in cement paste. Number of microcracks created in the paste specimen increases with decreasing temperature. The lower the nucleation temperature is, the more serious the internal damage is caused. This is attributed to a higher damage source (hydraulic pressure and crystallization pressure), listed in Table 4. Microcracks created at different nucleation temperatures of −0.3 °C, − 1.8 °C and − 18 °C by simulation are illustrated in Fig. 10. Randomly distributed microcracks are observed in paste specimens. Previous studies (Liu et al., 2011) have indicated that the crack pattern attributed to crystallization pressure in cement paste are randomly distributed microcracks. It can also be found that a transverse crack is formed at a lower nucleation temperature. This is in agreement with the simulation results in Liu et al. (2014) where transverse cracks are observed attributed to hydraulic pressure. Discussion about frost damage by crystallization pressure model and by this study will be presented in Section 5. The influence of hydration time on the frost damage of cement paste by simulation is given in Fig. 11, where microcracks in cement pastes at −18 °C are illustrated. By simulation, for a w/c 0.4 CEM I 42.5 N cement paste, the degrees of hydration are 0.71, 0.85 and 0.97 after 26, 166 and more than 365 days of hydration. The earlier the cement-based materials are subjected to frost action, the more serious the internal damage is caused.

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

147

Table 2 Parameters associated with the thermodynamics of ice and water. Symbols

Values

Property

References

Tm γCL ΔSm δ KL KC

273 K 0.0409 J/m2 1.2 MPa/K 1.0 ~ 1.2 nm 1.79 × 103 MPa 7.81 × 103 MPa

Freezing temperature of bulk ice Energy of the crystal/liquid interface Melting entropy per unit volume of crystal Thickness of the liquid film between crystals and pore wall Bulk modulus of water at 263 K Bulk modulus of ice at 263 K Molar volume of water at 273 K

Coussy and Monteiro (2008) Coussy and Monteiro (2008) (Scherer (1999), Scherer and Valenza (2005)) Sun and Scherer (2010b) Speedy (1987) Lide (2001)

Molar volume of ice at 273 K

Lide (2001)

0

VL

0

0

0

1 − V L / V C ≈ 0.083

VC

At a specific hydration time of 166 days, the simulated internal damage of cement paste with w/c of 0.3, 0.4 and 0.5 subjected to a temperature of − 18 °C is illustrated in Fig. 12. For w/c of 0.3, only randomly distributed microcracks are observed. For w/c of 0.4 and 0.5, both randomly distributed microcracks and transverse cracks are observed. The higher the w/c is, the more serious the internal damage is caused in cement-based materials. The above simulations indicate that the frost damage of cement-based systems is closely associated with its pore structure. For cement-based materials with a more porous structure, the internal damage caused by frost action would be more serious. 4.2. Mechanical property changes due to frost action Ignoring micro-damage caused in the thawing process, changes in mechanical properties of cement-based systems (i.e., elasticity, tensile strength, load–displacement curve) after freezing/thawing can also be obtained by 3D lattice analysis. By applying a prescribed displacement in a direction, the load applied on the specimen can be obtained by the kernel solver in 3D lattice analysis. Fig. 13 illustrates the load– displacement curve of cement paste before and after one freezing/ thawing cycle. From the load–displacement curve in the elastic stage, the Young's modulus of the specimen can be derived. Young's moduli of cement pastes before and after one freezing/thawing cycle are indicated in Fig. 14. Since there are transverse cracks are formed in cement pastes after freezing, as indicated in Figs. 10 to 12, the Young's modulus of cement paste is calculated in the direction vertical to the plane with transverse crack. Therefore, a maximum decrease of Young's modulus is obtained. Damage parameter Df usually utilized to evaluate the extent of damage, can be estimated using Eq. (17), D f ¼ 1−

Ef E0

in the elastic modulus of cement paste due to frost action. For cement pastes with w/c of 0.4, the damage parameter Df decreases with its increasing degree of hydration. For cement pastes with 166 days of hydration, the damage parameter Df increases with its increasing water-to-cement ratio. By testing the dynamic elastic modulus of mortar samples with w/c of 0.42, Li et al. (2012) has given the damage parameter Df as 0.15 after one freeze–thaw cycle. In this study, Df is about 0.10 for cement paste samples with w/c of 0.4 after freezing/ thawing. The discrepancy may be attributed to the differences of samples. In Li et al.'s experiments, the mortar samples is air-entrained and then vacuum water saturated, while in this study, the paste samples are non-air-entrainment. The dynamic elastic modulus tested during freeze–thaw cycles may be different from the Young's modulus simulated in this study. In addition, the internal damage induced in the thawing process is not taken into account in the simulation.

5. Discussion This study has considered the coupled effects of hydraulic pressure and crystallization pressure. Only the hydraulic pressure model and crystallization pressure model are special cases for the frost damage of porous solids. From the aspects of thermodynamic considerations, damage sources, volume changes and crack pattern, discussions about hydraulic pressure model, crystallization pressure model and this model will be described in this section.

ð17Þ

where E0 and Ef represent the Young's moduli of paste specimen before and after freezing/thawing, GPa. The obtained results show a reduction

Fig. 6. Crystal saturation degree SC and frozen porosity versus temperature.

Fig. 7. An overview of the lattice fracture analysis.

148

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152 Table 3 Mechanical properties of hydration products and ice crystals (Petrovic, 2003; Qian, 2012; Qian et al., 2010). Solid phase

Young's modulus (GPa)

Shear modulus (GPa)

Tensile strength (GPa)

Unhydrated cement Inner product Outer product Ice crystals

135 30 22 10.0

52.0 12.0 8.9 5.0

1.80 0.24 0.15 –

Fig. 8. Schematic diagram of constructing lattice structure and imposing internal loads. Table 4 Parameters of cement paste at different nucleation temperatures by simulation.

5.1. Comparison of frost damage by the crystallization pressure model, by the hydraulic pressure model and by this study A comparison of the hydraulic pressure model, the crystallization pressure model and this model is listed in Table 5. To a great extent, the specimen conditions (drainage or non-drainage) and nucleation temperatures determine the damage causes of cement-based materials. The hydraulic pressure model focuses on the expulsion effect, while the crystallization pressure model pays close attention to the interaction of in-pore crystals with pore walls attributed to surface energy effect and cryosuction effect. This model is more realistic to consider the coupled effects of expulsion, surface energy and cryosuction. In order to get the amount of ice, pore size distribution of cement-based materials is considered in the three models. Water mass balance is calculated in the hydraulic pressure model and this model to get the hydraulic pressure PL. That is not considered in the ice crystallization model where drainage condition is assumed. For cement paste with w/c of 0.4 at 0.71 degree of hydration, crack patterns at −0.6 °C and −18 °C simulated by the crystallization pressure model and by this study are illustrated in Fig. 15. According to the crystallization pressure model, random distributed microcracks are created and no transverse crack is formed. In the crystallization pressure model, the stress on pore walls are assumed as disjoining pressure PD, the expulsion of water is ignored. If the contribution of water expulsion is taken into account, local propagation of microcracks would take place, as illustrated in Fig. 15. In this study, the stresses on pore walls are assumed as two parts: one is applied by hydraulic pressure in unfrozen pores (PL N 0), and the other is applied by in-pore crystals of ice equals to the sum of PL and PD (see Eq. (7)).

5.2. Comparison of specimens' deformation during freezing by simulation and by experiment In the process of temperature decreasing, volume changes of cement paste simulated in this study and by the ice crystallization pressure model are given in Fig. 16. By comparison, it is found that the volume

Nucleation temperature

SC

Volume expansion (%)

Number of microcracks

PC (MPa)

PL,max (MPa)

PL,E (MPa)

0 −0.3 °C −0.6 °C −1.0 °C −1.8 °C −7.7 °C −18 °C

0 0.10 0.27 0.33 0.37 0.407 0.41

0 0.012 0.591 0.753 0.867 0.920 1.024

0 2122 8773 9122 9069 11,645 12,083

0 3.48 3.04 4.09 4.86 10.5 23.86

0 17.8 54.9 71.30 83.39 95.54 96.57

0 3.44 2.59 3.19 2.90 1.33 2.60

PL,max, maximum hydraulic pressure at nucleation by Eq. (15). PL,E, hydraulic pressure at equilibrium after nucleation where microcracks and volume change of specimens are taken into account.

expansion attributed to ice crystallization pressure is far less than that by this study considering the coupled effects of hydraulic pressure and ice crystallization pressure. The volume change of specimens during freezing process can be generally assumed as three stages: (I) before nucleation, the volume decreases with the decreasing temperature. A thermal expansion factor of 10.9 × 10−6/°C of cement-based materials is obtained by differential mechanical analyzer (Sun and Scherer, 2010a) and that is adopted in this simulation. (II) At nucleation, a sudden increase in the length of specimen is observed. Sun and Scherer (2010a) have found that a sudden strain increase of the mortar sample is about 1150 × 10−6 at nucleation temperature of − 10.6 °C. In this study, the strain increase for cement paste specimen is 3065 × 10− 6 at nucleation temperature of − 7.7 °C. The deformation of cementbased materials by simulation is larger than that by experiment. This may result from the restrain effect of aggregate which is not considered in simulation. (III) After nucleation, the specimens continuously expand as the temperature decreases. A slope of −34 × 10−6/°C is obtained in the experiments in Sun and Scherer (2010a) and − 35 × 10− 6/°C is obtained by simulation in this study after nucleation. Nevertheless, a slope of −3.5 × 10−6/°C is obtained if assuming that the damage source is only disjoining pressure according to the crystallization pressure model. This indicates that the model established in this study can provide an order of magnitude of volume change consistent with the experiments.

Fig. 9. Determination of the local mechanical properties of an element (Liu et al., 2014). Es, Gs and fs represent the Young's modulus, shear modulus and tensile strength of hydration products. Ei and Gi represent the Young's modulus and shear modulus of ice crystals. φn and φi represent the volume fraction of pores and ice crystals in the voxel.

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

a) T=0.3 ºC, SC=0.10

b) T=1.8 ºC, SC=0.37

149

c) T=18 ºC, SC=0.41

Fig. 10. Microcracks created in cement paste (CEM I 42.5 N, w/c = 0.4, α = 0.71) at nucleation temperatures of −0.3 °C, −1.8 °C and −18 °C.

a) 26 days specimen

b) 166 days specimen

c) 365 days specimen

Fig. 11. Influence of hydration time on the internal damage of cement paste with w/c of 0.4 due to frost action by simulation.

a) w/c=0.3,

b) w/c=0.4

c) w/c=0.5

Fig. 12. Influences of w/c on the internal damage of cement paste with 166 days of hydration due to frost action by simulation.

Before freezing

After freeze/thaw

Fig. 13. Load–displacement curve of cement paste (CEM I 42.5 N, w/c = 0.4, α = 0.71) before and after freezing/thawing.

150

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

a) w/c=0.4

b) 166 days specimen

Fig. 14. Young's moduli of cement pastes before and after freezing.

Table 5 Comparison of hydraulic pressure model, crystallization pressure model and this model.

Specimen conditions Nucleation temperatures Damage causes Stress sources on pore walls Pore size distribution consideration Water mass balance consideration

Hydraulic pressure model

Ice crystallization pressure model

This study

Under undrained conditions Lower than Tm Expulsion effect Hydraulic pressure PL Yes Yes

Under drained condition About at Tm Surface energy effect and cryosuction effect Disjoining pressure PD Yes No

Under undrained conditions Tm ≥ T ≥ −18 °C Coupled effects of expulsion, surface energy and cryosuction Hydraulic pressure PL and the sum of PL and PD (see Eq. (7)) Yes Yes

6. Conclusion This study has established a microstructure-based model to investigate the effects of freezing on the micro-internal damage and macromechanical property changes of cement-based systems. In this model, the coupled effects of water expulsion, in-pore crystallization and cryosuction are taken into account to assess the realistic internal

damage due to frost action. The newly established model is applied to cement pastes with different w/c and degrees of hydration and is validated by experiments from the Young's modulus and the deformation of specimens exposed to frost action. Changes of the microstructure and the mechanical properties of cement paste due to frost action have indicated that cement-base materials with a more porous microstructure would have more serious internal damage. In the discussion

T=0.6ºC, by crystallization pressure model

T=0.6ºC, by this study

T=18ºC, by crystallization pressure model

T=18ºC, by this study

Fig. 15. Cracks pattern simulated by crystallization pressure model and by this study (CEM I 42.5 N, w/c = 0.4, α = 0.71).

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

151

from State Key Laboratory of High Performance Civil Engineering Materials (No. 2012CEM005) and Jiangsu Key Laboratory of Construction Materials in Southeast University (No. 2012CEM01) are also greatly acknowledged. References

Fig. 16. Volume changes of specimens during freezing by simulation and by experiment.

about the hydraulic pressure model, the ice crystallization pressure model and the newly established model, it is found that only randomly distributed microcracks are created according to the ice crystallization pressure. While the propagation of microcracks takes place and a transverse crack is formed according to this model. The internal damage simulated by this model is more serious than those by the ice crystallization pressure model. The model established in this study can provide an order of magnitude of strains consistent with experiments. Due to the limitation of current experimental techniques, the microinternal damage of cement paste during freezing needs to be further validated. This study focused on water-saturated situation, and this model can be extended to consider the effects of air voids during the freezing process in the future. Moreover, the micro-damage caused in the thawing process, which has not been taken into account in this study, will be discussed in future research. From the damaged microstructure of cement-based materials, the reduction in mechanical properties after freezing/thawing has been estimated in this study. In future study, on the basis of the damaged microstructure, the increase in permeability can be evaluated by the use of a transport model, and these can be utilized as inputs for multiscale modeling of concrete properties. Along with the development of the micro-/meso-structure model, the microstructure-based model can also be utilized to other types of concrete. Acknowledgment The financial supports of National Natural Science Foundation of China via Grant No. 51308187, and Natural Science Foundation of Jiangsu Province via Grant No. BK20130837, are greatly acknowledged. The China Postdoctoral Science Foundation (No. 2013M531266), Jiangsu Postdoctoral Science Foundation (No. 1202022C), the foundation

Bazant, Z.P., Chern, J.-C., Rosenberg, A.M., Gaidis, J.M., 1988. Mathematical model for freeze–thaw durability of concrete. J. Am. Ceram. Soc. 71 (9), 776–783. Beaudoin, J.J., MacInnis, C., 1974. The mechanism of frost damage in hardened cement paste. Cem. Concr. Res. 4 (2), 139–147. Coussy, O., 2010. Mechanics and Physics of Porous SolidsJohn Wiley & Sons, Ltd, Publication. Coussy, O., Monteiro, P., 2008. Poroelastic model for concrete exposed to freezing temperatures. Cem. Concr. Res. 38 (1), 40–48. Dai, Q., Ng, K., Liu, Y., Yu, X., 2013. Investigation of internal frost damage in concrete with thermodynamic analysis, microdamage modeling and time-domain reflectometry sensor measurements. J. Mater. Civ. Eng. 25 (9), 1248–1259. Diao, D., Zhang, B., Ye, J., Cheng, S., 2013. Effects of freeze–thaw cycles and seawater corrosion on the behavior of reinforced air-entrained concrete beams with persistent loads. J. Cold Reg. Eng. 27 (1), 44–53. Duan, A., Chen, J., Jin, W., 2013. Numerical simulation of the freezing process of concrete. J. Mater. Civ. Eng. 25 (9), 1317–1325. Fabbri, A., Fen-Chong, T., 2013. Indirect measurement of the ice content curve of partially frozen cement based materials. Cold Reg. Sci. Technol. 90–91, 14–21. Hain, M., Wriggers, P., 2008. Computational homogenization of micro-structural damage due to frost in hardened cement paste. Finite Elem. Anal. Des. 44 (5), 233–244. Jennings, H.M., Thomas, J.J., Gevrenov, J.S., Constantinides, G., Ulm, F.-J., 2007. A multitechnique investigation of the nanoporosity of cement paste. Cem. Concr. Res. 37 (3), 329–336. Johannesson, B., 2010. Dimensional and ice content changes of hardened concrete at different freezing and thawing temperatures. Cem. Concr. Compos. 32 (1), 73–83. Li, W., Pour-Ghaz, M., Castro, J., Weiss, J., 2012. Water absorption and critical degree of saturation relating to freeze–thaw damage in concrete pavement joints. J. Mater. Civ. Eng. 24 (3), 299–307. Lide, D.R., 2001. CRC Handbook of Chemistry and Physics 2001–2002CRC Press. Liu, L., 2012. Modeling the Damaged Microstructure of Cement Paste and Studies on Its Mechanical and Diffusivity Properties(PhD Thesis) Southeast University, Nanjing, China. Liu, M., Wang, Y., 2012. Damage constitutive model of fly ash concrete under freeze–thaw cycles. J. Mater. Civ. Eng. 24 (9), 1165–1174. Liu, L., Ye, G., Schlangen, E., Chen, H.S., Qian, Z.W., Sun, W., van Breugel, K., 2011. Modeling of the internal damage of saturated cement paste due to ice crystallization pressure during freezing. Cem. Concr. Compos. 33 (5), 562–571. Liu, L., Sun, W., Ye, G., Chen, H.S., Qian, Z.W., 2012. Estimation of the ionic diffusivity of virtual cement paste by random walk algorithm. Constr. Build. Mater. 28 (1), 405–413. Liu, L., Chen, H.S., Sun, W., Ye, G., 2013. Microstructure-based modeling of the diffusivity of cement paste with micro-cracks. Constr. Build. Mater. 38, 1107–1116. Liu, L., Shen, D.J., Chen, H.S., Sun, W., Qian, Z.W., Zhao, H.T., Jiang, J.H., 53, 2014. Analysis of damage development in cement paste due to ice nucleation at different temperatures. Cem. Concr. Compos. 53, 1–9. Miao, C., Mu, R., Tian, Q., Sun, W., 2002. Effect of sulfate solution on the frost resistance of concrete with and without steel fiber reinforcement. Cem. Concr. Res. 32 (1), 31–34. Molero, M., Aparicio, S., Al-Assadi, G., Casati, M.J., Hernández, M.G., Anaya, J.J., 2012. Evaluation of freeze–thaw damage in concrete by ultrasonic imaging. NDT&E Int. 52, 86–94. Perron, S., Beaudoin, J.J., 2002. Freezing of water in portland cement paste — an ac impedance spectroscopy study. Cem. Concr. Compos. 24 (5), 467–475. Petrovic, J.J., 2003. Review mechanical properties of ice and snow. J. Mater. Sci. 38 (1), 1–6. Powers, T.C., 1949. The air requirement of frost-resistant concrete. Proc. Highway Res. Board 29, 184–211. Qian, Z., 2012. Multiscale Modeling of Fracture Process in Cementitious Materials (PhD Thesis) Delft University of Technology, Delft. Qian, Z., Schlangen, E., Ye, G., van Breugel, K., 2010. Prediction of mechanical properties of cement paste at microscale. Mater. Constr. 60 (297), 1–12. Scherer, G.W., 1999. Crystallization in pores. Cem. Concr. Res. 29 (8), 1347–1358. Scherer, G.W., Valenza, J.J., 2005. In: Young, F., Skalny, J. (Eds.), Mechanisms of Frost Damage. Materials Science of Concrete, VII. American Ceramic Society. Schlangen, E., 1993. Experimental and Numerical Analysis of Fracture Processes in Concrete(PhD Thesis) Delft University of Technology, Delft. Shang, H., Song, Y., 2013. Triaxial compressive strength of air-entrained concrete after freeze–thaw cycles. Cold Reg. Sci. Technol. 90–91, 33–37. Shi, X., Fay, L., Peterson, M.M., Yang, Z., 2010. Freeze–thaw damage and chemical change of a portland cement concrete in the presence of diluted deicers. Mater. Struct. 43 (7), 933–946. Skripkiūnas, G., Nagrockienė, D., Kerienė, J., Janavičius, E., Girskas, G., Špokauskas, A., 2013. Microstructure changes in hardened cement paste after freezing–thawing cycles. Mater. Sci. 19 (1), 108–114. Speedy, R.J., 1987. Thermodynamic properties of supercooled water at 1 atm. J. Phys. Chem. 91 (12), 3354–3358. Sun, Z., Scherer, G.W., 2010a. Effect of air voids on salt scaling and internal freezing. Cem. Concr. Res. 40 (2), 260–270.

152

L. Liu et al. / Cold Regions Science and Technology 106–107 (2014) 141–152

Sun, Z., Scherer, G.W., 2010b. Pore size and shape in mortar by thermoporometry. Cem. Concr. Res. 40 (5), 740–751. van Breugel, K., 1991. Simulation of Hydration and Formation of Structure in Hardening Cement-based Materials(PhD Thesis) Delft University of Technology, Delft. Wardeh, G., Perrin, B., 2008. Numerical modelling of the behavior of consolidated porous media exposed to frost action. Constr. Build. Mater. 22 (4), 600–608. Wardeh, G., Mohamed, A.S.M., Ghorbel, E., 2010. Analysis of concrete internal deterioration due to frost action. J. Build. Phys. 35 (1), 54–82. Yang, Y., Cao, K., Wang, T., 2013. Experimental behavior of CFST stub columns after being exposed to freezing and thawing. Cold Reg. Sci. Technol. 89, 7–21. Ye, G., 2003. Experimental Study and Numerical Simulation of the Development of the Microstructure and Permeability of Cementitious Materials(PhD Thesis) Delft University of Technology, Delft.

Zeng, Q., 2012. Poromechanical Behavior of Cement-based Materials Subjected to Freeze– Thaw Actions With Salts: Modeling and Experiments(PhD Thesis) Université ParisEst, Paris. Zeng, Q., Li, K., Fen-chong, T., Dangla, P., 2010. A study of the behaviors of cement-based materials subjected to freezing. International Conferences on Mechanic Automation and Control Engineering (MACE), Wuhan, 26–28 June, 2010. Zeng, Q., Fen-chong, T., Dangla, P., Li, K., 2011. A study of freezing behavior of cementitious materials by poromechanical approach. Int. J. Solids Struct. 48 (22–23), 3267–3273. Zuber, B., Marchand, J., 2000. Modeling the deterioration of hydrated cement systems exposed to frost action: part I: description of the mathematical model. Cem. Concr. Res. 30 (12), 1929–1939.