# A numerical microscopically informed upscale approach for analyzing the reliability of testing method for concrete resistance to freeze-thaw

Sara Al Haj Sleiman, Abderrahmane Rhardane, Syed Yasir Alam, Frédéric Grondin, Ahmed Loukili, Laurent Izoret

► To cite this version:

Sara Al Haj Sleiman, Abderrahmane Rhardane, Syed Yasir Alam, Frédéric Grondin, Ahmed Loukili, et al.. A numerical microscopically informed upscale approach for analyzing the reliability of testing method for concrete resistance to freeze-thaw. Construction and Building Materials, 2022, 317, pp.125772. 10.1016/j.conbuildmat.2021.125772. hal-04235565

HAL Id: hal-04235565

https://hal.science/hal-04235565v1

Submitted on 7 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A numerical microscopically informed upscale approach for analyzing the reliability of testing method for concrete resistance to freeze-thaw

Sara Al Haj Sleiman¹,², Abderrahmane Rhardane³,⁴, Syed Yasir Alam¹, Frédéric Grondin*¹,
Ahmed Loukili¹, Laurent Izoret²

¹Institut de Recherche en Génie Civil et Mécanique (GeM), UMR 6183, Centrale Nantes – Université de Nantes – CNRS, 1 rue de la Noë, BP 92101, 44321 Nantes cedex 3, France

²Syndicat Français de l'Industrie Cimentière (SFIC), 7, place de la défense
92974 Paris-la-Défense Cedex, France

³CNES, 2 place Maurice-Quentin, Paris, France

⁴Institut PPRIME, UMR CNRS 3346, Ecole Nationale Supérieure de Mécanique et Aéronautique,
2 Boulevard des Frères Lumière, 86360 Chasseneuil-du-Poitou, FRANCE

*frederic.grondin@ec-nantes.fr; tel: +33240371668 (Corresponding author)

## ABSTRACT:

Due to the complexity of the problem, experimental studies on the frost action in concrete cannot be validated directly, since the major challenge of such a complex phenomenon is to be able to isolate the influence of each parameter involved in freezing and thawing processes. Therefore, the development of numerical tools offers an opportunity to represent the complexity of the material under study and to investigate as closely as possible the local phenomena at small scales. In this work, we adopt a micromechanical approach aiming to integrate the effect of each parameter defining laboratory exposure conditions of concrete. A multiscale approach originating from the hydration products and goes up to concrete was presented and results are compared with the experimental measurements. The results show that this numerical approach offers a satisfactory tool to predict concrete strains when exposed to freeze and thaw damage, without any calibration, by taking into account the effects of the developed pressures in pores at the microscale. An application of the model is presented to predict the strains of concrete specimens exposed to freeze and thaw cycles placed at different positions inside a climatic chamber. The temperature was recorded at different specimens' positions inside the climatic chamber and these temperatures were applied as thermal loadings in the model. A clear correlation between the experimental degraded surfaces and the strain amplitudes recorded for each specimen is found. The study also found an inconsistency in the recorded temperature evolutions and deformation levels among the concrete specimens placed at different positions inside the climatic chamber, which calls into question the reliability of the existing normative methods in order to characterize concrete resistance to freeze and thaw.

**Keywords:** concrete, freeze-thaw, deterioration, durability, multi-scale modeling, micromechanical approach, thermo-poro-mechanical modeling.
# 1. Introduction

Concrete infrastructures, notably highways, bridges and dams, located in cold climates suffer from the effect of frost action during winter, where under extreme cold environments salts are also regularly used to de-ice surfaces and roads. Under these conditions, a superficial scaling of concrete can result from the freezing of the saline solution on the surface and an internal cracking caused by cyclic freezing/thawing cycles resulting in a reduction of strength may occurs. However, these two types of degradation may not occur at the same time. According to [1], scaling damage is observed without a significant loss of strength. It is dependent on the existence of a frozen saline solution on the surface. This can be seen in experimental measurements when no scaling is identified if the saline solution is evaporated from the surface during freeze-thaw tests. The main parameters influencing the deterioration of concrete under freeze and thaw laboratory cycles are the minimum temperature, the salt and moisture content, the freezing rate, and the time duration at minimum temperature [2]. Under the action of freeze and thaw cycles, water phase transformations and thermal dilation of solid phases within the heterogeneous and porous cement paste contribute to the gradual degradation of concrete, as these physical processes induce pressure build-up in the pore network.

There are numerous theories as well as models that can describe the effect of freeze-thaw cycles on damage in cementitious materials. None of them fully explain the problem, because of the coupling of different thermo-hydro-chemo-mechanical phenomena. But each one covers at least one area of this vast problem, so that by their complementarity, certain phenomena can be interpreted. Indeed, due to the complexity of the problem, experimental studies on frost action in concrete cannot be validated directly since the major challenge of such a complex phenomenon is to be able to isolate and understand the influence of each parameter.

The phase change of water into ice is accompanied with a volume expansion of 9%. Due to this transformation, water is expelled from frozen to unfrozen pores; a hydraulic pressure is therefore developed [3]. Furthermore, when ice forms, it releases the ions dissolved in the interstitial solution, originating from the de-icing salts. Therefore a gradient of concentration is formed that causes the development of osmotic pressure [4]. The third type of pressures developed is the crystallization pressure. This is the result of the repulsive forces between the formed ice crystal and the pore walls, these pressures causes tension that can reach 10 MPa [5].

Under cold exposure, as a consequence of the cited mechanisms and others, the exposed concrete may undergo two types of degradation: internal damage and surface scaling. Internal damage displays as an important micro-cracking that degrades the coupling of the cement paste with the aggregates. Pigeon and Regourd [6] showed that since the calcium hydroxide and the sulfoaluminates are dissolved under low temperature, their recrystallization in the form of fibrous secondary hydrates in the air voids weakens the concrete by increasing its capillary porosity. It results a loss in the tensile strength and in the Young's modulus as well as an increase in the permeability of the material to aggressive agents such as chloride ions. [6]–[8].
As its name suggests, surface scaling is a mode of degradation that only affects concrete surfaces exposed to an external environment with freeze-thaw cycles and de-icing salts. Scaling is manifested by the progressive detachment of small fragments or flakes of paste or mortar whose thickness does not exceed few millimeters. Consequently, not only an aesthetic problem is caused to concrete under scaling effect, but also a loss of concrete cover as well as an increase of permeability and porosity of the surface which can have great influence on the durability of the structure when exposed to aggressive agents. This type of degradation might be related to the local variability in the superficial concrete properties such as the air void system and its connectivity, surface finishing, and mostly depending on the exposure conditions as well such as factors related to the de-icing salts concentration, water saturation level, etc. [1], [9]. Even if scaling is measured by the mass of separated cement paste particles, after damage, it is difficult to define its origin. Due to the multiple causes of degradation, the development of numerical tools offers an opportunity to represent the complexity of the material under study and to investigate as closely as possible the local phenomena at small scales. In this study, a damage model is proposed to first of all link the thermo-chemical effects to damage. At this step it is too early to simulate the scaling which needs more statistical studies on the relationships between the local physical phenomena and the local separation of material.

Quantitative modeling of freeze-thaw mechanisms is relatively recent due to the latest theoretical developments and the development of sophisticated numerical tools. Based on the work elaborated by Coussy [10]–[12], poromechanical models are by far the most widely used. However, the classical poromechanical approach considers a homogeneous medium equivalent to the material studied and does not allow studying the effect of heterogeneities. Moreover, it seems interesting to make an explicit representation of the porous network by using virtual microstructures of the hydrated cement paste to model explicitly the freeze-thaw mechanisms in the porosity. Liu et al [13] have used HYMOSTRUC3D to simulate internal damage in a virtual hydrated cement paste. Internal cracks could be visualized as well as the ice propagation process in the porous network. However, an overestimation of the deformation and the ice saturation was observed. Several models are available for modeling the degradation of cement-based materials exposed to freeze-thaw, at meso or macroscopic scales. Each model is based on a different method, such as the discrete element method used in the work carried out by Si et al. [14] to simulate internal cracking, and the fractal method used by Jin et al. [15], for a quantitative analysis of the concrete pore size distribution. A damage model based on the variation of the pore structure on the microscale of concrete is then used to determine the freeze-thaw resistance of concrete. Another constitutive model formulated within the framework of continuum damage mechanics is presented in [16]. However these damage models do not take into account the different mechanical behaviour between the various hydrates and cement paste phases. Moreover, the different pressures induced by freeze action must be taken explicitly into account in a model that allows studying their effect locally at cement paste scale. On the other hand, an in-depth study of the parameters involved in the modelling of these mechanisms is absent in the literature and, to our knowledge, no study has exploited this type of modelling to study the effect of the saturation degree, the minimum temperature, the water-to-cement (w/c) ratio and the effect of supercooling on the deformation of the paste. These parameters define the concrete exposure to freeze-thaw
conditions. Therefore, this study proposes a multi-scale model in which the mechanisms attributed to freeze-thaw are explicitly modeled and applied on a virtual microstructure through a micromechanical approach developed in [17]. A specific attention has been given to the cross-effects between the salt presence and the water saturation degree on the frost degradation of cementitious materials. The first part of this study on the effect of each frost parameter only at the cement paste scale has already been done by the authors and details are given in a previous article by the same authors [18]. This paper, however, presents the development of a comprehensive model to consider all types of pressures observed during freezing, the thermal dilation of all phases, in particular that of water, and a mechanical damage model. The up-scaling was performed considering three scales: the cement paste scale, called the microscopic scale; the mortar scale, called the mesoscopic scale; the concrete scale, called the macroscopic scale. At the macroscopic scale, concrete is assumed to be constituted of two homogenous components: a solid matrix (mortar) surrounding the inclusions (aggregates). The same assumption is considered at the mesoscopic scale with a solid matrix (representing homogeneous cement paste) surrounding inclusions that are sand grains. At these two scales, thermo-mechanical behaviour laws are relevant to calculate strain and damage of the material under frost action [19]. However, since the input data at these two scales is related to the chemical changes and pressure effects which occur at the cement paste scale, an explicit representation of pores was necessary at the microscopic scale to distinguish all types of pressures, the mechanical behaviour of the different phases as well as their chemical reactions. In this study, an upgrade of the micromechanical approach was done aiming to model the strain at concrete and mortar scales under freeze-thaw cycles.

In the following, some remarks on different measurement setups and discrepancies in the previous studies are first provided, which motivates the present numerical study. Then the main mechanisms related to freezing as well as the micro-mechanical approach are recalled from the previous paper [18] to model the thermo-poro-mechanical behavior of cement paste under freeze-thaw cycles. The multi-scale approach is then given in details to simulate the behaviour of concrete based on simulations at the cement paste scale. The numerical results are compared to experimental measurements. These results show that numerical strain curves are in good agreement with strain obtained from experimental measurements. Also, an application of the developed model is presented to predict the strain amplitudes that concrete specimens will develop under freeze-thaw cycles placed at different positions inside a climate chamber. Here, a good correlation was also found between the levels of surface degradation detected in the test for the different specimens and the strain amplitudes given by the simulation results.

## 2. Some remarks on previous studies

In order to test the validity of the model, comparisons have been made with experiments on freeze-thaw of concrete and mortar found in literature. The expansion was measured using strain gauges in some researches [20], [21], [22], by means of length change measurement tools such as LVDT [23], [24] or using differential mechanical analyzers (DMA) as in [25]. Peak strains at -20°C of different literature studies are reported in
Table 1. An important discrepancy is noticed between the different peak strain values that may be first related to experimental arrangements that are different. Using strain gauges to measure displacement is one of the most commonly used methods. However, the results may not be representative of the overall expansion of the material and might not be very precise. This is because the gauge gives local strain measurements, and therefore cannot represent global concrete strain since the deformation of this latter happens in all directions under freeze and thaw cycles. Moreover, the gauge could be placed parallel to an aggregate which can significantly alter the results. Using LVDT might be a better choice for length change measurements.

<table><caption>Table 1. Peak strains from different works in the literature.</caption><thead><tr><th>Reference</th><th>Material</th><th>Size (mm)</th><th>Peak strain (μm/m)</th></tr></thead><tbody><tr><td>Johannesson [24]</td><td>Concrete<br>D<sub>max</sub> = 8mm</td><td>20×20×160</td><td>650</td></tr><tr><td>Sun and Scherer [25]</td><td>Mortar</td><td>5×5×15</td><td>900</td></tr><tr><td>Kaufmann [27]</td><td>Concrete<br>D<sub>max</sub> = 32 mm</td><td>200×200×50</td><td>450</td></tr><tr><td rowspan="2">Gong et al. [22]</td><td rowspan="2">Mortar</td><td rowspan="2">40×40×2</td><td>750</td></tr><tr><td>1000</td></tr></tbody></table>

Other reasons behind the discrepancy of the experimental strain measurements might be related to material properties and differences in composition (cf. [26], [22], [21]). In [22] for example, the maximum strain amplitude reached after one freeze-thaw cycle was about 1000 µm/m for a W/C =0.5 mortar with 1090 kg/m³ fine aggregates, and about 750 µm/m for a W/C =0.5 mortar with 754 kg/m³ fine aggregates. As for experiments run on concrete specimens, the peak strain at -20°C is equal to 450 µm/m in [27] for a concrete with maximum aggregate diameter D<sub>max</sub>= 32 mm and 750 µm/m in [24] measured on concrete prisms with D<sub>max</sub>=8 mm. More values are reported in [22].

Other reasons behind these large scattered experimental values might be related to test conditions: e.g. age of specimen during the test, storing and curing conditions as well as the risk of carbonation that concrete may undergo during the preconditioning phase which can greatly affect the results by modifying the composition of the exposed testing face. Furthermore, experimental devices and procedures used for measuring strain can greatly affect the results. The measured strains depend strongly on the experimental conditions and the physical problem that was considered by each researcher. In Kaufmann's experiments for example [100], [27], the displacement of concrete specimens was measured between the centers of two steel plates fixed to the specimen by means of two inductive transducers. In the experiment one steel plate was kept fixed to the invar frame, and the corresponding displacement transducer was used for control purpose only. The opposite concrete side was allowed to move freely so that the expansion of the whole plate was actually measured by one displacement transducer only [27]. However in [22], strain gauges were glued on different sides of the specimen. [21] showed that there is strain gradient along the height of the specimen due to the differences in
moisture at different heights; strain measured in the middle part of the prism was quite higher than the upper part. Therefore, the strain gauge measurement is a local measurement that depends on the differences in moisture along the height of the specimen and therefore, cannot be considered as a relevant parameter to evaluate the overall deformation of the specimen.

## 3. Problem formulation

### 3.1. Coupled nonlinear thermo-poro-mechanical model at the cement paste scale

At the cement paste scale, we consider a representative elementary volume $V_{REV}$ consisting of the volumes of the pore network $V_p$ and the one of the solid cementitious matrix $V_{CM}$. $V_{CM}$ is the binding of $V_{CSH}$ (the volume of the calcium silicate hydrates (C-S-H) phase) and $V_{Si}$ which is the volume of the other separated constituents of the cement matrix. The solution of the below described thermo-mechanical problem (eqs.(1)-(7)) gives the local displacement field ($\bar{y}$):

$$ \bar{\nabla} \cdot \bar{\bar{\sigma}}_I(\bar{y}) = 0 \qquad \bar{y} \in V_{REV} \quad (1) $$

$$ \bar{\bar{\epsilon}}_I(\bar{y}) = \frac{1}{2} [\bar{\nabla} u_I(\bar{y}) + t \bar{\nabla} u_I(\bar{y})] \qquad \bar{y} \in V_{REV} \quad (2) $$

$$ \bar{\bar{\sigma}}_I(\bar{y}) = \mathbb{C}^d_I(\bar{y}, \bar{\epsilon}(\bar{y})) : [\bar{\epsilon}_I(\bar{y}) - \bar{\alpha}_{Si}(\bar{y}, T) \Delta T] \qquad \bar{y} \in V_{Si} \quad (3) $$

$$ \bar{\bar{\sigma}}_I(\bar{y}) = \mathbb{C}^d_I(\bar{y}, \bar{\epsilon}(\bar{y})) : [\bar{\epsilon}_I(\bar{y}) - \bar{\alpha}_{CSH}(\bar{y}, T) \Delta T - p_{osm}(T, S_C) \bar{I}] \qquad \bar{y} \in V_{CSH} \quad (4) $$

$$ \bar{\bar{\sigma}}_I(\bar{y}) = \mathbb{C}^0_I(\bar{y}, \bar{\epsilon}(\bar{y})) : \{\bar{\epsilon}_I(\bar{y}) - \bar{\alpha}_P(\bar{y}, T) \Delta T - [\mathcal{H}(S-1)p_{hyd} + \mathcal{H}(1-S)p_{cap}] \bar{I}\} \qquad \bar{y} \in V_P \quad (5) $$

$$ \rho_I(\bar{y}, T) c_{pI}(\bar{y}, T) \frac{\partial T}{\partial t}(\bar{y}, t) + \nabla q_I(\bar{y}, T) = f(\bar{y}, t) \qquad \bar{y} \in V_{REV} \quad (6) $$

$$ q_I(\bar{y}, T) = -\lambda_I(\bar{y}, T) \cdot \nabla T(\bar{y}, t) \qquad \bar{y} \in V_{REV} \quad (7) $$

$$ \bar{T}(\bar{y}) = \bar{T}_{imp} \qquad \bar{y} \in \partial V \quad (8) $$

$$ \bar{u}_I(\bar{y}) = \bar{U}_{imp,I} \qquad \bar{y} \in \Gamma_u \quad (9) $$

$$ \| \bar{u}_I(\bar{y}) \| = \bar{0} \qquad \bar{y} \in \Gamma_i \quad (10) $$

$$ \| \bar{\bar{\sigma}}_I(\bar{y}) \| = \bar{0} \qquad \bar{y} \in \Gamma_i \quad (11) $$

where $\bar{\bar{\alpha}}_{Si}$, $\bar{\bar{\alpha}}_{CSH}$ and $\bar{\bar{\alpha}}_P$ are the thermal expansion coefficient (TEC) of the cement solid phases $i$ (without C-S-H), C-S-H, and the porosity, respectively. $\Gamma_i$ represents the interface between two different components and the brackets $||$. || the continuity condition of the field. The index “I” used in the above equations refers to the cement paste scale. The capillary pressure is activated when the pore is not totally saturated with water ($S<1$)
and in case S=1 the hydrostatic pressure is activated. This is done by introducing the Heaviside function $\mathcal{H}$. The boundary conditions (8) and (9) describe respectively the imposed temperature and the restrained displacement in the middle of the volume. Damaged and undamaged elastic stiffness tensors $\mathcal{C}^d$ and $\mathcal{C}^0$ are linked by the scalar damage value D (considering isotropic materials) as follows:

$$ \mathcal{C}^d = (\mathbb{I} - D\mathbb{I})\mathcal{C}^0 \qquad (12) $$

The loading surface f is given by following where only the tensile eigenvalues of the strain tensor are considered:

$$ f(\bar{n}) = \epsilon_{EQ} - \kappa_0 = \sqrt{\langle\epsilon_{11}\rangle_+^2 + \langle\epsilon_{22}\rangle_+^2 + \langle\epsilon_{33}\rangle_+^2 - \kappa_0} \qquad (13) $$

Where $\kappa_0$ represents the strain damage threshold, $\epsilon_{EQ}$ the equivalent strain and $\langle\cdot\rangle_+$ the Macaulay brackets.

The following exponential form gives the damage law evolution due to tensile microcracking in an isotropic model:

$$ \begin{cases} f = 0 \text{ and } \dot{\epsilon}_{EQ} > 0 \rightarrow \dot{D} = \frac{\kappa_0(1 + B_T\epsilon_{EQ})}{\epsilon_{EQ}^2}\exp[-B_T(\epsilon_{EQ} - \kappa_0)] \\ f < 0 \text{ or } \dot{\epsilon}_{EQ} = 0 \rightarrow \dot{D} = 0 \end{cases} \qquad (14) $$

with $B_T$ the strain softening slope parameter. Using the below equations, $\kappa_0$ and $B_T$ can be determined from mechanical properties of cement paste that can be measured:

$$ \kappa_0 = \frac{f_t}{E} \qquad B_T = \frac{f_t}{\frac{G_f}{h} - \frac{f_t \kappa_0}{2}} \qquad (15) $$

where E, $f_t$ and $G_f$ represent respectively the Young's modulus, tensile strength and fracture energy of the cement paste's phase. By introducing the characteristic size h, the energy dissipation in the finite element during damage becomes independent of the mesh, this can be seen in the equation of $B_T$.

Eq. (6) describes the thermal problem formulation, where transient heat conduction is considered, for any point $\bar{y}$ inside the volume of the material, at any scale. Note that $c_p$ is the specific heat capacity, t the time, q the heat flux given by eq. (7), f the heat source and $\lambda$ the heat conductivity coefficient. We consider a stationary fluid inside the porosity of the cement paste. The heat transport mechanisms in porous materials such as concrete are complex. Thermal radiation between particles can be involved, as well as convective exchange inside the porosity. However, simulations performed in the representative elementary volume consider it in the middle of the structure far away from the boundary where radiative effect can occurs. Also, internal convective transfer is supposed negligible due to the low size of pores.

The heat flow problem is coupled with the mechanical problem to have explicit interactions between phenomena. At the upper scales, the heat flow is considered to calculate the thermo-mechanical coupling. And the pressure effects are implicitly taken into account via numerical homogenization.
## 3.2. Calculation of the pore pressures at the cement paste scale

Under normal conditions of temperature and pressure (atmospheric pressure $p_0 = 1atm$, fusion temperature $T_0 = 273.15K$ at $p_0$), the thermodynamic equilibrium of the porosity is established when equation (16) is satisfied.

$$ \mu_L(T, p_L, a_L) = \mu_C(T, p_C) = \mu_V(T, p_V) \qquad (16) $$

The liquid, ice crystal and vapor phases are designated by the subscripts L, C and V respectively. We can write the chemical potential of water vapor in function of the relative humidity RH as following:

$$ \mu_V(T) = RT\ln(RH) \qquad (17) $$

where RH takes the following pore size-dependent form [28]:

$$ RH(r, T) = \exp\left(-\frac{2\gamma\cos\theta}{r - \delta}\frac{v_L}{RT}\right) \qquad (18) $$

where r represents the pore radius, θ the liquid water angle (the contact angle of the water surface with the pore walls, assumed to be zero), γ the interface energy between vapor and liquid and δ the thickness of the film of adsorbed water, also called the liquid-like layer (LLL), near the pore boundary, given by [28], [29]:

$$ \delta(T) = a + bRH(T) + \frac{c}{1.105 - RH(T)} \qquad (19) $$

with $a = 0.0626nm$, $b = 0.834nm$ and $c = 0.02309nm$.

Under -25°C, the freezing process is significantly affected by the LLL[30]. In fact, in addition to the water content, the amount of the frozen water in a cement paste is determined by the structure of the pore network (size and shape distributions). Moreover, an interface with a curvature $\kappa_{LC}$ maintains the equilibrium between the frozen and unfrozen sites. The curvature of this interface depends on the capillary pressure expressed as a pressure gradient $p_{cap} = p_C - p_L$. We can write this gradient in two forms:

$$ p_C - p_L = S_f(T_f - T) \qquad (20) $$

$$ p_C - p_L = \gamma_{LC} \frac{dA_{LC}}{dV_L} = \frac{\gamma_{LC}}{r_H} \qquad (21) $$

The first one is linked to the change in temperature, with $S_f \approx 1.2 MPa/K$ representing the molar fusion entropy per unit volume of crystal [31] and $T_f$ the fusion temperature.

The second one is the Laplace-Young equation and corresponds to the local equilibrium at the interface [31], [32], where $A_{LC}$ is the liquid-crystal interface area, $\gamma_{LC} = 40.9 \times 10^{-3} J/m^2$ the liquid-crystal interface energy [33] and $r_H$ the hydraulic radius of the interface.

These two relations lead to the equality:

$$ r_H = \frac{\gamma_{LC}}{S_f(T_f - T)} \qquad (22) $$
$r_H$ can be also linked to $r_C$ (known as the ice front) via a shape factor $\lambda$ (in case of an hemispherical ice front $\lambda = 0.5$, and it is equal to 1 for a spherical one [34]) as follows:

$$r_H = \lambda r_C = \lambda(r - \delta) \qquad (23)$$

Assuming hemispherical pore shapes, equations (22) and (23) give the Gibbs-Thomson equation linking the threshold radius of frozen sites $r = r_f$ to the frozen temperature $T_f$:

$$r_f(T) = -\frac{2\gamma_{LC}}{S_f(T - T_f)} + \delta(T) \qquad (24)$$

The partial volume $V_N(r_f)$ of frozen sites with pore radius greater than $r_f$ and smaller than the maximum capillary pore radius $r_{max}$ (also called crystal saturation) can be determined from the pore size distribution:

$$V_N(r_f) = \int_{r_f}^{r_{max}} \varphi(r)dr \qquad (25)$$

where $\varphi(r)$ represents the pore sizes distribution in cement paste.

The normalized cumulative porosity volume with pore radii bigger than $r_0$ has been calculated using the Brunauer-Skalny-Bodor (BSB) adsorption model [28], [29]:

$$V_N(r \ge r_0) = \frac{(1-k)}{\left[1 - k \cdot \exp\left(-\frac{2\gamma\cos\theta}{r_0 - \delta} \frac{v_L}{RT}\right)\right]} \times \frac{[1 + (C-1)k] \exp\left(-\frac{2\gamma\cos\theta}{r_0 - \delta} \frac{v_L}{RT}\right)}{\left[1 - (C-1)k \cdot \exp\left(-\frac{2\gamma\cos\theta}{r_0 - \delta} \frac{v_L}{RT}\right)\right]} \qquad (26)$$

where $C$ and $k$ are the parameters of the BSB model. These parameters depend on the w/c ratio, the cement type, the temperature and the hydration age of the cement paste [35]. According to that, the pore size distribution function $\varphi(r_0)$ can be deduced as:

$$\varphi(r_0) = \frac{dV_N(r \ge r_0)}{dr} \qquad (27)$$

Consequently, the degree of liquid water saturation $S_L$ can be calculated as follows:

$$S_L(T) = \frac{V_L}{V_C + V_L} = 1 - S_C(T) = 1 - \left[ V_N \left( r_f = -\frac{2\gamma_{LC}}{S_f(T - T_f)} + \delta \right) \right] \qquad (28)$$

where $V_L$ represents the amount of unfrozen water and $V_C$ the amount of crystal and $S_C$ the ice crystal saturation.

The liquid water saturation $S_L$ as well as the ice crystal saturation $S_C$ control the pressure applied on the pore walls. When the volume of the porosity is fulfilled, especially after the 9% increase in ice crystal volume during the water phase change, the excess water exert a hydrostatic pressure on the capillary walls, this causes the cement paste to expand. On the other hand, when a space is still available in the pore, a contracting pressure is developed that contributes to the overall shrinkage of the cement paste. The model takes into account the case of an air entrained concrete. The entrained air voids are connected to capillary pores, therefore, hydrostatic pressure cannot be imposed in these latters since are able to contain more water in their
volumes. This reduces the total expansion of the cement paste. An algorithm of the numerical calculations of the pressures is presented in authors' other publication [18].

The other internal pressure is the osmotic pressure that is developed in the presence of ion concentration differences. When water freezes, it rejects the ions it originally contained, thus locally increasing the concentration of unfrozen water [9]. However, the free energy of the solvent in solution is lower than that of pure solvent. Thus, the solvent passes from the medium where its free energy is high to the medium where its free energy is low and thus to the solvent in solution. The unfrozen water in the pores of the gel is therefore more energetic than the water frozen in the capillary pores. The concentration gradient creates a diffusive transport of the water molecules. To this mechanism is added the existence of a chemical potential gradient between, on the one hand, the water in equilibrium with the ice in the large pores affected by the frost front and, on the other hand, the water in the small pores (gel pores). Thermodynamic equilibrium is only restored by the migration of water from the gel pores to the gel formation site through semipermeable canals of nanopores of C-S-H. In this case, the chemical potentials are written at the interface between the capillary pores and the C-S-H nanoporosity as follows:

$$ \left\{ \begin{array}{l} \mu_L^{cap}(T, p_L, a_L^{cap}) = \mu_L^0(T, p_L) + RT\ln(a_L^{cap}) \\ \mu_L^{CSH}(T, p_L + p_{osm}, a_L) \approx \mu_L^0(T, p_L) + v_L p_{osm} + RT\ln(a_L^{CSH}) \end{array} \right. \qquad (29) $$

where R is the gas constant, p the pressure, a_L the water activity and p_{osm} the osmotic pressure.

Under the assumption of a thermodynamic equilibrium, chemical potentials are equal, therefore the osmotic pressure can be written as following:

$$ p_{osm} = \frac{RT}{v_L} \ln \left( \frac{a_L^{cap}}{a_L^{CSH}} \right) \qquad (30) $$

In C-S-H gel and capillary pores, the activity of a solvent (here liquid water) in the presence of an electrolyte is given according to an osmotic coefficient Φ:

$$ \ln(a_L) = \Phi b_T M_L \qquad (31) $$

where $M_L = 18.015g/mol$ is the molar mass of water and $b_T$ the total molality of ions:

$$ b_T = \sum_i b_i = \sum_i \frac{n_i}{m_L} \qquad (32) $$

where $n_i$ is the amount of substance of the ion i and $m_L$ the mass of the solvent. In the case of a 1:1 type electrolyte such as NaCl, we can determine the molality from the mass fraction $x_{salt}$ of the salt:

$$ b_i = M_{salt} \frac{x_{salt}}{1 - x_{salt}} \qquad (33) $$

The osmotic coefficient is given by the semi-empirical model of Lin and Lee [36]:

$$ \Phi = 1 - A_{\Phi} \left( \sum_i x_i z_i^2 \right) \frac{\sqrt{I}}{1 + B_i \sqrt{I}} + \frac{\alpha I^{\alpha}}{(\alpha + 1)I} \left( \sum_i C_i x_i z_i^2 \right) \qquad (34) $$
where $x_i$ and $z_i$ are the molar fraction and the absolute charge of the ion $i$ (for NaCl, $x_{Na^+} = x_{Cl^-} = 0.5$ and $z_{Na^+} = z_{Cl^-} = +1$). The parameters $\alpha$, $B_i$ and $C_i$ depend on the nature of the electrolyte, [36] gives the values of these parameters for several types of electrolytes. $A_\Phi$ is the Debye-Hückel constant (full expression is given in [36]) and $I$ is the ionic strength given as follow for a non-ideal solution as:

$$I = \frac{1}{2} \sum_{i} b_{i} z_{i}^{2} \qquad (35)$$

### 3.3. Coupled non-linear thermo-mechanical model at the mortar and concrete scales

Moving to the mortar and concrete scales, the cementitious material can be simulated using the thermo-mechanical behavior law while considering homogenous non-porous components:

$$\bar{\nabla} \cdot \bar{\sigma}_{II}(\bar{y}) = 0 \qquad \bar{y} \in V_{REV} \qquad (36)$$

$$\bar{\epsilon}_{II}(\bar{y}) = \frac{1}{2} [\bar{\nabla} u_{II}(\bar{y}) + t \bar{\nabla} u_{II}(\bar{y})] \qquad \bar{y} \in V_{REV} \qquad (37)$$

$$\bar{\sigma}_{II}(\bar{y}) = C^{d}_{II}(\bar{y}, \bar{\epsilon}(\bar{y})) : [\bar{\epsilon}_{II}(\bar{y}) - \bar{\alpha}_{i}(\bar{y}, T) \Delta T] \qquad \bar{y} \in V_{i} \qquad (38)$$

$$\rho_{II}(\bar{y}, T) c_{pII}(\bar{y}, T) \frac{\partial T}{\partial t}(\bar{y}, t) + \nabla q_{II}(\bar{y}, T) = f(\bar{y}, t) \qquad \bar{y} \in V_{REV} \qquad (39)$$

$$q_{II}(\bar{y}, T) = -\lambda_{II}(\bar{y}, T) \cdot \nabla T(\bar{y}, t) \qquad \bar{y} \in V_{REV} \qquad (40)$$

$$\bar{T}(\bar{y}) = \bar{T}_{imp} \qquad \bar{y} \in \partial V \qquad (41)$$

$$[\bar{\sigma}_{II}(\bar{y}), \bar{n}]. \bar{n} = 0 \qquad \bar{y} \in \partial V \qquad (42)$$

$$\bar{u}_{II}(\bar{y}) = \bar{U}_{imp,II} \qquad \bar{y} \in \Gamma_{u} \qquad (43)$$

$$\| \bar{u}_{II}(\bar{y}) \| = \bar{0} \qquad \bar{y} \in \Gamma_{i} \qquad (44)$$

$$\| \bar{\sigma}_{II}(\bar{y}) \| = \bar{0} \qquad \bar{y} \in \Gamma_{i} \qquad (45)$$

The index II used in the above equations refers to mortar and concrete scales. The thermal problem is given in equations (39)-(40). Equation (38) considers damage as formulated in the section 3.1. $V_{REV} = V_g \cup V_m$, where $V_g$ is the aggregates volume and $V_m$ is the volume of the matrix (homogeneous paste/mortar). When talking about aggregates (sand or gravel) $V_i=V_g$ otherwise $V_i=V_m$. The thermal dilation tensors of sand grains and aggregates, in mortar and concrete respectively, are found in the literature. That of the matrix is obtained from the calculation at the lower scale. The thermal dilation of the mortar's matrix is deduced from the calculation at the cement paste scale, and that of the concrete's matrix is then obtained from the calculation at the mortar scale. Because different values of strain are obtained at one temperature reached at each cycles (due to the freeze/thaw hysteresis), a function for the thermal dilation has been suggested, deduced from the total strain equation:

$$\bar{\alpha}'_{i} = \frac{\bar{\alpha}_{0} T_{0} + \bar{\epsilon}_{i}}{\bar{T}_{i}} \qquad (46)$$
where $\overline{\alpha}_0$ represents the initial thermal dilation at the initial temperature $T_0$, $\overline{\epsilon}_i$ tensor is the total strain calculated at the temperature $\overline{T}_i$.

The elastic properties of the matrix are obtained from the calculation at the lower scale as follows: the homogenized Young modulus and the Poisson's ratio are calculated from the slope of the linear part of the curve ($\sigma-\epsilon$) obtained by simulating the mechanical response of the heterogeneous matrix under compressive (or tensile) loading at the lower scale.

The thermal problem formulation is the same as described at micro-scale, as seen in equations (39) and (40). However, the thermal properties of the matrix ($\rho_{II}$, $c_{pII}$ and $\lambda_{II}$) are calculated by homogenization methods as follows:

- At mortar scale, the homogenized density of the cement paste is calculated by $\rho_{(hom)} = \sum_{i=1}^{56} \rho_i \cdot f_i$ where $\rho_i$ [kg/m³] is the density of the cement paste's phase i ($i=1... 56$, where 56 is the total number of phases in the microstructure generated by VCCTL (cf. section 4.1): hydrates and unhydrated cement phases) and $f_i$ the volume fraction of each phase. The same equation is used to calculate the homogenized properties of the mortar matrix considered at the concrete scale, considering only two components: cement paste and sand grains.

- The homogenized heat capacity is obtained from the following equation:

$$C_P(\text{hom}) = \frac{\sum_{i=1}^{56} \rho_i \cdot A_i \cdot f_i}{\rho(\text{hom})} + \frac{\sum_{i=1}^{56} \rho_i \cdot B_i \cdot f_i}{\rho(\text{hom})} \cdot T + \frac{\sum_{i=1}^{56} \rho_i \cdot C_i \cdot f_i}{\rho(\text{hom})} \cdot \frac{1}{T^2} \quad (47)$$

- The thermal conductivity is determined by inverse analysis as follows. First, a transient heat transfer simulation is conducted on the heterogeneous REV by imposing different temperatures on two opposite boundaries so as to simulate unidirectional flow. Second, the heat flow curve is obtained and is used as input for the inverse analysis. Third, a simulation with the same conditions is performed on a homogeneous REV with the same size as the heterogeneous volume, giving an initial value for the thermal conductivity. This last value is optimized using the Levenberg-Marquardt method in order to reproduce the same heat curve obtained for the heterogeneous volume.

## 3.4. Thermal properties of the cementitious materials at different scales

The assigned values of the thermal properties (density, thermal conductivity and specific heat capacity) of the cement paste phases are taken from the literature. All the properties are considered as temperature dependent only; the dependency to relative humidity is neglected since, to our knowledge, a thorough study of the influence of relative humidity of most of the cement phases on the thermal properties is currently unavailable. As for the cement phases' thermal conductivity, scalar values are assigned by considering all of them homogeneous and isotropic. Table 5 summarizes the assigned values of the thermal properties of the solid phases. For some phases, incomplete data are filled with those of the closet know mineral phase.

As for the porosity phases (water and ice), there is a strong dependency of their thermal properties on temperature and this is partly the factor controlling freeze-thaw behaviour of cement paste saturated with
water. The density of liquid water presents a maximum at 4°C (Figure 1-a). As water freezes, its density drops from 1 kg/m³ to just 0.916 kg/m³, or about a 9% decrease. The thermal dilation coefficient of liquid water presents an ‘anomaly’, as it changes sign at 4°C (Figure 1-d). On the other hand, both liquid and solid water phases show contrasting trends regarding their thermal conductivity and capacity (Figure 1-b, c). The thermal conductivity of water is low, it varies between 0.2 and 0.6 W.m⁻¹.K⁻¹[37]. As for ice, it is almost four times the conductivity of water (~2.5W.m⁻¹.K⁻¹), and decreases slightly with temperature [38]. Regarding the thermal capacity, for water it is high and equal to 4200 J.kg⁻¹.K⁻¹ and increases sharply at the temperature of homogeneous nucleation in the case of supercooled water [39]. While the ice has a thermal capacity of around 2000 J.kg⁻¹.K⁻¹ and increases linearly with temperature [40] (Figure 1-c).

The homogenized thermal properties of the equivalent fluid (the homogenized density $\rho_{hom}$, and specific thermal capacity $c_{p,hom}$) are given as following according to the volume fraction of solid ice crystal and liquid water:

$$ \rho_{hom}(T) = S_L(T)\rho_L(T) + S_C(T)\rho_C(T) \quad (48) $$

$$ \rho_{hom}(T)c_{p,hom}(T) = S_L(T)\rho_L(T)c_{p,L}(T) + S_C(T)\rho_C(T)c_{p,C}(T) \quad (49) $$

As for the homogenized thermal conductivity $\lambda_{hom}$, it is given using a spherical inclusion model [41], [42]:

$$ \lambda_{hom}(T) = \lambda_L(T) \left( 1 + \frac{S_C(T)}{\frac{1 - S_C(T)}{3} + \frac{\lambda_L(T)}{\lambda_C(T) - \lambda_L(T)}} \right) \quad (50) $$

We neglect the conductivity of water in vapor state.

At the upper scales, the thermal properties of the matrix are determined using homogenization from the lower scales, using the methodology described in section 3.3. The properties can therefore be completely predicted simply from the composition, as it is illustrated in the example studied in this work. Thermal properties are listed in Tables 2 and 3 for mortar and concrete scales respectively.

<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>Figure 1. *Temperature-dependence of the thermal properties of water and ice (solid curves are approximated from* [37], [39])
*(a,b,c) and Linear thermal expansion coefficient of liquid water and ice (d).*</figcaption></figure>

**3.5. Mechanical properties of cementitious materials at different scales**

We used the elastic moduli (*E* and *ν*), the tensile strength (*f*<sub>t</sub>), the fracture energy (*G*<sub>f</sub>) and the linear thermal expansion coefficient (TEC designated by α) as mechanical properties of the cement paste phases. The linear TEC and thermal properties of the cement phases have been determined from the literature and are reported in Table 5 in appendix. Mainly nano-indentation and molecular dynamics tests were found in literature giving the elastic moduli. [43]. For the most cement paste phases, the values of the fracture properties, such as *f*<sub>t</sub> and *G*<sub>f</sub>, are missing because of the difficulty to determine experimentally the exact values for individual phases. Table 5 contains the individual properties of the cement paste phases that are used in the model, this was found by reviewing the corresponding most recently published works found in [43].

The mechanical properties at upper scales are determined by simulating simple tensile and compressive tests on representative volumes of the heterogeneous matrix. The global response is obtained and the homogenized properties are directly extracted from the results, without any inverse calibration. Mechanical properties are listed in Tables 2 and 3 for mortar and concrete scales respectively.

**4. Application of the model to predict the deformation of concrete specimen
subjected to freeze-thaw cycles**

**4.1. Digital representation of the cement paste microstructure**

In this study, we used the VCCTL [44] software based on the CEMHYD3D model because its explicit representation of cement hydrates as well as its pixel architecture can be easily integrated into the Cast3M software used for thermo-mechanical simulation. It is adapted to represent all the phases of cement and allows calculating some elastic properties of the cement paste. However, its limitation lies in the voxel size which is equal to 1 µm. Therefore, it does not allow including all pore size distributions. It is supposed that a 1µm³ voxels of porosity represent a “homogenized” pore distribution. The pore size distribution is estimated using the model described in Section 3.2. This allows a precise estimation of the amount of ice in each voxel, since
the water and ice saturations in each voxel are representative of the water and ice saturations of all capillary pores around that voxel. Figure 2 shows the curves of the pore size distributions of two cement pastes with different w/c ratios (w/c = 0.4 and 0.5). We noticed that 40% of the pore volume is smaller than 2.5 µm. This is an approximate estimation of the existing limit between the strongly adsorbed-evaporable water in capillary pores and the non-evaporable water in the C-S-H gel pores [45]. While the water in the capillary pores is explicitly represented in the microstructure, the water in the C-S-H gel micropores is considered as a part of the homogenized C-S-H voxel. Therefore, in order to account only for the freezable water in the capillary porosity, the cumulated pore volume is adjusted as shown in Figure 2. In the temperature range [+5°C; +20°C], the material undergoes a linear thermal contraction since no phase change of water is occurring. Calculations were limited at 5°C for the highest positive temperature in calculations at cement paste, mortar and concrete scales. So, it was assumed an initial uniform temperature $T_{in} = 5$°C. Dirichlet-type boundary conditions corresponding to the freeze-thaw cycles are applied on the outer boundary of the REV (Figure 3) with temperature cycles ranging from 5°C to -20°C (Figure 4). The theta-method, implemented in the TRANSNON [46] procedure in the finite element software Cast3M, is used for the numerical resolution of the heat equations. At each step, displacements are smaller than the mesh size, so the effect of the depth in deformation can be considered negligible for this problem; therefore all simulations were performed in 2D. At cement sale, a squared REV of 200 µm is used.

<figure><img src="image_6.png" /><figcaption>Figure 2. Pore size distribution curve $\phi(r)$ (dotted), normalized cumulated pore volume curve $V_N(r)$ (dashed line) and corrected curve $V_N(r)$ for the capillary porosity (solid line), for two different w/c ratios.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 3. Thermo-mechanical boundary conditions imposed on the cement paste volume (red color represents CSH phase, blue for porosity, white for other phases)</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 4. Applied thermal cycle +5/-20°C in the numerical simulations.</figcaption></figure>
## 4.2. Digital representation of mortar and concrete samples

For concrete and mortar scales, the mesh generation was conducted with the algorithm developed by Mounajed [47], [48] and used in [49]. It consists of generating a volume by introducing a random distribution of the various inclusions of different sizes in a solid matrix. The distribution algorithm generates random positions of the inclusion gravity centers in the interested area of the volume. At each step of the algorithm, an inclusion i is added. The unit volume of this inclusion of type i is calculated as follows:

$$V_{inclusion} = \pi R_i^2 \qquad (51)$$

where $R_i$ is the inclusion radius. A test is conducted to check that no two inclusions overlap. The distance between the gravity center of an existing inclusion with a radius of $R_1$ and that of a candidate inclusion with a radius of $R_2$ should be higher than the distance $d_{min}$ defined as:

$$d_{min}=R_1 + R_2 \qquad (52)$$

After the validation of the inclusions' positions, the elements were assigned with material characteristics. First, all the elements were defined as having the matrix properties. Then, the model assigned different material properties to the elements by checking the distance d between the gravity center of the element and that of the inclusion. If $d \le R_i$ the inclusions properties are used; else it is the matrix properties that are applied. Note that the interfacial zone between the matrix and the aggregates was not considered in this study. An example of a mortar generated mesh is given in (Figure 5). It is often recommended to take a REV size that is 3 to 4 times larger than the diameter of the largest inclusion [50]. Therefore, at mortar scale for sand inclusions with a maximum diameter of 4 mm, a REV of 1.5 cm is used. At concrete scale, aggregates 4/16 mm are used and the squared REV is size is equal to 5 cm.

<figure><img src="image_9.png" /><figcaption>Figure 5. Example of a mortar microstructure</figcaption></figure>

## 4.3 Strain curves calculations at different scales under freeze and thaw cycles

Figures 6 and 7 show the calculated total strain evolution as a function of time and temperature during 10 freeze-thaw cycles for a W/C = 0.5 and fully saturated cement paste. The total macroscopic strain curves
represent the average strain along the longitudinal and transversal directions calculated over the total surface of the REV, which differ by less than 5%. The upscaling was done using eq. (46) that gives the mean value of $\alpha'_i$ from the calculated mean strain value and the mean temperature of the system at each time step $t$. The presence of salt increases the degradation of the material which can be expressed by increased strain amplitude; damage evolution is illustrated in Figure 8. The presence of salt also delays the ice formation due to supercooling induced by the salt presence. Nucleation starts at -1°C in absence of salts. This temperature drops to -2°C in the case with 3% NaCl.

After these calculations, homogenized parameters of the cement pastes have been calculated to use for the mortar's matrix. These properties, along with the aggregates' properties, are given in Tables 2 and 3. Figures 9 and 10 show the strain evolutions calculated at mortar scales, respectively for 0% NaCl and 3% NaCl. Comparisons with strain curves between cement paste and mortar are also shown on these figures.

<figure><img src="image_10.png" /><figcaption>Figure 6. Evolution of strain in cement paste after 10 freeze-thaw cycles, with 0% NaCl</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Figure 7. Evolution of strain in cement paste after 10 freeze-thaw cycles, with 3% NaCl</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 8. Evolution of damage after 5 cycles in cement paste for different mass fractions of NaCl.</figcaption></figure>
<table><caption>Table 2. *Properties of the mortar's components.*</caption><thead><tr><th>Matrix properties (Homogenized cement paste)</th><th>Siliceous sand inclusions properties<br>(0/4 mm)</th></tr></thead><tbody><tr><td>E = 17.2 (GPa)</td><td>E= 78 (GPa)</td></tr><tr><td>v = 0.227</td><td>v = 0.2</td></tr><tr><td>f<sub>t</sub> = 4 (MPa)</td><td>f<sub>t</sub> = 8 (MPa)</td></tr><tr><td>ρ = 1911 (kg.m<sup>-3</sup>)</td><td>ρ = 2600 (kg.m<sup>-3</sup>)</td></tr><tr><td>λ = 1.2 (W.m<sup>-1</sup>.K<sup>-1</sup>)</td><td>λ = 3.75 (W.m<sup>-1</sup>.K<sup>-1</sup>)</td></tr><tr><td>C<sub>p</sub> = 1647,5 - 0,017 T - <sup>2445000</sup>/<sub>T<sup>2</sup></sub> (J.kg<sup>-1</sup>.K<sup>-1</sup>)</td><td>C<sub>p</sub> = 740 (J.kg<sup>-1</sup>.K<sup>-1</sup>)</td></tr><tr><td>α<sub>p</sub> calculated from eq. (42)</td><td>α<sub>p</sub> = 15x 10<sup>-6</sup> (/K)</td></tr></tbody></table>

<figure><img src="image_13.png" /><figcaption>Figure 9 *Mortar vs heterogeneous cement paste strains as a function of temperature in case of 0% salt concentration*</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Figure 10 *Mortar vs heterogeneous cement paste strains as a function of temperature in case of 3% salt concentration*</figcaption></figure>

As observed, the presence of inclusions leads to a decrease in the maximal and residual strain amplitudes. In fact, the presence of aggregates in the cementitious material gives more rigidity to the resulting two-phase mix due to the rigidity of the aggregates. The result is a denser and stiffer material which can be expressed by lower strain amplitude compared to the cement paste scale. After one freeze-thaw cycle, the strain amplitude reached at -20°C is decreased from 4687µm/m at cement paste scale to 2616µm/m at mortar scale in case of 0% NaCl. When the NaCl concentration was equal to 3% the strain amplitude at -20°C decreases from 5391µm/m for the cement paste to 3033µm/m for the mortar. Moreover, the repetition of freeze-thaw cycles leads to cumulative residual deformations that persist but decrease rapidly, more than in the case of heterogeneous cement paste. In the same way, the expansion rate, after reaching the nucleation point, decreases due to the contribution of aggregates in the thermal contraction of the bi-phasic material when the temperature decreases. Therefore, the beginning of hydraulic expansion is delayed and begins at -1.95°C in mortar while it starts at around -1°C in the cement paste, in case of 0% NaCl. It is also the same situation for cement paste with 3% NaCl, where the nucleation starts at around -2°C due to the presence of salts in the interstitial solution that delays the ice formation. At the corresponding mortar scale, nucleation starts at -2.9°C.
486

<table><caption>Table 3. *Properties of the concrete's components.*</caption><thead><tr><th>Matrix properties (Homogenized mortar)</th><th>Gravels inclusions' properties</th></tr></thead><tbody><tr><td><b>E</b> = 30 (GPa)</td><td><b>E</b> = 85 (GPa)</td></tr><tr><td><b>v</b> = 0,2</td><td><b>v</b> = 0,2</td></tr><tr><td><b>f<sub>t</sub></b> = 4 (MPa)</td><td><b>f<sub>t</sub></b> = 8 (MPa)</td></tr><tr><td><b>ρ</b> = 2218 (kg·m<sup>-3</sup>)</td><td><b>ρ</b> = 2800 (kg·m<sup>-3</sup>)</td></tr><tr><td><b>λ</b> = 1.63 (W·m<sup>-1</sup>·K<sup>-1</sup>)</td><td><b>λ</b> = 3.5 (W·m<sup>-1</sup>·K<sup>-1</sup>)</td></tr><tr><td><b>C<sub>p</sub></b> = 880 (J·kg<sup>-1</sup>·K<sup>-1</sup>)</td><td><b>C<sub>p</sub></b> = 1000 (J·kg<sup>-1</sup>·K<sup>-1</sup>)</td></tr><tr><td><b>α<sub>p</sub></b> calculated from the eq. (42)</td><td><b>α<sub>p</sub></b> = 10x 10<sup>-6</sup> (/K)</td></tr></tbody></table>

487

488

489

490

## 4.4 Experimental testing

491

492

493

Concrete mixture was prepared with Portland CEM II type cement, whose mineral composition is listed in Table 4, and with sand aggregates of 0/4 mm and gravels of 4/20 mm. The water-to-cement ratio was taken equal to 0.5.

494

495

496

<table><caption>Table 4. *Mineral composition of CEM II (CEM II/A-LL 42,5N CE NF).*
(Qualitative analysis by X-ray diffraction).</caption><thead><tr><th rowspan="3">Volume fractions</th><th colspan="5">Primary constituents</th><th colspan="2">Secondary constituents</th></tr><tr><th colspan="4">Clinker phase</th><th>Gypsum</th><th>Calcite</th><th>Slag</th></tr><tr><th>C<sub>3</sub>S</th><th>C<sub>2</sub>S</th><th>C<sub>3</sub>A</th><th>C<sub>4</sub>AF</th><th>C<sub>3</sub>S<sub>2</sub></th><th>C<sub>3</sub>C</th><th>-</th></tr></thead><tbody><tr><td>CEM II</td><td>51.1 %</td><td>12.5 %</td><td>2.9%</td><td>9.6%</td><td>1.1%</td><td>13.8%</td><td>0%</td></tr></tbody></table>

497

498

499

500

501

502

503

504

505

506

507

508

After mixing, concrete specimens were cast into cubic molds of size 15×15×15cm³ for 24h. After demolding they were kept into lime saturated water during 7 days, and then specimens were stored until age of 28 days in 65% RH and 20°C. At the age of 28 days, concrete specimens were configured according to the slab test of the CEN/TS 12390-9 [51] (which is a normative document recognized at the European level) and saturated by maintaining 3mm layer of demineralized water for 72 hours. The concrete specimens were tested during 56 freeze-thaw cycles with thermal amplitude equal to 40°C, ranging from +20 to -20°C, with each cycle lasting 24 hours according to the slab test. The thermal cycle imposed by the technical specification report is shown in Figure 11 with the corresponding spindle. The temperature was recorded on the top of the testing surfaces. The specimens are considered at complete saturation state since the saline solution is maintained on the surface of the specimen throughout the test.

509

510

511

512

513

514

515

Specimens were placed on different shelves in the climatic chamber Binder MKF-240 at GeM laboratory of Centrale Nantes. A good simulation of complex phenomena needs precise data of the boundary conditions. In most cases, the comparison between experimental measurements and simulations considers applied loading conditions, but not the real loading values which can be different due to several variables. Indeed, in the case of freeze-thaw tests the applied temperature condition inside the climatic chamber is not really the same as the one measured on the surface of the specimens. Figure 12 shows the measured temperature in the saline solution of different specimens placed at different locations in the chamber and the maximum and minimum
temperature of the spindle of the CEN/TS 12390-9 normative cycle. The clear discrepancy could have an effect on the heat flux into the specimen, and therefore an effect on the mechanical response. That is why the simulations performed in this study have considered the thermal loading issued from the measurements given in Figure 12.

<figure><img src="image_15.png" /><figcaption>Figure 11. Thermal cycle control spindle according to the slab test of the TS report 12390-9</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Figure 12. Temperature cycles registered in the saline solution of each specimen (code input)</figcaption></figure>

## 4.5 Simulation of the deformation of concrete specimen inside a climatic chamber

Simulations were run at concrete scale on a REV of 5×5 cm² with sand aggregates of 0/4 mm and gravels of 4/16 mm. The matrix properties were calculated by homogenization of the mortar properties at the above scale (Table 3). The results of numerical modeling of the total macroscopic strains of concrete specimens exposed to temperature cycles (given in Figure 12) are presented in Figure 13. Only the first cycle is shown here, as the most of the residual strain is developed during this cycle. The total deformation quickly reaches its highest amplitude, and stabilizes for the rest of the freeze-thaw cycles (true for both cement paste and mortar).

The numerical results highlight a good correlation to be made between the amplitude of the deformation simulated by the model, and the level of degradation to which a surface is exposed to. Figure 14 shows the deformed concrete surfaces noted (S1-S5). A correlation between the degradation level of the test surface and the strain amplitude reached at the minimum temperature can be extrapolated. The most degraded specimen is (S2) that the strain curve shows the maximum amplitude at Tmin. This could reach 3.27 mm/m at the minimum temperature of -21.23°C. Moreover, the residual deformation observed for this same specimen after thawing is the largest, equal to 3 mm/m. This residual deformation appears as a result of damage (microcracks) in the material. However, for the least degraded specimen (S5), the temperature cycle recorded in the brine of this specimen is the smallest in amplitude, going from 10°C to -18°C. This result in a strain curve whose maximum amplitudes at -18°C (2.3 mm/m) and residual after thawing (1.8 mm/m) are small compared to the strain curves of other specimens whose thermal cycles reach lower temperatures. This difference in strains according to the thermal loading cycle imposed in the code could be explained by the Gibbs-Thomson equation which links the freezing temperature to the threshold gel pore radius. As the temperature decreases,
the finer radius pore freezes, which causes higher pressures and thus greater damage. This, in turn, results in maximum strains and significant residual strains. Based on this study, it is possible to say that the cycle that causes significant deformation in the simulation results may cause pronounced surface deformation in the experimental tests. Figures 15 (a and b) shows 3D surface profiles of 2 specimens, S1 and S5 respectively. The depth of degradation can reach 2.5 mm and form a cavity in the surface. This is the case of specimen S1 presented in Figure 15-a, as shown in its 3D cartography where an aggregate is detached from the surface under freeze and thaw cycles. The roughness of the surface is also highly affected after degradation. Figure 16 shows the pore structure in different sections of a concrete specimen of the same formulation presented in this article after being subjected to 56 [+20/-20°C] freeze and thaw cycles. These images are obtained using X-ray microtomography, which is a non-destructive 3D imaging technique that allows reconstructing in 3 dimensions the volume of a sample. An X-ray passes through the sample and is more or less absorbed depending on the density of the material. A series of images is taken under different angles. The system used is an XRadia Micro OCT-400 tomograph equipped with a Hamamatsu 150kV X-ray source and is located in Ecole Centrale de Nantes. A resolution of 8.77 µm is used. Images where then treated using ImageJ. From these sections, the damage of the upper surface can be clearly seen. Moreover, some cracks can be visualized far away from the exposed surface. These cracks are generated because of the movement of water inside the paste inducing pressure development that will damage the material. After exposure to freeze-thaw, the porosity of the concrete becomes connected due to the opening of some new pores and the increase in volume of some others due to ice crystallization and induced water movements. The result is cracks and deformation of the material because of the developed strains, which is very visible in the sections presented in Figure 16. The interesting point of this numerical model lies in its representativeness of the complex phenomena involved in the material subjected to freezing, but also in the possibility of quickly predicting the damage of the exposed material.

<figure><img src="image_17.png" /><figcaption>Figure 13. Total strain of concrete specimens subjected to the freeze-thaw cycles.</figcaption></figure>

566
<figure><img src="image_18.png" /></figure>

S1

<figure><img src="image_19.png" /></figure>

S2

<figure><img src="image_20.png" /></figure>

S3

<figure><img src="image_21.png" /></figure>

S4

<figure><img src="image_22.png" /></figure>

S5

<figure><img src="image_23.png" /><figcaption>Figure 14. Deformed concrete surfaces (S1-S5) after 56 [+20/-20°C] F/T cycles.</figcaption></figure>

(a)

<figure><img src="image_24.png" /></figure>

(b)

<figure><img src="image_25.png" /><figcaption>Figure 16. X-ray microtomography of a concrete specimen after 56 [+20/-20°C] F/T cycles.</figcaption></figure>

## 5. Discussion

Numerical results at the mortar and concrete scales show that inclusions only act as restraint or reinforcement and do not influence the poro-chemical behaviour of the material. This was also seen in the experimental measurements of strains in concrete [20], [52] where the specimen continues to shrink until the nucleation point is reached. This model can be a useful tool for predicting strains of concrete under freeze and thaw cycles. The effect of delayed ice formation due to supercooling in case of saturation with 3% NaCl salt solution was observed. A slight difference in the strain amplitude at minimum temperature is observed as well between concrete and mortars saturated with pure water and others saturated with a 3% NaCl salt solution. This was also found in experimental measurements in previous works [27]. However, it is not possible to
analyze the salt saturation effect at concrete scale. This is because the phenomenon takes place in the pore solution, so at the cement paste scale. An in-depth study is presented in authors' previous publication [18]. In fact, at the nucleation point (-1.79 °C) the hydrostatic pressure increases quickly in the porosity of the cement paste from zero to 17.8 MPa in presence of salt and 14 % of ice are approximately formed, while the osmotic pressure effect increases very slightly from 0 to -0.16 MPa only (the negative sign represents a contracting pressure). At the minimum temperature of -20 °C the hydrostatic pressure rises rapidly and conquers the osmotic pressure; which explains the global and continuous expansion of the paste. The maximum pressure is reached at -20°C.

Figure 17 shows the evolution of the maximum internal pressures (hydrostatic and osmotic pressures) during each cycle of the freeze-thaw process in the presence of 3% NaCl de-icing salts, over 50 cycles. We can notice that the hydrostatic pressure stands initially higher than the osmotic pressure. However, the osmotic pressure undergoes a slight linear increase during freeze-thaw cycles. This increase could be explained by the decrease in the molality of the ions present in the solution of the capillary pores, where the hydrostatic pressure develops and causes pores expansion. Maximum hydrostatic pressure, on the other hand, stands constant for few cycles, until the pores are large enough to contain more water, and then gradually decreases. The magnitude of the hydrostatic pressure only gradually decreases and never reaches zero before 50 cycles.

<figure><img src="image_26.png" /><figcaption>Figure 17. Evolution of the maximum internal pressures during each cycle of the freeze-thaw process in the presence of de-icing salts over 50 cycles</figcaption></figure>

A slight variation in the saturation state of the specimen can highly affect the strain amplitude, as shown in the previous authors' article [18].

On the other part, a lack of reproducibility and repeatability of frost tests on concrete are noticed due to complexity of the material and many mechanisms that may be dependent on the experimental conditions. As observed in this work, the imposed cycle is not applied exactly with the same amplitude and freezing and thawing rates in all the positions inside the climatic chamber. The position of the specimen in the chamber determines the thermal cycle that it is in fact undergoing, which can be largely different from the imposed one, such as the normative cycle. This experimental condition can possibly explain the discrepancy that is found in the different literature studies. Another important point to mention here is that the applied thermal
cycles imposed according to the standardized methods are in fact questionable whether they are representative of real exposure conditions or not. An in-depth study of the real temperature evolutions in reputedly severe climates is presented authors' previous publication [53]. This latter calls into question the reliability of laboratory tests for the characterization of the resistance of concrete to freeze and thaw.

## 6. Conclusions

In this work, a numerical study of the effect of minimum temperature, saturation degree and de-icing salt concentration on cement-based material under freeze and thaw cycles was done using a micro-mechanical approach developed by [54]. This approach takes basic variables involved in freeze and thaw mechanisms: temperature, internal developed pressure, saturation degree curve which is a governing parameter, pore size distribution calculated according to water content of the microstructure used in the study. An upscaling approach was also performed, where the deformations as well as all the mechanical and thermal properties of the heterogeneous cement paste, then mortar, were homogenized for the upper scales.

An application of this tool in modeling the deformation of concrete specimens subjected to real freeze-thaw cycles inside a climatic chamber showed very good correlation between the degradation level of the real concrete test surface and the amplitudes of maximum and residual deformations. The most deteriorated specimens showed a curve of deformation with the highest amplitude and the highest value of residual deformations. The effect of the presence of aggregates was noticeable. The maximum amplitude was reduced as well as the velocity of increment of the residual deformation since the aggregates restrict the strains and stiffens the matrix. Many interesting issues were obtained by isolating each parameter involved in the freeze and thaw action, which is not always possible to accomplish in the laboratory experiments because of the complexity of the involved mechanisms. The effect of important parameters such as the initial saturation degree was previously examined by the authors in a previous study, and its major influence on the global freeze-thaw behaviour was noticeable. The authors confirmed the existence of a critical saturation degree already found by Fagerlund, above which an important expansion of the material is observed that leads to damage.

Therefore, we can see how this numerical model offers a reliable tool for predicting the deformation of concrete exposed to freeze-thaw by taking into account the effects of pressures in the porosity developing at the microscopic scale by involving the properties of all phases of the cement paste. The different laboratory exposure parameters can be taken into account into the model and representative results can be obtained without any need of calibration.

This model is developed with the objective to be used to drive and control frost tests in view of the discrepancies that are found in the existing experimental methods for testing concrete resistance to freeze-thaw. Significant differences in the thermal conditions applied to the specimens placed at different positions in the chamber during a freeze-thaw test were noticed. The temperatures in the brines of the tested specimens are not uniform. The modeling of the strains of the concrete REVs subjected to these different thermal loadings
also shows a difference in the order of magnitude of the reached modeled strains. This model is used within the framework of the proposition of a new protocol for predicting the degradation of concrete in the winter period by replacing the existing standardized cycle of [+20/-20°C] and improving the existing testing method. This further development will be the topic for a forthcoming paper. Finally, in a future work a statistical study using this model will be done to assess the scaling risk in comparison with experiments.

## 7. Acknowledgments

The financial support of this work by ATILH (Technical Association of French Cement Industry) now absorbed by SFIC (Professional Syndicate of French Cement Industry) and the grant by ANRT by mean of CIFRE convention n°2017/891 for Sara AL HAJ SLEIMAN are greatly acknowledged.

## 8. Appendix

<table><caption>Table 5. Mechanical and thermal properties of cement phases used in this study. This table is taken from [43].</caption><thead><tr><th colspan="2">Cement phases</th><th>Young's modulus (GPa)</th><th>Poisson's ratio (-)</th><th>Tensile strength (MPa)</th><th>Fracture energy (J/m²)</th><th colspan="3">Thermal expansion coefficients</th><th>Density (kg/m³)</th><th>Thermal conductivity (Wm⁻¹K⁻¹)</th><th colspan="3">Specific heat capacity</th></tr><tr><th></th><th></th><th></th><th></th><th></th><th></th><th>α₀ (10⁻⁶K⁻¹)</th><th>α₁ (10⁻⁸K⁻²)</th><th>α₂ (K)</th><th></th><th></th><th>cₚ₁ (Jkg⁻¹K⁻¹)</th><th>cₚ₂ (10⁻³Jkg⁻¹K⁻²)</th><th>cₚ₃ (10⁵Jkg⁻¹K)</th></tr></thead><tbody><tr><td rowspan="7">Clinker phases</td><td>Alite</td><td>C₃S</td><td>137.4</td><td>0.299</td><td>430.7</td><td>65.8</td><td>6.2</td><td>0.8</td><td>0</td><td>3.2</td><td>3.4</td><td>710</td><td>0</td><td>0</td></tr><tr><td>Belite</td><td>C₂S</td><td>135.5</td><td>0.297</td><td>396.0</td><td>56.3</td><td>15.3</td><td>0</td><td>-0.4</td><td>3.3</td><td>3.5</td><td>720</td><td>0</td><td>0</td></tr><tr><td>Celite</td><td>C₃A</td><td>145.2</td><td>0.278</td><td>534.6</td><td>94.3</td><td>8.5</td><td>0.3</td><td>-0.3</td><td>3</td><td>3.5</td><td>720</td><td>0</td><td>0</td></tr><tr><td>Ferrite</td><td>C₄AF</td><td>150.8</td><td>0.318</td><td>470.3</td><td>72.5</td><td>8.2</td><td>0.1</td><td>-0.</td><td>3.7</td><td>3.5</td><td>720</td><td>0</td><td>0</td></tr><tr><td>Arcanite</td><td>KS</td><td>44.2</td><td>0.269</td><td>29.7</td><td>1.0</td><td>10.8</td><td>6.6</td><td>0</td><td>2.7</td><td>1</td><td>754</td><td>0</td><td>0</td></tr><tr><td>Thenardite</td><td>N'S</td><td>57.1</td><td>0.281</td><td>62.9</td><td>3.3</td><td>10.8</td><td>6.6</td><td>0</td><td>2.7</td><td>1</td><td>741</td><td>815</td><td>-72.4</td></tr><tr><td>Gypsum</td><td>C₅H₅</td><td>44.5</td><td>0.33</td><td>29.6</td><td>1.0</td><td>51.3</td><td>-3.4</td><td>0.7</td><td>2.3</td><td>3</td><td>1081</td><td>0</td><td>0</td></tr><tr><td rowspan="3">Gypsum</td><td>Bassanite</td><td>C₅H₆</td><td>132.0</td><td>0.25</td><td>61.7</td><td>1.4</td><td>51.3</td><td>-3.4</td><td>0.7</td><td>2.7</td><td>4.3</td><td>823</td><td>0</td><td>0</td></tr><tr><td>Anhydrite</td><td>C'S</td><td>88.8</td><td>0.233</td><td>93.6</td><td>4.6</td><td>51.3</td><td>-3.4</td><td>0.7</td><td>3</td><td>5.7</td><td>753</td><td>462</td><td>-131</td></tr><tr><td>Slags</td><td>§</td><td>137.4</td><td>0.299</td><td>430.7</td><td>65.8</td><td>6.2</td><td>0.8</td><td>0</td><td>2.3</td><td>3.4</td><td>710</td><td>0</td><td>0</td></tr><tr><td rowspan="3">Admixtures</td><td>Silica fume</td><td>§</td><td>71.7</td><td>0.169</td><td>346.5</td><td>73.4</td><td>6</td><td>0.8</td><td>-0.1</td><td>2.2</td><td>1.4</td><td>754</td><td>218</td><td>318</td></tr><tr><td>AS glass</td><td>C₅S₃A₂</td><td>71.7</td><td>0.169</td><td>346.5</td><td>73.4</td><td>-</td><td>-</td><td>-</td><td>3.3</td><td>0.1</td><td>730</td><td>0</td><td>0</td></tr><tr><td>Inert</td><td>§</td><td>71.0</td><td>0.296</td><td>69.3</td><td>3.3</td><td>7.1</td><td>3.4</td><td>-1.2</td><td>3.2</td><td>3.4</td><td>996</td><td>269</td><td>-216</td></tr><tr><td rowspan="12">Main hydrates</td><td>Portlandite</td><td>CH</td><td>43.5</td><td>0.294</td><td>73.8</td><td>6.1</td><td>20</td><td>1.1</td><td>0</td><td>2.2</td><td>1.3</td><td>1205</td><td>447</td><td>-140</td></tr><tr><td>C-S-H</td><td>C₁₇SH₄</td><td>23.8*</td><td>0.24*</td><td>55.0</td><td>5.9</td><td>15</td><td>0</td><td>0</td><td>2.1</td><td>1</td><td>930</td><td>0</td><td>0</td></tr><tr><td>Hydro-garnet</td><td>C₃AH₆</td><td>93.8</td><td>0.32</td><td>324.7</td><td>55.6</td><td>6.5</td><td>0.3</td><td>-0.2</td><td>2.5</td><td>5.5</td><td>772</td><td>1483</td><td>0</td></tr><tr><td>Ettringite</td><td>C₄A'S₄H₆</td><td>24.1</td><td>0.321</td><td>39.6</td><td>3.2</td><td>35.3</td><td>0</td><td>0</td><td>1.8</td><td>1</td><td>1545</td><td>629</td><td>0</td></tr><tr><td>Monosulfate</td><td>C₄A'S₄H₁₂</td><td>43.2</td><td>0.292</td><td>262.4</td><td>77.3</td><td>35.3</td><td>0</td><td>0</td><td>2</td><td>1</td><td>954</td><td>1876</td><td>0</td></tr><tr><td>Iron hydroxide</td><td>FH₃</td><td>22.4</td><td>0.25</td><td>292.1</td><td>178.5</td><td>20</td><td>1.1</td><td>0</td><td>3.1</td><td>1</td><td>951</td><td>0</td><td>0</td></tr><tr><td>Calcite</td><td>C'C</td><td>97.0</td><td>0.296</td><td>69.3</td><td>2.4</td><td>2.4</td><td>1.1</td><td>-0.4</td><td>2.7</td><td>3.4</td><td>996</td><td>269</td><td>-216</td></tr><tr><td>Limestone</td><td>C</td><td>183.9</td><td>0.21</td><td>111.2</td><td>3.1</td><td>10.1</td><td>0.3</td><td>0</td><td>3.4</td><td>2.6</td><td>750</td><td>0</td><td>0</td></tr><tr><td>Iron-rich ettringite</td><td>C₆(A·F)</td><td>24.1</td><td>0.321</td><td>39.6</td><td>3.2</td><td>35.3</td><td>0</td><td>0</td><td>1.8</td><td>1</td><td>1545</td><td>629</td><td>0</td></tr><tr><td>Pozzolanic / slag C-S-H</td><td>§</td><td>23.8</td><td>0.24</td><td>55.0</td><td>5.9</td><td>15</td><td>0</td><td>0</td><td>2.1</td><td>1</td><td>930</td><td>0</td><td>0</td></tr><tr><td>Friedel's salt</td><td>C₃A(C₆)</td><td>22.4</td><td>0.25</td><td>76.2</td><td>12.2</td><td>15</td><td>0</td><td>0</td><td>1.9</td><td>1</td><td>930</td><td>0</td><td>0</td></tr><tr><td>Strütingite</td><td>C₂ASH₈</td><td>22.4</td><td>0.25</td><td>328.2</td><td>225.4</td><td>7.7</td><td>0.09</td><td>0</td><td>1.9</td><td>1</td><td>1274</td><td>0</td><td>0</td></tr><tr><td>Calcium mono-carbonate</td><td>C₄A'C₄H₃</td><td>79.0</td><td>0.25</td><td>262.3</td><td>40.9</td><td>-</td><td>-</td><td>-</td><td>2.1</td><td>1</td><td>1299</td><td>0</td><td>0</td></tr><tr><td>Brucite</td><td>MH</td><td>83.3</td><td>0.183</td><td>49.5</td><td>1.3</td><td>26.7</td><td>0</td><td>0</td><td>2.4</td><td>1</td><td>1752</td><td>259</td><td>-449</td></tr><tr><td>Magnesium sulphate</td><td>M'S</td><td>10.5</td><td>0.1</td><td>41.6</td><td>6.8</td><td>51.3</td><td>-3.4</td><td>0.6</td><td>2.7</td><td>1</td><td>799</td><td>0</td><td>0</td></tr><tr><td rowspan="3">Secondary hydrates</td><td>Porosity (wet)</td><td>H</td><td>1.9**</td><td>0.499**</td><td>0.001†</td><td>-0†</td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Porosity (dry)</td><td>Ø</td><td>1.9**</td><td>0.001**</td><td>0.001</td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Porosity (ice)</td><td>I<sub>h</sub></td><td>9.4‡</td><td>0.42‡</td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr></tbody></table>

<footnote>* Elastic properties of C-S-H are homogenized, assuming that C-S-H at the microscale is a mixture of 70% LD C-S-H and 30% HD C-S-H</footnote>

<footnote>** Calibrated using the experimentally-measured Young's modulus and Poisson's ratio of cement paste CP1-0.45 as a basis for numerical simulations of elastic behaviour.</footnote>

<footnote>† Calibrated value for tension simulations using the experimentally-determined tensile strength as a basis for calibration.</footnote>

<footnote>†† Calibrated value for compression simulations using the experimentally-determined compressive strength as a basis for calibration.</footnote>

<footnote>‡ Values given at 0°C and atmospheric pressure.</footnote>
[1] J. J. Valenza and G. W. Scherer, "A review of salt scaling: I. Phenomenology," Cem. Concr. Res., vol. 37, no. 7, pp. 1007-1021, 2007, doi: 10.1016/j.cemconres.2007.03.005.

[2] S. Lindmark, "Mechanisms of Salt Frost Scaling on Portland Cement-bound Materials: Studies and Hypothesis," Lund, Sweden, 1998.

[3] T. Powers, "The air requirement of frost resistant concrete," in Proceedings of the Twenty-Ninth Annual Meeting of the Highway Research Board, 1949, pp. 184-211.

[4] R. A. Powers, T. C., and Helmuth, "Theory of volume changes in hardened Portland cement paste during freezing," in Highway Research Board Proceedings 32, 1953, vol. 32, pp. 285-297.

[5] G. W. Scherer, "Stress from crystallization of salt," Cem. Concr. Res., vol. 34, no. 9, pp. 1613-1624, 2004, doi: 10.1016/j.cemconres.2003.12.034.

[6] M. Pigeon and M. Regourd, "Effects of Freeze-Thaw Cycles on the Microstructure of Hydration Products.," in Durability of building materials, vol. 4, no. 1, 1986, pp. 1-19.

[7] A. Müller, C. Fuhr, and D. Knöfel, "Frost resistance of cement mortars with different lime contents," Cem. Concr. Res., vol. 25, no. 4, pp. 809-818, 1995, doi: 10.1016/0008-8846(95)00071-J.

[8] C. Girodet, "Endommagement des mortiers sous sollicitations thermo-mécaniques: Influence des caractéristiques des matériaux sur la résistance au gel interne," Institut National des sciences appliquées de Lyon, 1996.

[9] J. J. Valenza and G. W. Scherer, "A review of salt scaling: II. Mechanisms," Cem. Concr. Res., vol. 37, no. 7, pp. 1022-1034, 2007, doi: 10.1016/j.cemconres.2007.03.003.

[26] V. Penttala and F. Al-Neshawy, "Stress and strain state of concrete during freezing and thawing cycles," Cem.
720 Concr. Res., vol. 32, no. 9, pp. 1407-1420, 2002, doi: 10.1016/S0008-8846(02)00785-8.

721 J. Kaufmann, "Experimental identification of damage mechanisms in cementitious porous materials on phase
722 transition of pore solution under frost deicing salt attack," Swiss Federal Institute of Technology Lausanne
723 (EPFL), 2000.

724 [28] F. Gong, D. Zhang, E. Sicat, and T. Ueda, "Empirical estimation of pore size distribution in cement, mortar, and concrete," *J. Mater. Civ. Eng.*, vol. 26, no. 7, pp. 1-11, 2014, doi: 10.1061/(ASCE)MT.1943-5533.0000945.

726 [29] Q. Zeng, "Discussion of 'Empirical estimation of pore size distribution in cement, mortar, and concrete' by Fuyuan Gong, Dawei Zhang, Evdon Sicat, and Tamon Ueda," *J. Mater. Civ. Eng.*, vol. 27, no. 11, pp. 1-3, 2015, doi: 10.1061/(ASCE)MT.1943-5533.0001298.

729 [30] B. Zuber and J. Marchand, "Modeling the deterioration of hydrated cement systems exposed to frost action: Part 1: Description of the mathematical model," *Cem. Concr. Res.*, vol. 30, no. 12, pp. 1929-1939, Dec. 2000, doi: 10.1016/S0008-8846(00)00405-1.

732 [31] G. W. Scherer, "Crystallization in pores," *Cem. Concr. Res.*, vol. 29, no. 8, pp. 1347-1358, 1999, doi: 10.1016/S0008-8846(99)00002-2.

734 [32] M. Brun, A. Lallemand, J.-F. Quinson, and C. Eyraud, "A new method for the simultaneous determination of the size and shape of pores: the thermoporometry," *Thermochim. Acta*, vol. 21, no. 1, pp. 59-88, Oct. 1977, doi: 10.1016/0040-6031(77)85122-8.

737 [33] O. Coussy and P. J. M. Monteiro, "Poroelastic model for concrete exposed to freezing temperatures," *Cem. Concr. Res.*, vol. 38, no. 1, pp. 40-48, Jan. 2008, doi: 10.1016/J.CEMCONRES.2007.06.006.

739 [34] J. Kaufmann, "Experimental identification of damage mechanisms in cementitious porous materials on phase transition of pore solution under frost deicing salt attack," 1999.

741 [35] Y. Xi, Z. P. Bažant, and H. M. Jennings, "Moisture diffusion in cementitious materials - Adsorption isotherms," *Adv. Cem. Based Mater.*, vol. 1, no. 6, pp. 248-257, Nov. 1994, doi: 10.1016/1065-7355(94)90033-7.

743 [36] H. Lin and L. Lee, "Estimations of activity coefficients of constituent ions in aqueous electrolyte solutions with the two-ionic-parameter approach," vol. 237, pp. 1-8, 2005, doi: 10.1016/j.fluid.2005.08.005.

745 [37] J. W. Biddle, V. Holten, J. V Sengers, and M. A. Anisimov, "Thermal conductivity of supercooled water," 2013.

746 [38] L. J. Bonales, A. C. Rodriguez, and P. D. Sanz, "Thermal conductivity of ice prepared under different conditions," *Int. J. Food Prop.*, vol. 20, no. sup1, pp. 610-619, Dec. 2017, doi: 10.1080/10942912.2017.1306551.

749 [39] V. Holten, C. E. Bertrand, M. A. Anisimov, and J. V. Sengers, "Thermodynamics of supercooled water," *J. Chem. Phys.*, vol. 136, no. 9, p. 094507, Mar. 2012, doi: 10.1063/1.3690497.

751 [40] W. F. Giauque and J. W. Stout, "The Entropy of Water and the Third Law of Thermodynamics. The Heat Capacity of Ice from 15 to 273°K.," *J. Am. Chem. Soc.*, vol. 58, no. 7, pp. 1144-1150, Jul. 1936, doi: 10.1021/ja01298a023.

754 [41] E. Herve, "Thermal and thermoelastic behaviour of multiply coated inclusion-reinforced composites," *Int. J. Solids Struct.*, vol. 39, no. 4, pp. 1041-1058, Feb. 2002, doi: 10.1016/S0020-7683(01)00257-8.

756 [42] A. Fabbri, "Physico-mécanique des matériaux cimentaires soumis au gel-dégel," Université de Marne-la-Vallée, 2006.

758 [43] A. Rhardane, F. Grondin, and S. Y. Alam, "Development of a micro-mechanical model for the determination of damage properties of cement pastes," *Constr. Build. Mater.*, Nov. 2020, doi: https://doi.org/10.1016/j.conbuildmat.2020.120514.

761 [44] J. W. Bullard and P. E. Stutzman, "Analysis of CCRL Portland Cement Proficiency Samples Number 151 and Number 152 Using the Virtual Cement and Concrete Reference Laboratory," *Cem. Concr. Res.*, vol. 36, pp. 1548-1555, 2006.

764 [45] A. B. Poole and I. Sims, *Concrete petrography : a handbook of investigative techniques*. 2016.

765 [46] T. Charras and F. Di Paola, "Utilisez Cast3M, &quot;La procédure PASAPAS&quot;," 2011.

766 [47] G. Moujaned, "Exploitation du nouveau modèle Béton numérique dans symphonie: Concept, homogénéisation du comportement thermomécanique des BHP et simulation de l'endommagement thermique, Centre Scientifique et Technique du Bâtiment," 2002.

769 [48] F. Grondin, H. Dumontet, A. Ben Hamida, G. Mounajed, and H. Boussa, "Multi-scales modelling for the behaviour of damaged concrete," *Cem. Concr. Res.*, vol. 37, no. 10, pp. 1453-1462, 2007, doi: 10.1016/j.cemconres.2007.05.012.

772 [49] M. Guo, F. Grondin, and A. Loukili, "Numerical method to model the creep of recycled aggregate concrete by considering the old attached mortar," *Cem. Concr. Res.*, vol. 118, no. January, pp. 14-24, 2019, doi: 10.1016/j.cemconres.2019.01.008.

775 [50] A. A. Gusev, "Representative volume element size for elastic composites: A numerical study," *J. Mech. Phys. Solids*, vol. 45, no. 9, pp. 1449-1459, 1997, doi: 10.1016/S0022-5096(97)00016-1.

777 [51] C. 12390-9 Document, "CEN/TS 12390-9, "Testing hardened concrete. Freeze-thaw resistance- Scaling",," 2017.

778 [52] J. P. Kaufmann, "Experimental identification of ice formation in small concrete pores," *Cem. Conc. Res.*, pp. 1421-1427, 2004.
780 [53] N. et al. Izoret, L., Al Haj Sleiman, S., Matoiri-Chaibati, "Concrete in a severe freezing environment: a meteorological characterization," *Mater. Struct.*, vol. 54, no. 1, 2021, doi: 10.1617/s11527-020-01603-8.

781
782 [54] A. Rhardane, "Élaboration D'Une Approche Micromécanique Pour Modéliser L'Endommagement Des Matériaux Cimentaires Sous Fluage Et Cycles De Gel-Dégel," p. 250, 2018.

784