# Diffusion and physicochemical behavior of chloride and sodium ions in C-S-H gel pores: A molecular dynamics investigation

Tarek Ihaddadene, Jérôme Claverie, François Bignonnet, Ouali Amiri

► To cite this version:

Tarek Ihaddadene, Jérôme Claverie, François Bignonnet, Ouali Amiri. Diffusion and physicochemical behavior of chloride and sodium ions in C-S-H gel pores: A molecular dynamics investigation. Cement and Concrete Research, 2026, 199, pp.108053. 10.1016/j.cemconres.2025.108053. hal-05329948

HAL Id: hal-05329948

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

Submitted on 24 Oct 2025

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
Diffusion and physicochemical behavior of chloride and sodium ions in C-S-H
gel pores : A molecular dynamics investigation

Tarek Ihaddadene<sup>a</sup>, Jérôme Claverie<sup>a</sup>, François Bignonnet<sup>a</sup>, Ouali Amiri<sup>a</sup>

<sup>a</sup>Nantes Université, École Centrale Nantes, CNRS, GeM, UMR 6183, F-44600 Saint-Nazaire, France

**Abstract**

Calcium silicate hydrate (C-S-H) gel is an amorphous material with a complex, disordered structure that complicates the understanding of its atomic-scale properties. While many studies have investigated ionic diffusion in C-S-H pores, the behavior of sodium (Na<sup>+</sup>) and chloride (Cl<sup>-</sup>) ions, particularly in relation to the Ca/Si ratio, is not well understood. In this study, atomistic models of C-S-H were developed with varying Ca/Si ratios, pore sizes, and NaCl concentrations. Molecular dynamics simulations were used to calculate the self-diffusion profiles of Na<sup>+</sup> and Cl<sup>-</sup> ions. Results indicate that pore size strongly affects diffusivity and adsorption through confinement and electrical double layer effects. The solid surface influences ionic mobility up to approximately 1.2 nm. The Ca/Si ratio has minimal impact on the diffusion profiles of non-adsorbed ions, but its rise enhances Cl<sup>-</sup> mobility near the surface and increases chloride binding capacity, consistent with prior studies.

***Keywords:*** Molecular dynamics, Calcium silicate hydrates (C-S-H), Chloride Diffusion, Adsorption, Concrete durability

**1. Introduction**

Corrosion of steel due to the diffusion of carbon dioxide and/or chlorides into reinforced concrete struc-
tures is the main cause of their degradation, leading to significant repair costs. Predicting the service life is
crucial to optimize maintenance operations for these structures. In the case of chlorides, numerous predic-
tion models of the diffusion coefficient have been developed at various microscopic and macroscopic scales
[1–5]. These models demonstrate that the low ratio between the diffusion coefficient of an ionic species in
concrete and its value in bulk solution cannot solely be attributed to a tortuosity factor depending on the
geometry of the porous medium. Most of the models proposed for predicting the effective diffusivity in
concrete require the calibration of parameters, particularly those related to the diffusivity in gel pores (\~
1-10 nm) of calcium silicate hydrates (C-S-H), the primary hydrate in Portland cement paste [6]. This shows
that interactions between ions and the pore surface of C-S-H gels must be considered, including through the
phenomenon of the electrical double layer (EDL) [7]. Molecular dynamics (MD) is a suitable method for
studying such interactions and has been widely used to explore ion behavior in the C-S-H gel pores.

The first MD-based study of the C-S-H structure was carried out in 1996 by Faucon et al. [8, 9]. Their main goal was to find the causes of the structural instability of C-S-H chains and look into the mechanisms that lead to the failure. They built two C-S-H models based on tobermorite structure. Starting from the 1.1 nm tobermorite structure with Ca/Si=0.66 [10], the Ca/Si ratio of 0.83 was attained by replacing two protons with a calcium ion. Although jennite was first considered to be a potential basis to build atomistic C-S-H models, all recent studies use tobermorite because of its structural similarity with C-S-H. Since the method proposed by Pellenq in 2009 to build an atomistic C-S-H model [11], consisting in silicate chain modification of tobermorite and insertion of water molecules with Grand Canonical Monte Carlo (GCMG) some improvements were put forward, in particular to reproduce the absence of silicate monomers, supported by experimental studies [12, 13]. More recently, a brick C-S-H model was proposed to build C-S-H model from simple building blocks, allowing to introduce defect in a 14 Å tobermorite model to obtain structural properties in good agreement with experimental measurements [14, 15].
The structure of C-S-H, and especially the Ca/Si ratio, influence its physical properties. Smaller Ca/Si ratio corresponds to C-S-H with longer silicate chains, resulting in greater rigidity and mechanical strength while larger Ca/Si ratio corresponds to a more disordered, glassy structure with softer mechanical behavior [12, 16]. Regarding the influence of the Ca/Si ratio on the retention and diffusion of ions, an MD investigation has proven that the adsorption capacity of cesium ions is decreasing and the diffusivity increasing with the Ca/Si ratio [17]. The increase in adsorption capacity with the Ca/Si ratio was attributed to a better accessibility of the pairing silicate to the cations.

Little is known about the diffusion of NaCl solutions at different molar concentrations and most of the existing MD investigations on chlorides only considered tobermorite as a model for C-S-H [18–24]. Using the 0.9 nm tobermorite with fully protonated surface silicates as a model for C-S-H gel, Kalinichev and Kirkpatrick [18] focused their study on the chloride binding at the surface of several minerals. They described metal-oxygen interactions using the ClayFF force field [25], which allowed flexibility of surface hydroxyl groups. No binding of chlorides ions was observed at the surface of tobermorite in their simulations. The authors stated that the chloride adsorption capacity on silicates should be even lower in real C-S-H, where the Si-O sites should be deprotonated at high pH values expected for most cements.

To our knowledge, only two studies investigating chlorides diffusion considered a more realistic C-S-H structure with varying Ca/Si ratios [26, 27]. One of these studies considered the saturated (ion-transport) in pores of few nanometers with non-equilibrium MD simulations. Simulations of unsaturated transport in gel pores are also reported in the literature [23, 24]. However, this type of transport may only occur under extreme drying conditions and is absent in structures located in coastal environments exposed to high relative humidity. According to GCMG simulations, it may only occur at low relative humidity since a pore of one nanometer is completely filled down to HR=20% at 300 K [28].

Other molecular simulations of the pores in C-S-H gel, focused on the interfaces of tobermorite-tobermorite, tobermorite-jennite, and jennite-jennite [29]. Although the surfaces were hydroxylated, the Ca/Si ratio was not adjusted to match experimental characterization of real C-S-H structures. Additionally, the impact of sodium chloride concentration on adsorption and diffusion processes was not addressed in this study. Other researchers proposed an extensive study on chloride binding in cement paste including both experimental measurements and molecular dynamics simulations [30]. Defective C-S-H and calcium aluminosilicate hydrates (C-A-S-H) structures were considered and the number of chemical bonds created between chlorides ions and Ca²⁺ was found to increase with the Ca/Si ratio. However, the diffusion of chloride and sodium ions was not studied in this article.

Furthermore, in many studies considering the motion of chemical species at the interface with a solid, the diffusion coefficient is calculated for the whole pore, independently of the position of the species. This is a limitation, as a thorough analysis of the trajectory enables the determination of a diffusion coefficient profile as a function of the distance from the surface. Such a type of analysis may enhance the interpretation and the discussion on the origin of the low diffusion coefficient of chlorides measured experimentally in concrete. To address the aforementioned limitations, the adsorption and diffusivity of chloride and sodium ions in saturated C-S-H gel pore are systematically analyzed based on MD simulations. The main objective is to investigate the impact of the molecular properties of C-S-H on the structural and dynamic properties of the species present in the gel pore. This paper makes several original contributions to the literature, notably: (a) atomistic C-S-H structures in agreement with experimental data are constructed following a detailed methodology, which, to our knowledge, is lacking in existing MD studies on chloride diffusion; (b) the diffusion of Cl⁻ and Na⁺ ions is thoroughly investigated with respect to directionality and the different diffusion regimes are discussed; and (c) the diffusion analysis is position-dependent, offering detailed diffusion profiles across the pore and explicitly differentiating between interfacial and bulk regions

To our knowledge, these aspects have not been addressed, or have only been partially explored, in previous studies on the interaction and diffusion of chloride and sodium ions in C-S-H.

In Section 2, the process of creating the C-S-H gel pore models with Ca/Si ratios ranging from 1.2 to 1.8 is explained, involving modifications to the crystalline structure of 11 Å tobermorite. The molecular models under study have Ca/Si, H₂O/Si ratios, and silicate chain length in agreement with experimental measurements. Gel pores from 5 Å to 50 Å are considered with idealized pore solution constituted of inter-layer calcium ions and NaCl solution with concentration of 0.5 M to 4 M. The Section 3 presents a
detailed study of the results on the adsorption mechanism, focusing on the structure of ions at the interface, bond lifetime, and the quantification of free and bound ions. Finally the Section 5 provides an analysis of the different diffusion regime and anisotropy within the pore. In addition, the influence of the pore size, concentration and Ca/Si ratios on the self-diffusion coefficient profiles is investigated. Finally, a conclusion and an outlook are presented in Section 6.

## 2. Methods

### 2.1. *Creation of C-S-H models*

Based on a review of the state-of-the-art [31], atomistic C-S-H models are created from a defective 11 Å tobermorite structure [10] (a = 6.69 Å, b = 7.39 Å, c = 22.77 Å). The steps to create the models are depicted in Fig. 1 and presented below.

<figure><img src="image_2.png" /><figcaption>Figure 1: Workflow of C-S-H model creation steps, with varying Ca/Si ratios, starting from the Hamid 11 Å tobermorite structure [10]. Calcium atoms are represented in lime green, silicon in yellow, oxygen in red, hydrogen in white, chloride in cyan and sodium in purple.</figcaption></figure>

All the water molecules are removed from the interlaminar space to adjust the water content afterwards. The unit cell is scaled by a factor of 3 along the “a” direction (crystal base vector) and 5 along the “b” direction. The length in the “c” direction remains unchanged. Then, bridging SiO₂ are first randomly removed to reach the targeted Ca/Si ratio. The scaling of the system allows to reach more accurately the desired Ca/Si ratios, however the dimensions are kept small enough to limit the computation time with the reactive force field. The mean chain length of silicates (MCL) are kept larger than 3, in accordance with the experimental findings of Richardson [32]. The value of MCL in this study aligns with that reported by Qomi et al. [12]. Then, water molecules are inserted into the interlayer pore of the initial model with Packmol [33] to reach the desired H₂O/Si ratio. When the removal of SiO₂ alone is insufficient to reach the targeted Ca/Si ratio, Ca²⁺ ions are introduced into the interlayer pore, and H⁺ ions are removed from random water molecules to maintain charge balance, thus creating OH⁻ ions.
The H₂O/Si ratio is calculated as a function of the Ca/Si ratio by Eq. (1), determined through interpolation on a collection of experimental results [32]:

$$ \frac{\text{H}_2\text{O}}{\text{Si}} = \frac{19 \text{Ca}}{17 \text{Si}} - \frac{7}{17} \qquad (1) $$

This work builds on an established method involving the random removal of SiO₂ groups and the addition of Ca²⁺ ions to adjust the Ca/Si ratio. To better match experimental observations, a minimal MCL is also imposed. In some studies, dangling oxygen atoms resulting from silicate removal are compensated by extra Ca²⁺ ions [17, 34]. In contrast, no compensation is applied upon each SiO₂ group removal, as these units are charge neutral and formal charge balance is therefore preserved [12]. Instead, the Ca²⁺ are added once the MCL as reached a minimum value of 3. Instead, Ca²⁺ ions are added only after the MCL reaches the target minimum. These differences in methodology may result in small differences between our structures and those reported in the literature, particularly in the evolution of the MCL as a function of the Ca/Si ratio.

Once the base models are created, several MD simulations are conducted using the LAMMPS code [35] to prepare the final models. Initially, the ReaxFF force field [36] is applied with an optimized set of parameters for Ca/H/O/Si elements [37–39] in the NPT ensemble for a simulation duration of 500 ps with a time step dt = 0.2 fs. The application of this reactive force field aims to dissociate water molecules and form hydroxyl groups in the interlayer. Also, as already reported in the literature, a minor condensation of silicate groups is observed [12].

The ReaxFF uses a polarizable charge calculation approach, where the force field parameters are applied to each element, regardless of the atomic environment. However, during subsequent equilibration using a classical force field, the parameters for each atomic type must be defined. The local environment of atoms is thus detected from topological analysis based on geometrical criteria. Atom types are assigned according to the following nomenclature: Osih and Hsi are oxygen and hydrogen atoms of silanol groups, Ob corresponds to bridging oxygen atoms of silicate chains, while Osi are terminal oxygen atoms of inner and outer silicate chain sheets. Oh and Hh are oxygen and hydrogen atoms in hydroxide ions, while Ow and Hw are oxygen and hydrogen atoms of water molecules. Finally, Ca and Cw correspond to intra-laminar and inter-laminar calcium ions, respectively.

The parameters of the classical force field potentials used in the simulation of C-S-H systems are derived from CSH-FF [40] and ClayFF with the flexible SPC water model [25]. Detailed parameter values are available in Section S1 of the online Supplementary Information (SI). A geometric mixing rule is adopted to calculate the Lennard-Jones parameters between different atom types. The interactions between atoms are modeled taking into account the bonded harmonic contribution for H-O bonds and the angular contribution for H-O-H, allowing flexibility in the molecular structure. To account for different Ca/Si ratios, the number of atoms of each type vary from one system to another. After assignment of the initial charges, the excess charge is equally distributed among Ob, Osi and Osih atoms. This method was already employed in the literature to maintain the electroneutrality [12, 17, 34]. However, in contrast with the work of Duque-Redondo et al., the excess charge is only distributed on oxygen atoms of the skeleton and not on oxygen in hydroxides (Oh) to avoid charge unbalance in the bulk of the pore. The excess charge is increasing with the Ca/Si and H₂O/Si and is related to the relaxation with the ReaxFF and topological analysis. This results in different charge distribution within each atom types and thus in different charges for Ob, Osi and Osih atoms that raise up to +25% compared to the initial values for high Ca/Si ratios. Although the distribution of the charge does not affect the structural stability of the systems under study, one should note that this method may affect the interaction of chemical species within the pore. Hence, the formation of stronger bonds between cations and surface oxygen atoms is expected at high Ca/Si ratios. The modified partial charges are given in the Section S1 of the SI.

After the topological analysis and the application of the classical force field, the structures are relaxed in the isobaric NPT ensemble (constant number of particles, pressure, temperature) for a duration of 500 ps with a timestep of 1 fs. This step aims to equilibrate the system pressure and obtain an energetically stable initial configuration. An appropriate relaxation time is defined for the simulations to allow the system to reach an equilibrium state. This relaxation duration aims to ensure the convergence of the total energy,
thereby ensuring the global stability of the models. The stability of the geometry of the box during the different steps of creation of the model is ensured by tracking it (see Section S2 of the SI). Moreover, since the main goal of this study is to assess the diffusivity of chlorides ion in C-S-H, a sensitivity study is carried out to choose the water model and the mixing rules for Lennard-Jones parameters the most appropriate to obtain a diffusion coefficient close to experimental observations. Additional details can be found in the Section S3 of the SI.

Table 1 presents the chemical characteristics of the C-S-H models. The proportion of molecular water %<sub>mH2O</sub> and the proportion of chemically bonded water %<sub>bH2O</sub> is calculated as

$$ \%_{mH_2O} = \frac{N_{Ow}}{N_H/2} \cdot 100 \qquad (2) $$

$$ \%_{bH_2O} = \frac{N_H/2 - N_{Ow}}{N_H/2} \cdot 100 \qquad (3) $$

where $N_H$ is the total number of hydrogen atoms and $N_{Ow}$ is the total number of oxygen atoms in water molecules.

$Q^i$ is an indicator for silicate structures where $i$ is the connectivity, which refers to the number of silicates group in the surrounding of a given silicate group. $Q^0$ are thus monomers, $Q^1$ are dimers, $Q^2$ corresponds to pentamers, and $Q^3$ are silicate groups with three neighbors, that is bridging silicates between two chains. Here, the MCL is calculated approximately from the proportion of $Q^1$ and $Q^2$ silicates (see Eq. (4)) [41]:

$$ MCL \approx 2\left(1 + \frac{Q^1}{Q^2}\right) \qquad (4) $$

<table><caption>Table 1: Chemical properties of the C-S-H models</caption><thead><tr><th>Ca/Si</th><th>H<sub>2</sub>O/Si</th><th>ρ (g/cm<sup>3</sup>)</th><th>Molecular H<sub>2</sub>O (%)</th><th>Chemically bonded H<sub>2</sub>O (%)</th><th>Q<sup>0</sup> (%)</th><th>Q<sup>1</sup> (%)</th><th>Q<sup>2</sup> (%)</th><th>Q<sup>3</sup> (%)</th><th>MCL</th></tr></thead><tbody><tr><td>1.2</td><td>0.93</td><td>2.51</td><td>57.2</td><td>42.8</td><td>0.00</td><td>39.67</td><td>60.00</td><td>0.33</td><td>5.03</td></tr><tr><td>1.4</td><td>1.15</td><td>2.57</td><td>53.6</td><td>46.4</td><td>0.00</td><td>66.30</td><td>33.33</td><td>0.37</td><td>3.01</td></tr><tr><td>1.6</td><td>1.35</td><td>2.61</td><td>50.6</td><td>49.4</td><td>0.00</td><td>63.70</td><td>36.30</td><td>0.00</td><td>3.14</td></tr><tr><td>1.8</td><td>1.55</td><td>2.65</td><td>47.4</td><td>52.6</td><td>1.11</td><td>54.81</td><td>40.37</td><td>3.70</td><td>3.47</td></tr></tbody></table>

Because of the method employed to create the models, where the bridging groups are first removed, the MCL decreases significantly from Ca/Si=1.2 to 1.4. For Ca/Si=1.6 and 1.8, the MCL increases moderately because of the condensation (formation of $Q^2$ and $Q^3$ groups) during the relaxation step with ReaxFF. Marginally, few monomers are formed for Ca/Si=1.8. Nevertheless, the $Q^i$ and MCL values obtained from the molecular models fall within the range reported by nuclear magnetic resonance studies on synthesized C-S-H [42].

## 2.2. Creation of gel pore models and simulation details

The creation of the gel pore involves opening the C-S-H interlayer to establish a specific distance, allowing for the insertion of a solution box with varying widths: 5 Å, 10 Å, 20 Å, 30 Å, or 50 Å. NaCl solution boxes with concentrations of 0.5, 1.5, 2, 3, and 4 M are then introduced. One should note that numerical simulation on a representative elementary volume impose a minimal concentration in order to obtain statistically significant results while keeping the simulation box, and thus the simulation time small. For the size of the systems under study, a concentration lower than 0.5 M would lead to too few ions to observe and to ensure the sampling of physical properties. The step 5 of the Fig. 1 depicts C-S-H model with a Ca/Si ratio of 1.2 featuring a 30 Å nanopore containing a NaCl solution with a concentration of 0.5 M. The systems are relaxed in the isothermal NVT ensemble for a duration of 500 ps with a timestep of 0.5 fs to ensure temperature equilibration and energetically stable configuration of atoms. An equilibrium is performed in
the NPT ensemble for the same duration and timestep for the regulation of the pressure and adjustment of the simulation box volume. Once the system is equilibrated in terms of temperature and pressure, the production run is performed in the NVT ensemble for a duration of 30 ns with a timestep of 1 fs to collect statistically significant data. To determine the optimal simulation time and see its influence on the diffusion profile, a parametric study was conducted (see Section S4 of the online SI). After a thorough analysis of the data, it is concluded that a simulation time of 30 ns is most suitable to ensure reliable and representative results, and to keep a reasonable balance between accuracy and computational time. In all the simulations, the Verlet algorithm is used to integrate the equation of motion and a cutoff of 12 Å is used for the Lennard-Jones interactions. The Coulombic interactions are computed directly for pairwise interactions within this distance and the Ewald summation method is used to calculate interactions outside this distance in the reciprocal space.

The Table 2 summarizes all the simulations performed in this study. The full range of pore sizes was only simulated for Ca/Si=1.8. This choice is supported by experimental measurements and literature reviews, which indicate that although the molar Ca/Si ratio in the C-S-H phase of original Portland cement pastes typically ranges from approximately 1.4 to 2.0, values around 1.7 to 1.8 are generally considered the most representative [43–46]. It is worth noting that the inclusion of ground granulated blast furnace slag has been reported to lower the Ca/Si ratio, though it remains close to 1.8 at a 25% replacement level for both inner and outer product C-S-H [47].

<table><caption>Table 2: Simulated system. Code for concentrations:</caption><thead><tr><th>Pore Size<br>Ca/Si</th><th>5 Å</th><th>10 Å</th><th>20 Å</th><th>30 Å</th><th>50 Å</th></tr></thead><tbody><tr><td>1.2</td><td></td><td></td><td></td><td>♣♣♠</td><td></td></tr><tr><td>1.4</td><td></td><td></td><td></td><td>♣♣</td><td></td></tr><tr><td>1.6</td><td></td><td></td><td></td><td>♣♣</td><td></td></tr><tr><td>1.8</td><td>★</td><td>♣♣</td><td>♣♣</td><td>♣♣</td><td>♣♣</td></tr></tbody></table>

♣=0.5 M, ★=2 M, ♠=1.5 M, 3 M and 4 M.

# 3. Structure of the interface and ion adsorption

## 3.1. *Atomic density profiles*

Silicates in C-S-H can structure as the so-called drierketten chains with a bridging group pointing toward the interlayer region. With such a configuration, the solid surface, defined as the silicate chains linked by intralayer calcium, is not plane and the liquid and solid phase coexists within channels delimited by bridging silicates. For this reason, two interfaces can be defined:

* The outer interface which marks the end of the atomic density of the solid part of the C-S-H skeleton, including Si, Osi, Osi, Ob, Hsi and Ca atoms. This interface is defined as the end of the envelope of this set of atoms and marks the transition between the start of fluid phase and the interphase region, where the solid and liquid phases coexist.

* The inner interface, which marks the beginning of the solution phase that is composed of the following atoms: Ow, Hw, Oh, Hh, Cw, Cl, and Na. In other words, the inner interface defines the boundary where the atomic density of these atoms reaches zero, thereby indicating the transition between the interphase region and the solid region.
<figure><img src="image_3.png" /><figcaption>Figure 2: Atomic density profile for a gel pore of 30 Å. Subfigures a, b, c, and d (on the left) illustrate atomic densities for a concentration of 0.5 M, while subfigures e, f, g, and h (on the right) correspond to a concentration of 2 M. The left axis in each figure represents the overall atomic density of the solution and the solid (envelope of atomic density profiles), whereas the right axis provides a specific zoom on Cl, Na, and Cw atoms. The dashed blue lines indicate the inner interfaces while the dashed black lines indicate the outer interfaces. The origin (z = 0) is set at the left outer interface.</figcaption></figure>

The study of atomic density profiles (see Fig. 2) highlights the impact of the Ca/Si ratio on the structure of the interface. Since symmetry is not enforced in the C-S-H models, where bridging SiO₂ is removed randomly, the ion distribution within the pore is asymmetric. As reported for tobermorite in a pioneer molecular dynamics study [18], only cations can coordinate to the solid surface. While Cw and Na⁺ mainly adsorb chemically to the surface, Cl⁻ ions constitute a diffuse layer, together with part of the Na⁺. Such structure at the interface was already reported in previous MD investigations [22, 26, 31]. Cl⁻ ions thus adsorb onto the surface through interaction with cations. Two coordination peaks for Na-Cl and Cw-Cl pairs are reported on the radial distribution functions (RDF) available in the Section S5 of the online SI. The first peak is sharp and corresponds to the first coordination sphere at ~ 2.7 Å. A slight decrease of the coordination distance is observed with increasing Ca/Si ratio. The second peak is wide and corresponds to the second coordination sphere at ~ 4-6 Å. The overlapping of density profiles highlights a competition between Cw and Na⁺ ions which becomes more pronounced at high concentration and high Ca/Si ratio. Na⁺ can substitute for Cw ions at specific sites, particularly in calcium-rich environments where silicate sites are less available to stabilize the network. This gives rise to complex atomic structures, which may result
from interference between the hydration shells, making it harder to define the inner and outer coordination
spheres. However, Cw ions are much more strongly adsorbed to the silicates on the surface than Na⁺ ions
because of their higher charge, and no Cw ions were reported in the bulk of the pore. With an increasing
Ca/Si ratio, this adsorption becomes stronger and Cw closer to the surface. This behavior can be attributed
to a better accessibility to the pairing silicate with increasing Ca/Si ratio and leads to a tightening of the
Stern layer. To further illustrate the above, three peaks are reported on the RDF for Cw-Osi and Na-Osi
(see SI), at 2.3 Å (sharp), 3.5-5.5 Å and 6-8 Å (wide). The outer coordination spheres of Na-Osi is very
unstable compared to that of Cw-Osi, in particular at Ca/Si=1.6 where the second coordination sphere is
bi-modal.

For low Ca/Si ratios, a water layered structure is observed at the solid surface. This structure collapse
with increasing Ca/Si ratio, because the H₂O/Si ratio and the proportion of chemically bound water become
higher (see Table 1). This result suggests that the surface becomes more hydrophobic with the increase
in protonation of the surface silicates with the Ca/Si ratio. For this reason the H-bonds between the
water molecules and the silicates are weakened, resulting in the decay of the water H-bond network as the
Ca/Si ratio rises. These observations are consistent with simulations of C₃S/water interface with increasing
protonation of the surface [48]. Previous atomistic simulations of C-S-H gel pore with high concentration
of Cs ions (3.5-3.7 M) have shown that the water molecules were more structured at high Ca/Si ratio,
revealing an increase in hydrophilicity [17]. This observation was confirmed by a larger residence time of
water molecules on the surface and an upshifting of the dipole moment with increasing Ca/Si ratios. The
difference between this former study and the present one could indicate that the nature and concentration
of the gel pore solution influence how variations in the Ca/Si ratio affect the structure of water molecules at
the interface. The atomic density profiles as function of pore size for Ca/Si=1.8 are provided in the Section
S6 of the online SI. It is visible that the number of free ions is decreasing when the pore size decreases.
Furthermore, for pore sizes of 5 Å and 10 Å, the density profile of water shows a dip at the center of the
pore due to confinement effects.

3.2. *Adsorption*

The physicochemical adsorption of Cl⁻ and Na⁺ is evaluated based on the integration of the density
profile. Details on the method are provided in the Section S7 of the online SI. It should be noted that
both chemical and physical adsorption were considered hereafter without distinction. Although the density
profiles provide a good qualitative analysis of the EDL, they are not sufficiently clear to provide alone a
detailed quantified analysis of these two adsorption types. For this purpose, a geometrical-based analysis,
including pair distances or RDF as criteria could be used and may allow the characterization of inner-sphere
and outer-spheres surface complexes. The left-top part of Fig. 3 illustrates the fraction of bound ions as
a function of the Ca/Si ratio at 0.5 M and shows that more than 80% of Na⁺ and approximately 65-70%
of Cl⁻ are adsorbed onto the surface. At this concentration, a large fraction of ions are bound because
the number of available adsorption sites is high relative to the number of ions. At 2 M (see left-bottom of
Fig. 3), the fraction of bound ions decreases to approximately 60% due to the reduced availability of binding
sites.
<figure><img src="image_4.png" /><figcaption>Figure 3: Left: Fraction of bound chloride and sodium ions as a function of the Ca/Si ratio in a pore size of 30 Å (top: 0.5 M, bottom: 2 M). Right: Number of bound chloride and sodium ions per unit area as a function of the concentration (0.5 M; 1.5 M; 2 M; 3 M and 4 M) with a Ca/Si ratio of 1.2 and a pore size of 30 Å fitted to a Freundlich isotherm.</figcaption></figure>

A very slight increase in chloride adsorption with Ca/Si is observed at both 0.5 M and 2 M, which may be explained by the stronger attraction to the surface due to the increasing number of positive charges. The increase in the fraction of bound chloride is even smaller at 2 M compared to 0.5 M, but the absolute number of bound ions is higher. A previous study combining MD simulations and immersion experiments [49] also reported a higher chloride binding capacity of C-S-H at greater Ca/Si ratios. Using the same techniques, another study found that $Cl^-$ adsorption in C-S-H is greater when immersed in a $CaCl_2$ solution compared to NaCl [20], suggesting that calcium ions in the pore solution promote chloride binding. Similarly, a molecular investigation found that $Cl^-$ adsorption increases with the Ca/Si ratio due to stronger interactions between chlorides and calcium ions [27]. This was further supported by SEM-EDX image analysis [50] and chloride sorption isotherms [51]. Other experimental studies have also reported greater chloride binding capacity at higher Ca/Si ratios [52, 53], in agreement with our findings. However, a separate MD investigation found that chloride binding in C-S-H reaches a maximum at $Ca/Si = 1.2$ [26]. It is important to note that their numerical method differs from the present study, as it is based on transport models and considers the total amount of $Cl^-$ migrating from a bulk solution to a C-S-H gel pore. The results of that study align with chemical analyses of synthesized C-S-H immersed in a $CaCl_2$ solution [54]. This experimental study also reported a maximum chloride adsorption at $Ca/Si = 1.2$ for NaCl solution. However, the variation in bound chloride per unit mass of C-S-H with the Ca/Si ratio was found to be minimal. Regarding sodium adsorption, a drop is observed at $Ca/Si = 1.6$ at both 0.5 M and 2 M (see Fig. 3, left). This decrease may be attributed to the instability of the outer coordination, as indicated by the RDF, particularly at $Ca/Si = 1.6$. Indeed, due to their higher charge, Cw ions are strongly bonded to oxygen in silicates. In contrast, the interactions of $Na^+$ ions with the surface oxygen atoms are less stable and they can swap from the Stern layer to the diffuse layer. This instability of the second coordination sphere is visible in the RDF of Na-Osi pairs (see Section S5 of the online SI). A well-defined bimodal second coordination sphere is apparent for $Ca/Si = 1.6$, which might explain the drop in sodium adsorption observed at $Ca/Si = 1.6$ in Fig. 3 (left). Discrepancies between the results and conclusions of different numerical investigations are likely related to variations in simulation methods: some considering transport, while others calculate adsorption at equilibrium using atomic density profiles. Furthermore, differences in the methods used to create the C-S-H model and the choice of force field could impact the results. Also importantly, the methodology to determine the quantity of bound ions can be further improved with distance-based criteria based on the RDF to discriminate between chemically adsorbed, physically adsorbed and free ions. The right part of Fig. 3 shows that the total amount of bound ions follows a Freundlich isotherm and increases with the total concentration in the pore solution.
A gradual saturation of available binding sites is observed for Na⁺. However, this saturation effect is less pronounced for Cl⁻ within the concentration range under study. Moreover, the pore size affects significantly the proportion of bound species (see Section S8 of the online SI). The number of free ions decreases markedly in small pores due to confinement effects and the overlap of EDL. Across all pore sizes, more Na⁺ ions are bound than Cl⁻, reflecting their respective positions in the EDL: Na⁺ mainly accumulates in the Stern layer, whereas Cl⁻ remains in the diffuse layer.

## 4. Bond lifetime

The stability of Cw-Osi, Cw-Cl, Na-Osi, and Na-Cl bonds are analyzed in Fig. 4 using the following time autocorrelation function [55]:

$$C(t) = \left\langle \frac{\sum h_{ij}(t_0)h_{ij}(t_0 + t)}{\sum h_{ij}(t_0)^2} \right\rangle \quad (5)$$

where $h_{ij} = 1$ if the bond exists and $h_{ij} = 0$ if the bond is broken.

<figure><img src="image_5.png" /><figcaption>Figure 4: Analysis of the time autocorrelation function for interactions between Cw-Osi, Cw-Cl, Na-Osi, and Na-Cl pairs for a concentration of 2 M and a pore size of 30 Å</figcaption></figure>

The critical distance at which the bond is broken is determined from the first minima of the RDF (see SI). A continuous description is used, meaning the number of existing bonds is computed at the start of the correlation time and monitored throughout, without the possibility of bond reformation – once broken, a bond remains so. In order to facilitate the analysis and avoid notched curved, the autocorrelation functions are fitted with a sum of two exponential decays:

$$C(t) \approx A_1 \exp\left(-\frac{t}{\tau_1}\right) + (1 - A_1) \exp\left(-\frac{t}{\tau_2}\right). \quad (6)$$
From Fig. 4, it can be observed that Cw-Osi and Cw-Cl bonds are significantly more stable than Na-Osi and Na-Cl bonds. This can be explained by the higher charge of Cw ions, which makes their adsorption to the surface and their bonds with Cl⁻ much stronger than those of Na⁺. The relaxation time of the Cw-Osi increased for Ca/Si=1.8 and that of Na-Osi bonds increases for Ca/Si=1.6 and Ca/Si=1.8, reflecting greater stability of cations-surface interactions under high Ca/Si ratio. It may be counter-intuitive since a saturation effect is expected and may weaken this interaction. However, the removal of bridging SiO₂ may have increased the accessibility of cations to more adsorption sites and strengthened their binding to silicates. In contrast, the relaxation time for Cw-Cl bonds is shorter for high Ca/Si ratio, indicating that the Cw-Cl are weakening with the increase of Ca/Si. In particular, the relaxation time for Ca/Si = 1.2 is significantly longer than for other Ca/Si ratio. For the Na-Cl pair, no significant influence of the Ca/Si ratio on the relaxation time is observed, suggesting that this interaction remains relatively independent on the Ca/Si ratio. However, this trend is observed for the whole set of atoms presents in the pore but may be different if only the ions at the interface were observed. A trend similar to that of the Cw-Cl pair could be expected, but with higher bond flexibility.

## 5. Diffusivity of chloride and sodium ions

### 5.1. Diffusion regime and anisotropy

The diffusion coefficient of Cl⁻ and Na⁺ ions are calculated from the Einstein equation (Eq. (7)), which establishes a fundamental relationship between particle diffusion and their mobility based on the temperature and viscosity of the medium

$$D = \frac{1}{2d} \lim_{\tau \to \infty} \frac{MSD}{\tau} \qquad (7)$$

where d = 3 is the dimensionality of the system. The mean squared displacement (MSD), is calculated as:

$$MSD = \frac{1}{N} \sum_{i}^{N} |\mathbf{r}_{i}(t_{0} + \tau) - \mathbf{r}_{i}(t_{0})|^{2} \qquad (8)$$

where $\mathbf{r}_i$ is the position vector of the particle i, N is the number of particles in the block, and $t_0$ denotes the initial time of a given observation window of duration $\tau$. Herein, the diffusion coefficient is plotted in the gel pore as a function of the z coordinate (perpendicular to the surface). For this purpose, for each observation window, the MSD and the averaged z coordinate of each ions are recorded. Finally, the MSDs are grouped and averaged into spatial blocks. In practice, the MSD should be observed over a time $\tau$ sufficiently long to capture a stable diffusion regime. However, this time must be short enough to maintain the local nature of the measurement in each spatial block [56, 57]. In this study, a time $\tau = 30$ ps is chosen following a sensitivity analysis (see Section S4 of the SI). The plot of MSD over time reveals two distinct regimes independent of the distance from the surface (Fig. 5). The first regime, lasting less than a picosecond, is characterized by a super-diffusive, ballistic motion, as previously reported in the literature [58–60]. In contrast, the second regime corresponds to Fickian diffusion ($MSD \propto t$) for ions at the center of the pore, while ions closer to the surface exhibit a subdiffusive regime ($MSD \propto t^{\alpha}$, $\alpha < 1$). In this study, diffusion was calculated within the second regime (from 1 ps) using Einstein's equation, thereby assuming Fickian behavior to obtain a diffusion coefficient in square meters per second. A generalized self-diffusion approach can be employed for a deeper analysis of anomalous diffusion [61, 62].
<figure><img src="image_6.png" /><figcaption>Figure 5: MSD of chloride ions along the *y* direction as a function of the *z*-coordinate within a 30 Å nanopore in a C-S-H model with a Ca/Si ratio of 1.8 and a molar concentration of 0.5 M. The dashed line indicates the left outer interface, where the origin (*z* = 0) is set.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6: MSD of chloride ions in a C-S-H model with a Ca/Si ratio of 1.8, a pore size of 30 Å, and a concentration of 0.5 M, in the three directions (*x*, *y*, *z*), for ions close to the surface (left) and ions far from the surface (right).</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 7: Averaged, parallel, and perpendicular diffusion coefficient profiles for a C-S-H model with a Ca/Si ratio of 1.2, a pore size of 30 Å, and a chloride ion concentration of 0.5 M. The dashed line indicates the left outer interface, where the origin (z = 0) is set.</figcaption></figure>

The translational motion in the three directions for chloride ions close and far from the surface is depicted in Fig. 6. The *x* and *y* coordinates are parallel to the surface, while the *z* coordinate is perpendicular to the surface. An anisotropic motion is observed for ions close to the surface (*MSD*<sub>xx</sub> > *MSD*<sub>yy</sub> > *MSD*<sub>zz</sub>) where the MSD exhibits a pseudo-plateau behavior with a steep step at the very beginning of the observation window, followed by a slow evolution, which characterize almost motionless ions (see Fig. 6 on the left). As expected and already reported for particle motion between walls, the motion perpendicularly to the surface is smaller than the motion parallel to the surface [59]. Interestingly, the motion of Cl<sup>-</sup> in the *x* direction is by order of twice the motion in the *y* direction. This difference can be attributed to the topology and charge distribution of the atoms of the surface. In particular the direction of the silicate chains and may influence the difference in mobility observed between the *x* and *y* directions. At the center of the pore, the ions adopt the behavior of the bulk and the motion can be described as isotropic (*MSD*<sub>xx</sub> ≈ *MSD*<sub>yy</sub> ≈ *MSD*<sub>zz</sub>) (see right of Fig. 6).

The Fig. 7 represents the diffusion profile coefficient of Cl⁻ within the pore. It illustrates particularly the difference between the averaged diffusion coefficient over all direction x, y and z D<sub>avg</sub> (or D<sub>iso</sub>), the in-plane or parallel D<sub>‖</sub> (along x and y) and the perpendicular diffusion coefficient D<sub>⊥</sub> (along z). As expected from Fig. 6, D<sub>⊥</sub> is smaller than D<sub>‖</sub>, in particular at the vicinity of the surface where ions diffuse in an anisotropic manner. However, at the scale of the diffusion profile, there is no significant difference between D<sub>iso</sub> and D<sub>‖</sub>. The error on the diffusion coefficient is represented as a bright-colored area. It is calculated with a 95% confidence interval from a weighted linear regression on the MSDs and their standard deviation, where the weight is the inverse of the variance (see details in the Section S9 of the SI). In what follows, the results are only provided in terms of averaged diffusion coefficient D<sub>iso</sub>.

## 5.2. Influence of pore size, concentration, and Ca/Si ratio

The Fig. 8 presents the averaged diffusion coefficients of chloride ions at 300 K as a function of distance from the C-S-H outer interface.
<figure><img src="image_9.png" /><figcaption>Figure 8: Influence of the pore size on the diffusion profile of chloride ions. Left: Ca/Si = 1.8 at 0.5 M, right: Ca/Si = 1.8 at 2 M. The dashed lines indicate the left outer interface, where the origin (z = 0) is set.</figcaption></figure>

At 0.5 M, for nanopores of 20 Å and 30 Å, the maximum diffusion coefficient approaches that of Cl-ions in an unconfined solution, which has been measured at $2.03 \times 10^{-9} \text{ m}^2\text{s}^{-1}$ in seawater at 298 K at infinite dilution [63, 64]. The low value of the diffusion coefficient near the interface can be explained by the effect of the EDL, which hinders the free movement of ions. The range of influence of the surface is about 12 Å. Hence, for nanopores smaller than twice this distance, the bulk diffusion is never reached. For example, the maximum diffusion coefficient for a nanopore size of 10 Å is $9.2 \times 10^{-10} \text{ m}^2\text{s}^{-1}$. In such pores, the effect of both surface are coupled and the EDL overlaps. The mobility of ions may also be affected by the overlapping of hydration shells, especially at high concentration [58, 65, 66]. The influence of these effects on the diffusion profile can be further observed in the Section S10 of the online SI. As expected the concentration as a significant impact on the diffusion coefficient (see Fig. 9).

<figure><img src="image_10.png" /><figcaption>Figure 9: Influence of concentration on the diffusion profile of chloride and sodium ions for a pore size of 30 Å and a Ca/Si ratio of 1.8. a. for chloride ions, b. for sodium ions, c. variation of the diffusion coefficient at the center of the pore as a function of the concentration (red dotted line in a and b). The results at the center are fitted to a Kohlrausch law in plain line for our simulation and in dashed line for a previous experimental study. The dashed lines indicate the left outer interface, where the origin (z = 0) is set.</figcaption></figure>

The relationship between the diffusion coefficient *D* and the concentration *C* at the center of the pore
are fitted to a Kohlrausch law:

$$D/D_0 = 1 - A\sqrt{C} \qquad (9)$$

where $D_0$ is the diffusion coefficient in the bulk solution at infinite dissolution and A is an empirical constant that depends on the interionic interactions. The decrease in the diffusion coefficient with the concentration for NaCl solution was already reported in both MD studies of bulk solutions [67–69] and experimental investigations [70]. All the aforementioned MD studies underestimate the experimental values. The results of our simulations align closely with the experimental findings. However, based on the fits to the Kohlrausch law, the simulated coefficients are slightly overestimated at low concentrations (1.5-2 M) and underestimated at higher concentrations. Our results at 0.5 M are much closer to the experimental data at infinite dilution than the prediction of the Kohlrausch law (~2.5 × 10⁻⁹ m² s⁻¹ at 0 M for Cl⁻). According to the experimental measurements by Vitagliano and Lyons [70], the relationship between the diffusion coefficient and concentration for NaCl solution is approximately linear, suggesting that the Kohlrausch law may not be the most appropriate model to describe this behavior, especially close to 0 M. Globally, no major influence of the Ca/Si ratio is reported for a concentration of 0.5 M, except for a small shift in the diffusion coefficient profile at the vicinity of the surface (see Fig. 10).

The Fig. 11 examines closely the MSD of chloride ions at the vicinity of the surface. The spatial blocks correspond to the first three regions (each a 1 Å slice perpendicular to the surface), starting from the point where chloride ions are first detected. These regions correspond to the first three points of the diffusion profile in Fig. 10. In the first spatial block (solid lines), the MSD for Ca/Si=1.8 is much higher than for other ratios. In the second spatial block (dashed lines), 1 Å farther, the MSD increase gradually with the Ca/Si ratio. From the third spatial block (dotted lines), the Ca/Si ratio is no more influencing the MSD. This detailed analysis shows that, although the Ca/Si ratio has no impact on the overall diffusion coefficient profile, the diffusion of ions within the first ångströms increases with the Ca/Si ratio. At 2 M, the Fig. 10 (left) shows a slight increase of the coefficient diffusion of Cl⁻ at the center of the pore for Ca/Si=1.8. This might be explained by the weakening of Cw-Cl interactions at high Ca/Si ratios observed on Fig. 4. Additionally, the error margin becomes smoother with an increase in concentration, thus improving the sampling.

<figure><img src="image_11.png" /><figcaption>Figure 10: Influence of the Ca/Si ratio on the diffusion profile of chloride ions for a pore size of 30 Å. Left: 0.5 M, right: 2 M. The dashed lines indicate the left outer interface, where the origin (z = 0) is set.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>Figure 11: MSD of chloride ions in the three closest spatial blocks from the C-S-H surface at 0.5 M and for a pore size of 30 Å.</figcaption></figure>

In Fig. 12, the diffusion profiles of sodium ions are influenced by their strong adsorption to the solid surface, leading to a more restricted diffusion coefficient near the pore surface compared to $Cl^-$ . The results are similar but may overestimate the experimental value obtained at infinite dilution (1.33 × 10⁻⁹ m² s⁻¹) [63, 70].

<figure><img src="image_13.png" /><figcaption>Figure 12: Influence of the Ca/Si ratio of a pore size of 30 Å (left) and of the pore size (right) on the diffusion profiles of sodium ions at a concentration of 0.5 M. The dashed lines indicate the left outer interface, where the origin (z = 0) is set.</figcaption></figure>

$Na^+$ are less mobile than $Cl^-$ and very attracted by the surfaces, for this reason the diffusion coefficient could not be sampled in some places of the pore. This is particularly the case for pores of 5 Å where
the number of ions is low and their free movement is particularly restricted by the confinement and the overlapping of EDL (see Section S10 of the SI).

## 6. Conclusions and limitations

In this study, a molecular simulation approach was used to assess the diffusion and adsorption of Cl⁻ and Na⁺ ions in C-S-H gel pore models with varying Ca/Si ratios and different pore sizes. The models were created based on a 11 Å tobermorite crystal structure, and a NaCl solution was introduced, considering various molar concentrations (0.5 to 4 M). Through several MD simulations, the main conclusions that can be drawn are the following:

(a) The proposed C-S-H model is based on MCL of silicate chains in good agreement with experimental measurements across different Ca/Si ratios. The current study proposed a complete overview of the diffusion mechanism including a discussion on the diffusion regimes within C-S-H gel pores. More importantly, the current study proposed extensive results and discussion on the change of the diffusivity of chloride and sodium atoms with the distance from the pore surface.

(b) The structure of the EDL was highlighted, involving cations bound to the oxygen atoms in silicates in the Stern layer and chloride ions in the diffuse layer with part of sodium ions. The former tend to compete with calcium ions at the interface. However, Na-Osi bonds exhibit shorter lifetimes compared to the more polar Cw-Osi interaction.

(c) The adsorption isotherm of Na⁺ across the considered concentration range (0.5 to 4 M) indicates a gradual saturation of the binding sites at the surface. This can be explained by the competition between cations in the Stern layer, which favors Ca²⁺ over Na⁺ due to the much higher stability of the Cw-Osi bond.

(d) The proportion of bound Na⁺ is much higher than that of Cl⁻, and it does not follow a particular trend as a function of the Ca/Si ratio. In contrast, a slight increase in the binding capacity with the Ca/Si ratio was observed for Cl⁻ at 0.5 M.

(e) The study of ionic mobility in C-S-H: gel pores highlighted an anisotropic and sub-diffusive diffusion regime near the surface. The diffusion becomes isotropic and eventually reaches a Fickian regime at 12 Å from the surface, in sufficiently large gel pores. The pore size plays a critical role in governing the diffusivity and adsorption of Na⁺ and Cl⁻ ions through confinement effects and the overlapping of the EDL.

(f) The overall diffusion coefficient profile is not significantly affected by the Ca/Si ratio. However, an increase in diffusivity near the surface was observed with increasing Ca/Si ratio, particularly for Cl⁻ ions within the first coordination sphere of cations adsorbed on silicate chains.

The simulated diffusion coefficients at the center of the pore were validated against experimental measurements. However, the results may be significantly affected by the choice of force field, which influences ion mobility and their interactions with the pore surface. Specifically, the distribution of excess charge over the oxygen atoms of the C-S-H framework may lead to an overestimation of cation binding to silicate sites, especially at high Ca/Si ratios. Additionally, the force field employed in this study is non-polarizable and does not explicitly account for electronic polarization effects, which are known to impact ion-surface interactions. In particular, chloride ions, being relatively large and highly polarizable, may exhibit stronger binding to positively charged surface sites in reality than predicted by the simulations. This limitation could result in an underestimation of chloride surface affinity or changes in the local hydration structure and stability near the interface, thereby potentially affecting the diffusivity. Moreover, more sophisticated C-S-H models can be developed by incorporating additional experimental data, thereby enhancing the reliability of their predicted structural and chemical properties. [15] It is also clear that this initial work could be improved by considering a multi-species pore solution that is more representative of the real conditions of structures
exposed to seawater containing other ions than Na⁺ and Cl⁻ (SO₄²⁻, K⁺, Mg⁺...). Furthermore, blended cements with ground granulated blast furnace slag and metakaolin are known to refine the pore structure during hydration due to the pozzolanic reaction and to enhance chloride binding capacity due to the increase in alumina content [71, 72]. Future investigations could consider C-A-S-H to clarify the impact of the Al/Si ratio on the diffusion and sorption of chloride ions at the molecular scale.

## 7. Outlook

This work lays a foundation for further exploration of diffusion profiles and a deeper understanding of ion adsorption phenomena, with the ultimate goal of integrating these insights into multiscale homogenization models. These models will facilitate the determination of diffusion coefficients at the cement paste scale and eventually concrete. Multiscale model outputs can subsequently be implemented in predictive models for the durability of reinforced concrete structures. A forthcoming study will utilize the diffusion profiles as input for multiscale modeling. The present MD results have important implications for the development of such a model. The position-dependent diffusion coefficient profiles clearly indicate that the reduction in diffusivity occurs within the first 12 Å from the surface and is significant only in the first few ångströms. Therefore, in multiscale models, a distinction between capillary and gel pores is critical. In particular, the diffusion coefficient in capillary pores should be identified with the diffusion coefficient in the bulk $D_{bulk}$, which suggests a revision of the assumptions underlying some existing models [1, 5, 73]. Further, modeling the diffusion coefficient in gel pores with a uniform value much smaller than in the bulk, for example 0.025$D_{bulk}$ [4] is incompatible with the present results. The diffusion profiles confirm the effective diffusion coefficient of gel pores depends on pore size [3, 59]. Yet, the effective diffusion coefficients in C-S-H gels reported in the literature for water are one to two order of magnitude lower than $D_{bulk}$ [59, 74–76]. The present MD diffusion profiles thus suggest that the diffusion is controlled by interfaces or constrictions rather than the bulk of gel pores. Finally, the anisotropy of the diffusion in slit gel pores near interfaces can be accounted for [59], although the improvement over an isotropic model is not expected to be of leading order [5].

## 8. Acknowledgements

Financial support from the French Ministère de l'enseignement supérieur et de la recherche for the Ph.D. thesis of Tarek Ihaddadene is gratefully acknowledged.

For the purpose of Open Access, a CC-BY public copyright licence has been applied by the authors to the present document and will be applied to all subsequent versions up to the Author Accepted Manuscript arising from this submission.
References

[1] M. Zhang, D. Lin, Z. He, R. Yang, Modeling the ionic diffusion coefficient of unsaturated hardened cement paste: A micromechanical approach, Cement and Concrete Research 177 (2024) 107415. doi: 10.1016/j.cemconres.2023.107415.

[2] L.-y. Tong, Q. Xiang Xiong, M. Zhang, Z. Meng, F. Meftah, Q.-f. Liu, Multi-scale modelling and statistical analysis of heterogeneous characteristics effect on chloride transport properties in concrete, Construction and Building Materials 367 (2023) 130096. doi:10.1016/j.conbuildmat.2022.130096.

[3] C. Liu, M. Zhang, Microstructure-based modelling of chloride diffusivity in non-saturated cement paste accounting for capillary and gel pores, Cement and Concrete Research 168 (2023) 107153. doi:10.1016/j.cemconres.2023.107153.

[4] M. Achour, F. Bignonnet, J.-F. Barthélémy, E. Rozière, O. Amiri, Multi-scale modeling of the chloride diffusivity and the elasticity of Portland cement paste, Construction and Building Materials 234 (2020) 117124. doi:10.1016/j.conbuildmat.2019.117124.

[5] N. Damrongwiriyanupap, S. Scheiner, B. Pichler, C. Hellmich, Self-Consistent Channel Approach for Upscaling Chloride Diffusivity in Cement Pastes, Transport in Porous Media 118 (3) (2017) 495-518. doi:10.1007/s11242-017-0867-3.

[6] R. A. Patel, Q. T. Phung, S. C. Seetharam, J. Perko, D. Jacques, N. Maes, G. De Schutter, G. Ye, K. Van Breugel, Diffusivity of saturated ordinary Portland cement-based materials: A critical review of experimental and analytical modelling approaches, Cement and Concrete Research 90 (2016) 52-72. doi:10.1016/j.cemconres.2016.09.015.

[7] H. Friedmann, O. Amiri, A. Aït-Mokhtar, Shortcomings of geometrical approach in multi-species mod-
elling of chloride migration in cement-based materials, Magazine of Concrete Research 60 (2) (2008)
119-124. doi:10.1680/macr.2008.60.2.119.

[8] P. Faucon, J. Delaye, J. Virlet, Molecular Dynamics Simulation of the Structure of Calcium Silicate Hydrates, Journal of Solid State Chemistry 127 (1) (1996) 92-97. doi:10.1006/jssc.1996.0361.

[9] P. Faucon, J. M. Delaye, J. Virlet, J. F. Jacquinot, F. Adenot, Study of the structural properties of the C-S-H(I) BY molecular dynamics simulation, Cement and Concrete Research 27 (10) (1997) 1581-1590. doi:10.1016/S0008-8846(97)00161-0.

[10] S. A. Hamid, The crystal structure of the 11 Å natural tobermorite Ca₂.₂₅[Si₃ O 7.₅ (OH) 1.₅]·1H₂O, Zeitschrift für Kristallographie - Crystalline Materials 154 (1-4) (1981) 189-198. doi:10.1524/ zkri.1981.154.14.189.

[11] R. J.-M. Pellenq, A. Kushima, R. Shahsavari, K. J. V. Vliet, M. J. Buehler, S. Yip, F.-J. Ulm, A realistic molecular model of cement hydrates, Proceedings of the National Academy of Sciences 106 (38) (2009) 16102-16107. doi:10.1073/pnas.0902180106.

[12] M. J. A. Qomi, K. J. Krakowiak, M. Bauchy, K. L. Stewart, R. Shahsavari, D. Jagannathan, D. B. Brommer, A. Baronnet, M. J. Buehler, S. Yip, F.-J. Ulm, K. J. V. Vliet, R. J.-M. Pellenq, Combinatorial molecular optimization of cement hydrates, Nature Communications 5 (1) (2014) 4960. doi:10.1038/ncomms5960.

[13] G. Kovačević, B. Persson, L. Nicoleau, A. Nonat, V. Veryazov, Atomistic modeling of crystal structure of Ca1.67SiHx, Cement and Concrete Research 67 (2015) 197-203. doi:10.1016/j.cemconres.2014.09.003.
[14] A. Kunhi Mohamed, S. C. Parker, P. Bowen, S. Galmarini, An atomistic building block description of C-S-H - Towards a realistic C-S-H model, Cement and Concrete Research 107 (2018) 221-235. doi:10.1016/j.cemconres.2018.01.007.

[15] Z. Casar, J. López-Zorrilla, H. Manzano, E. Duque-Redondo, A. K. Mohamed, K. Scrivener, P. Bowen, pyCSH: Automated atomic-level structure generation of bulk C-S-H and investigation of their intrinsic properties, Cement and Concrete Research 183 (2024) 107593. doi:10.1016/j.cemconres.2024.107593.

[16] M. Bauchy, M. J. A. Qomi, C. Bichara, F.-J. Ulm, R. J.-M. Pellenq, Rigidity transition in materials: Hardness is driven by weak atomic constraints, Physical Review Letters 114 (12) (Mar. 2015). doi:10.1103/physrevlett.114.125502.

[17] E. Duque-Redondo, K. Yamada, H. Manzano, Cs retention and diffusion in C-S-H at different Ca/Si ratio, Cement and Concrete Research 140 (2021) 106294. doi:10.1016/j.cemconres.2020.106294.

[18] A. G. Kalinichev, R. J. Kirkpatrick, Molecular dynamics modeling of chloride binding to the surfaces of calcium hydroxide, hydrated calcium aluminate, and calcium silicate phases, Chemistry of Materials 14 (8) (2002) 3539-3549. doi:10.1021/cm0107070.

[19] T. Pan, K. Xia, L. Wang, Chloride binding to calcium silicate hydrates (C-S-H) in cement paste: A molecular dynamics analysis, International Journal of Pavement Engineering 11 (5) (2010) 367-379. doi:10.1080/10298436.2010.488732.

[20] Y. Zhou, D. Hou, J. Jiang, P. Wang, Chloride ions transport and adsorption in the nano-pores of silicate calcium hydrate: Experimental and molecular dynamics studies, Construction and Building Materials 126 (2016) 991-1001. doi:10.1016/j.conbuildmat.2016.09.110.

[21] B. Zehtab, A. Tarighat, Diffusion study for chloride ions and water molecules in C-S-H gel in nano-scale using molecular dynamics: Case study of tobermorite, Advances in concrete construction 4 (4) (2016) 305-317. doi:10.12989/ACC.2016.4.4.305.

[22] D. Hou, Z. Li, Molecular dynamics study of water and ions transport in nano-pore of layered structure: A case study of tobermorite, Microporous and Mesoporous Materials 195 (2014) 9-20. doi:10.1016/j.micromeso.2014.04.011.

[23] P. Wang, Q. Zhang, M. Wang, B. Yin, D. Hou, Y. Zhang, Atomistic insights into cesium chloride solution transport through the ultra-confined calcium-silicate-hydrate channel, Physical Chemistry Chemical Physics 21 (22) (2019) 11892-11902. doi:10.1039/C8CP07676F.

[24] D. Hou, T. Li, P. Wang, Molecular Dynamics Study on the Structure and Dynamics of NaCl Solution Transport in the Nanometer Channel of CASH Gel, ACS Sustainable Chemistry & Engineering 6 (7) (2018) 9498-9509. doi:10.1021/acssuschemeng.8b02126.

[25] R. T. Cygan, J.-J. Liang, A. G. Kalinichev, Molecular models of hydroxide, oxyhydroxide, and clay phases and the development of a general force field, The Journal of Physical Chemistry B 108 (4) (2004) 1255-1266. doi:10.1021/jp0363287.

[26] Z. Liu, D. Xu, S. Gao, Y. Zhang, J. Jiang, Assessing the Adsorption and Diffusion Behavior of Multicomponent Ions in Saturated Calcium Silicate Hydrate Gel Pores Using Molecular Dynamics, ACS Sustainable Chemistry & Engineering 8 (9) (2020) 3718-3727. doi:10.1021/acssuschemeng.9b06817.

[27] N. Wang, R. Zhao, L. Zhang, X. Guan, Molecular insights into the adsorption of chloride ions in calcium silicate hydrate gels: The synergistic effect of calcium to silicon ratio and sulfate ion, Microporous and Mesoporous Materials 345 (2022) 112248. doi:10.1016/j.micromeso.2022.112248.
[28] P. A. Bonnaud, Q. Ji, B. Coasne, R. J.-M. Pellenq, K. J. V. Vliet, Thermodynamics of water confined in porous calcium-silicate-hydrates, Langmuir: the ACS journal of surfaces and colloids 28 (31) (2012) 11422-11432. doi:10.1021/la301738p.

[29] H. Deng, Z. He, Interactions of sodium chloride solution and calcium silicate hydrate with different calcium to silicon ratios: A molecular dynamics study, Construction and Building Materials 268 (2021) 121067. doi:10.1016/j.conbuildmat.2020.121067.

[30] W. Zhang, D. Ding, M. Li, T. Wang, H. Ma, B. Chen, H. Hu, J. Chen, X. Liu, D. Hou, Chloride binding mechanism in seawater-mixed UHPC, Construction and Building Materials 427 (2024) 136191. doi:10.1016/j.conbuildmat.2024.136191.

[31] E. Duque-Redondo, P. A. Bonnaud, H. Manzano, A comprehensive review of C-S-H empirical and computational models, their applications, and practical aspects, Cement and Concrete Research 156 (2022) 106784. doi:10.1016/j.cemconres.2022.106784.

[32] I. G. Richardson, Model structures for C-(A)-S-H(I), Acta Crystallographica Section B: Structural Science, Crystal Engineering and Materials 70 (6) (2014) 903-923. doi:10.1107/S2052520614021982.

[33] L. Martínez, R. Andrade, E. G. Birgin, J. M. Martínez, P ACKMOL: A package for building initial configurations for molecular dynamics simulations, Journal of Computational Chemistry 30 (13) (2009) 2157-2164. doi:10.1002/jcc.21224.

[34] E. Duque-Redondo, Y. Kazuo, I. López-Arbeloa, H. Manzano, Cs-137 immobilization in C-S-H gel nanopores, Physical Chemistry Chemical Physics 20 (14) (2018) 9289-9297. doi:10.1039/c8cp00654g.

[35] A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown, P. S. Crozier, P. J. in 't Veld, A. Kohlmeyer, S. G. Moore, T. D. Nguyen, R. Shan, M. J. Stevens, J. Tranchida, C. Trott, S. J. Plimpton, LAMMPS - a flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales, Computer Physics Communications 271 (2022) 108171. doi: 10.1016/j.cpc.2021.108171.

[36] T. P. Senftle, S. Hong, M. M. Islam, S. B. Kylasa, Y. Zheng, Y. K. Shin, C. Junkermeier, R. Engel-Herbert, M. J. Janik, H. M. Aktulga, T. Verstraelen, A. Grama, A. C. T. van Duin, The ReaxFF reactive force-field: Development, applications and future directions, npj Computational Materials 2 (2016). doi:10.1038/npjcompumats.2015.11.

[37] J. C. Fogarty, H. M. Aktulga, A. Y. Grama, A. C. T. van Duin, S. A. Pandit, A reactive molecular dynamics simulation of the silica-water interface, The Journal of Chemical Physics 132 (17) (2010) 174704. doi:10.1063/1.3407433.

[38] H. Manzano, R. J. M. Pellenq, F.-J. Ulm, M. J. Buehler, A. C. T. van Duin, Hydration of calcium oxide surface predicted by reactive force field molecular dynamics, Langmuir: the ACS journal of surfaces and colloids 28 (9) (2012) 4187-4197. doi:10.1021/la204338m.

[39] H. Manzano, S. Moeini, F. Marinelli, A. C. T. van Duin, F.-J. Ulm, R. J.-M. Pellenq, Con- fined water dissociation in microporous defective silicates: Mechanism, dipole distribution, and im- pact on substrate properties, Journal of the American Chemical Society 134 (4) (2012) 2208-2215. doi:10.1021/ja209152n.

[40] R. Shahsavari, R. J.-M. Pellenq, F.-J. Ulm, Empirical force fields for complex hydrated calcio- silicate layered materials, Physical Chemistry Chemical Physics 13 (3) (2011) 1002-1011. doi: 10.1039/c0cp00516a.
[41] I. Richardson, Tobermorite/jennite- and tobermorite/calcium hydroxide-based models for the structure of C-S-H: Applicability to hardened pastes of tricalcium silicate, -dicalcium silicate, Portland cement, and blends of Portland cement with blast-furnace slag, metakaolin, or silica fume, Cement and Concrete Research 34 (9) (2004) 1733-1777. doi:10.1016/j.cemconres.2004.05.034.

[42] X. Cong, R. Kirkpatrick, 29Si MAS NMR Study of the Structure of Calcium Silicate Hydrate, Advanced Cement Based Materials 3 (3-4) (1996) 144-156. doi:10.1016/s1065-7355(96)90046-2.

[43] H. Taylor, Cement Chemistry, 2nd Edition, Thomas Telford Publishing, 1997. doi:10.1680/cc.25929.

[44] J. Beaudoin, I. Odler, Hydration, Setting and Hardening of Portland Cement, in: Lea's Chemistry of Cement and Concrete, Elsevier, 2019, pp. 157-250. doi:10.1016/B978-0-08-100773-0.00005-8.

[45] A. J. Allen, J. J. Thomas, H. M. Jennings, Composition and density of nanoscale calcium-silicate-hydrate in cement, Nature Materials 6 (4) (2007) 311-316. doi:10.1038/nmat1871.

[46] I. Richardson, The nature of C-S-H in hardened cements, Cement and Concrete Research 29 (8) (1999) 1131-1147. doi:10.1016/S0008-8846(99)00168-4.

[47] I. G. Richardson, G. W. Groves, Microstructure and microanalysis of hardened cement pastes involving ground granulated blast-furnace slag, Journal of Materials Science 27 (22) (1992) 6204-6212. doi: 10.1007/BF01133772.

[48] J. Claverie, F. Bernard, J. M. M. Cordeiro, S. Kamali-Bernard, Waters behaviour on Ca-rich tricalcium silicate surfaces for various degrees of hydration: A molecular dynamics investigation, Journal of Physics and Chemistry of Solids 132 (2019) 48-55. doi:10.1016/j.jpcs.2019.03.020.

[49] Y. Zhou, D. Hou, J. Jiang, L. Liu, W. She, J. Yu, Experimental and molecular dynamics studies on the transport and adsorption of chloride ions in the nano-pores of calcium silicate phase: The influence of calcium to silicate ratios, Microporous and Mesoporous Materials 255 (2018) 23-35. doi: 10.1016/j.micromeso.2017.07.024.

[50] F. Georget, C. Bénier, W. Wilson, K. L. Scrivener, Chloride sorption by C-S-H quantified by SEM-EDX image analysis, Cement and Concrete Research 152 (2022) 106656. doi:10.1016/j.cemconres.2021. 106656.

[51] Z. Shi, M. R. Geiker, K. De Weerdt, T. A. Østnor, B. Lothenbach, F. Winnefeld, J. Skibsted, Role of calcium on chloride binding in hydrated Portland cement–metakaolin–limestone blends, Cement and Concrete Research 95 (2017) 205–216. doi:10.1016/j.cemconres.2017.02.003.

[52] J. J. Beaudoin, V. S. Ramachandran, R. F. Feldman, Interaction of chloride and C-S-H, Cement and Concrete Research 20 (6) (1990) 875-883. doi:10.1016/0008-8846(90)90049-4.

[53] W. Cai, Z. Xu, Z. Zhang, J. Hu, H. Huang, Y. Ma, Z. Zhang, H. Wang, S. Yin, J. Wei, C. Shi, Q. Yu, Chloride binding behavior of synthesized reaction products in alkali-activated slag, Composites Part B: Engineering 238 (2022) 109919. doi:10.1016/j.compositesb.2022.109919.

[54] M. F. Gou, X. Guan, Q. Sun, Adsorption of chloride ion by calcium silicate hydrate, Journal of Building Materials 18 (3) (2015) 363-368. doi:10.3969/j.issn.1007-9629.2015.03.002.

[55] R. J. Gowers, P. Carbone, A multiscale approach to model hydrogen bonding: The case of polyamide, The Journal of Chemical Physics 142 (22) (2015) 224907. doi:10.1063/1.4922445.

[56] I. C. Bourg, G. Sposito, Molecular dynamics simulations of the electrical double layer on smectite surfaces contacting concentrated mixed electrolyte (NaCl-CaCl2) solutions, Journal of Colloid and Interface Science 360 (2) (2011) 701-715. doi:10.1016/j.jcis.2011.04.063.
[57] I. C. Bourg, C. I. Steefel, Molecular dynamics simulations of water structure and diffusion in silica nanopores, The Journal of Physical Chemistry C 116 (21) (2012) 11556-11564. doi:10.1021/jp301299a.

[58] T. Honorio, F. Benboudjema, T. Bore, M. Ferhat, E. Vourch, The pore solution of cement-based materials: Structure and dynamics of water and ions from molecular simulations, Physical Chemistry Chemical Physics (2019). doi:10.1039/c9cp01577a.

[59] T. Honorio, H. Carasek, O. Cascudo, Water self-diffusion in C-S-H: Effect of confinement and temperature studied by molecular dynamics, Cement and Concrete Research 155 (2022) 106775. doi:10.1016/j.cemconres.2022.106775.

[60] S. M. Mutisya, J. M. De Almeida, C. R. Miranda, Probing the dynamics of water over multiple pore scales in cement by atomistic simulations, Applied Surface Science 565 (2021) 150426. doi:10.1016/j.apsusc.2021.150426.

[61] R. Metzler, J.-H. Jeon, A. G. Cherstvy, E. Barkai, Anomalous diffusion models and their properties: Non-stationarity, non-ergodicity, and ageing at the centenary of single particle tracking, Physical Chemistry Chemical Physics 16 (44) (2014) 24128-24164. doi:10.1039/C4CP03465A.

[62] T. Honorio, M. Trifa, T. Herin, S. Le Caër, H₂ Anisotropic Subdiffusion and Induced Expansion in Portlandite, Gibbsite, and Boehmite, The Journal of Physical Chemistry C 128 (45) (2024) 19085-19097. doi:10.1021/acs.jpcc.4c04681.

[63] L. Yuan-Hui, S. Gregory, Diffusion of ions in sea water and in deep-sea sediments, Geochimica et Cosmochimica Acta 38 (5) (1974) 703-714. doi:10.1016/0016-7037(74)90145-8.

[64] P. Vanýsek, Ionic Conductivity and Diffusion, in: CRC Handbook of Chemistry and Physics, 86th Edition, Taylor and Francis, Boca Raton, FL, 2006, pp. 5-76-5-78.

[65] A. V. Verde, M. Santer, R. Lipowsky, Solvent-shared pairs of densely charged ions induce intense but short-range supra-additive slowdown of water rotation, Physical Chemistry Chemical Physics 18 (3) (2016) 1918-1930. doi:10.1039/C5CP05726D.

[66] G. Stirnemann, E. Wernersson, P. Jungwirth, D. Laage, Mechanisms of Acceleration and Retardation of Water Dynamics by Ions, Journal of the American Chemical Society 135 (32) (2013) 11824-11831. doi:10.1021/ja405201s.

[67] A. P. Lyubartsev, A. Laaksonen, Concentration Effects in Aqueous NaCl Solutions. A Molecular Dynamics Simulation, The Journal of Physical Chemistry 100 (40) (1996) 16410-16418. doi:10.1021/jp961317h.

[68] S. Bouazizi, S. Nasr, Self-diffusion coefficients and orientational correlation times in aqueous NaCl solutions, Journal of Molecular Liquids 162 (2) (2011) 78-83. doi:10.1016/j.molliq.2011.06.004.

[69] A. Ghaffari, A. Rahbar-Kelishami, MD simulation and evaluation of the self-diffusion coefficients in aqueous NaCl solutions at different temperatures and concentrations, Journal of Molecular Liquids 187 (2013) 238-245. doi:10.1016/j.molliq.2013.08.004.

[70] V. Vitagliano, P. A. Lyons, Diffusion Coefficients for Aqueous Solutions of Sodium Chloride and Barium Chloride, Journal of the American Chemical Society 78 (8) (1956) 1549-1552. doi:10.1021/ja01589a011.

[71] D. G. Mapa, H. Zhu, F. Nosouhian, N. Shanahan, K. A. Riding, A. Zayed, Chloride binding and diffusion of slag blended concrete mixtures, Construction and Building Materials 388 (2023) 131584. doi:10.1016/j.conbuildmat.2023.131584.
[72] M. D. A. Thomas, R. D. Hooton, A. Scott, H. Zibara, The effect of supplementary cementitious materials on chloride binding in hardened cement paste, Cement and Concrete Research 42 (1) (2012) 1-7. doi:10.1016/j.cemconres.2011.01.001.

[73] P. Pivonka, C. Hellmich, D. Smith, Microscopic effects on chloride diffusivity of cement pastes—a scale-transition analysis, Cement and Concrete Research 34 (12) (2004) 2251-2260. doi:10.1016/j.cemconres.2004.04.010.

[74] J.-P. Korb, L. Monteilhet, P. McDonald, J. Mitchell, Microstructure and texture of hydrated cement-based materials: A proton field cycling relaxometry approach, Cement and Concrete Research 37 (3) (2007) 295-302. doi:10.1016/j.cemconres.2006.08.002.

[75] H. Li, E. Fratini, W.-S. Chiang, P. Baglioni, E. Mamontov, S.-H. Chen, Dynamic behavior of hydration water in calcium-silicate-hydrate gel: A quasielastic neutron scattering spectroscopy investigation, Physical Review E 86 (6) (2012) 061505. doi:10.1103/PhysRevE.86.061505.

[76] E. Fratini, A. Faraone, F. Ridi, S.-H. Chen, P. Baglioni, Hydration Water Dynamics in Tricalcium Silicate Pastes by Time-Resolved Incoherent Elastic Neutron Scattering, The Journal of Physical Chemistry C 117 (14) (2013) 7358-7364. doi:10.1021/jp312684p.
# Supplementary Information

Diffusion and physicochemical behavior of chloride and sodium ions in C-S-H gel pores : A molecular dynamics investigation

Tarek Ihaddadeneª, Jérôme Claverieª, François Bignonnetª, Ouali Amiriª

ªNantes Université, École Centrale Nantes, CNRS, GeM, UMR 6183, F-44600 Saint-Nazaire, France

## S1. Parameters of the potential energy function

### S1.1. Lennard-Jones parameters

<table><thead><tr><th>Pairs</th><th>ε (kcal mol<sup>-1</sup>)</th><th>σ (Å)</th></tr></thead><tbody><tr><td>Ca-Ca</td><td>0.000005</td><td>5.567000</td></tr><tr><td>Cw-Cw</td><td>0.100000</td><td>2.872000</td></tr><tr><td>Hsi-Hsi</td><td>0.000000</td><td>0.000000</td></tr><tr><td>Hh-Hh</td><td>0.000000</td><td>0.000000</td></tr><tr><td>Hw-Hw</td><td>0.000000</td><td>0.000000</td></tr><tr><td>Si-Si</td><td>0.000002</td><td>3.302000</td></tr><tr><td>Na-Na</td><td>0.544300</td><td>2.350000</td></tr><tr><td>Osih-Osih</td><td>0.061800</td><td>3.448000</td></tr><tr><td>Ob-Ob</td><td>0.061800</td><td>3.448000</td></tr><tr><td>Osi-Osi</td><td>1.243000</td><td>2.735000</td></tr><tr><td>Oh-Oh</td><td>0.155400</td><td>3.166000</td></tr><tr><td>Ow-Ow</td><td>0.155400</td><td>3.166000</td></tr><tr><td>Cl-Cl</td><td>0.100100</td><td>4.400000</td></tr></tbody></table>

In the classical force field used, the geometric mixing rule is employed for the Lennard-Jones parameters based on the sensitivity study at Section S3. These parameters are calculated according to the following formulas:

$$ \sigma_{ij} = \sqrt{\sigma_{ii} + \sigma_{jj}} $$

$$ \epsilon_{ij} = \sqrt{\epsilon_{ii} + \epsilon_{jj}} $$

### S1.2. Harmonic stretching and bending parameters

<table><thead><tr><th>Pairs</th><th>k (kcal mol<sup>-1</sup> Å<sup>-2</sup>)</th><th>σ (Å)</th></tr></thead><tbody><tr><td>Hsi-Osih</td><td>554.135</td><td>1.000</td></tr><tr><td>Hh-Oh</td><td>554.135</td><td>1.000</td></tr><tr><td>Hw-Ow</td><td>554.135</td><td>1.000</td></tr></tbody></table>

<table><thead><tr><th>Triples</th><th>k (kcal mol<sup>-1</sup> °<sup>-2</sup>)</th><th>θ<sub>0</sub> (°)</th></tr></thead><tbody><tr><td>Hw-Ow-Hw</td><td>45.770</td><td>109.470</td></tr></tbody></table>
### S1.3. Charge parameters

Initially, the charges are assigned to each atom types as follow.

<table><thead><tr><th>Atom type</th><th>Charge</th></tr></thead><tbody><tr><td>Ca</td><td>1.430</td></tr><tr><td>Cl</td><td>-1.000</td></tr><tr><td>Cw</td><td>2.000</td></tr><tr><td>Hh</td><td>0.425</td></tr><tr><td>Hsi</td><td>0.290</td></tr><tr><td>Hw</td><td>0.410</td></tr><tr><td>Na</td><td>1.000</td></tr><tr><td>Ob</td><td>-1.140</td></tr><tr><td>Oh</td><td>-0.950</td></tr><tr><td>Osi</td><td>-1.000</td></tr><tr><td>Osih</td><td>-1.000</td></tr><tr><td>Ow</td><td>-0.820</td></tr><tr><td>Si</td><td>1.720</td></tr></tbody></table>

To account for different Ca/Si ratios, the number of atoms of each type vary from one system to another. After assignment of the initial charges, the excess charge was equally distributed among Ob, Osi and Osih atoms.

<table><caption>Table 1: Summary table of charge distribution on Ob, Osi, and Osih sites after the final relaxation using the CSH-FF force field, for different Ca/Si ratios (1.2, 1.4, 1.6, and 1.8).</caption><thead><tr><th>Atoms</th><th>Ob</th><th>Osi</th><th>Osih</th></tr></thead><tbody><tr><th>Ca/Si</th><td></td><td></td><td></td></tr><tr><th>1.2</th><td>-1.209</td><td>-1.069</td><td>-1.069</td></tr><tr><th>1.4</th><td>-1.248</td><td>-1.108</td><td>-1.108</td></tr><tr><th>1.6</th><td>-1.309</td><td>-1.169</td><td>-1.169</td></tr><tr><th>1.8</th><td>-1.393</td><td>-1.253</td><td>-1.253</td></tr></tbody></table>
## S2. Evolution of the box parameters during the creation of C-S-H models

The following tables provides a detailed comparison of the box parameters, starting with the initial values, followed by the changes observed after applying the ReaxFF force field. The final results are then compared with the parameters obtained after applying the specific CSH-FF force field, all within the NPT ensemble. This analysis helps to understand the impact of these force fields on the structure and dynamics of the simulated boxes.

### S2.1. Initial box parameters

<table><thead><tr><th>Box Parameters</th><th>Ca/Si</th><th>1.2</th><th>1.4</th><th>1.6</th><th>1.8</th></tr></thead><tbody><tr><td>a (Å)</td><td>20.07</td><td>20.07</td><td>20.07</td><td>20.07</td><td>20.07</td></tr><tr><td>b (Å)</td><td>30.82</td><td>30.82</td><td>30.82</td><td>30.82</td><td>30.82</td></tr><tr><td>c (Å)</td><td>22.77</td><td>22.77</td><td>22.77</td><td>22.77</td><td>22.77</td></tr><tr><td>α (°)</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td></tr><tr><td>β (°)</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td></tr><tr><td>γ (°)</td><td>90.59</td><td>90.59</td><td>90.59</td><td>90.59</td><td>90.59</td></tr></tbody></table>

### S2.2. Box parameters after the application ReaxFF

<table><thead><tr><th>Box Parameters</th><th>Ca/Si</th><th>1.2</th><th>1.4</th><th>1.6</th><th>1.8</th></tr></thead><tbody><tr><td>a (Å)</td><td>20.11</td><td>20.09</td><td>21.14</td><td>21.81</td><td>21.81</td></tr><tr><td>b (Å)</td><td>30.87</td><td>30.85</td><td>32.47</td><td>33.49</td><td>33.49</td></tr><tr><td>c (Å)</td><td>22.81</td><td>22.79</td><td>23.99</td><td>24.44</td><td>24.44</td></tr><tr><td>α (°)</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td></tr><tr><td>β (°)</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td><td>90.00</td></tr><tr><td>γ (°)</td><td>90.59</td><td>90.59</td><td>90.59</td><td>90.59</td><td>90.59</td></tr></tbody></table>

### S2.3. Box parameters after duplicating the cell along the a-direction and application of the CSH-FF force field

<table><thead><tr><th>Box Parameters</th><th>Ca/Si</th><th>1.2</th><th>1.4</th><th>1.6</th><th>1.8</th></tr></thead><tbody><tr><td>a (Å)</td><td>40.51</td><td>39.83</td><td>39.70</td><td>39.58</td><td>39.58</td></tr><tr><td>b (Å)</td><td>31.12</td><td>31.02</td><td>30.82</td><td>30.61</td><td>30.61</td></tr><tr><td>c (Å)</td><td>22.74</td><td>22.53</td><td>24.44</td><td>26.45</td><td>26.45</td></tr><tr><td>α (°)</td><td>90.45</td><td>89.65</td><td>89.44</td><td>88.44</td><td>88.44</td></tr><tr><td>β (°)</td><td>90.21</td><td>89.71</td><td>88.21</td><td>90.69</td><td>90.69</td></tr><tr><td>γ (°)</td><td>91.65</td><td>91.92</td><td>91.77</td><td>93.11</td><td>93.11</td></tr></tbody></table>
S3. Sensitivity study for the choice of the water model and the Lennard-Jones mixing rule

In this study, the calculation of diffusion coefficients is carried out with a simulation time of 10 ns by combining several parameters. On the left of Fig. 2, the analysis focuses on a C-S-H model with a Ca/Si ratio of 1.2, a pore size of 30 Å, and a concentration of 0.5 M at 300 K, using SPC and SPC/E force fields for water, along with two mixing rules (geometric and arithmetic) for the Lennard-Jones parameters.

In addition, a comparison of diffusion profiles is presented for two different temperatures with an SPC water model and a geometric mixing rule. Diffusion coefficient profiles are plotted for a C-S-H model with a Ca/Si ratio of 1.2, a 30 Å nanopore, and a molar concentration of 0.5 M at temperatures of 300K and 291K (see Fig. 1, right). As reported in the literature, a raise in temperature increases the diffusion coefficient. The results are in good agreement with experimental measurements of the diffusion coefficient of Cl⁻ at infinite dilution: 2.03 × 10⁻⁹ m²/s at 298 K and 1.71 × 10⁻⁹ m²/s at 291 K [1, 2].

<figure><img src="image_14.png" /><figcaption>Figure 1: Sensitivity study and comparison of diffusion coefficients based on SPC and SPC/E force fields and temperature, for Ca/Si=1.2, a pore size of 30 Å and a concentration 0.5M.</figcaption></figure>

Following this study, we choose an SPC force field for water molecules and a geometric mixing rule for the Lennard-Jones parameters for all atoms in our simulations.
S4. Sensitivity study on the influence of simulation time and observation time on the diffusion coefficient profile

<figure><img src="image_15.png" /><figcaption>Figure 2: Sensitivity analysis for a C-S-H model with a Ca/Si ratio of 1.2, a pore size of 30 Å, and a NaCl concentration of 0.5 M. (a) Effect of simulation time on the chloride diffusion coefficient. (b) Effect of the observation time $\tau$ on the chloride diffusion coefficient.</figcaption></figure>
S5. Radial distribution functions

<figure><img src="image_16.png" /><figcaption>Figure 3: Radial Distribution Functions (RDF) for Cw-Osi, Cw-Cl, Na-Osi and Na-Cl pairs at different Ca/Si ratios for a pore size of 30 Å and a concentration of 0.5 M.</figcaption></figure>

The function used to calculate the radial distribution function (RDF) is a discrete function, tailored to the discrete nature of our simulated systems. In a discrete system, particles are represented by fixed positions or coordinates calculated at specific intervals in space. Consequently, the RDF is computed by counting the number of particle pairs located at a given distance $r$ within defined intervals ($\Delta r$). This approach allows for the estimation of the relative density of particles around a given particle, while taking into account the inherent discretization of the simulated model.

$$g(r) = \frac{1}{4\pi r^2 \Delta r \rho_0} \sum_{i,j} \delta(r - |r_i - r_j|) \quad (1)$$

* $\Delta r$ is the thickness of a spherical shell around the reference particle.

* $|r_i - r_j|$ is the distance between particle i and particle j.

* The summation is performed over all pairs of particles i, j.

In the radial distribution function (RDF) for each pair of atoms, the first observed peak generally corresponds to the internal coordination sphere, indicating the most probable distance between a reference atom and its nearest neighbors. This first coordination sphere is crucial as it defines the immediate environment around each atom. Subsequent peaks in the RDF represent the external coordination spheres, which can manifest as one or more additional layers around the reference atom.
<figure><img src="image_17.png" /><figcaption>Figure 4: Atomic density profile for a gel pore of Ca/Si=1.8 for different pore sizes (5 Å, 10 Å, 20 Å, 30 Å, or 50 Å) and concentrations (0.5 M and 2 M). The left axis in each figure represents the overall atomic density of the solution and the solid (envelope of atomic density profiles), whereas the right axis provides a specific zoom on Cl, Na, and Cw atoms. The dashed blue lines represent the inner interfaces while the dashed black lines represent the outer interfaces.</figcaption></figure>
## S7. Bound and free ions calculations

The methodology used to determine the proportions of bound and free ions involves the following steps: First, it is assumed that the prominent peaks observed on both surfaces represent the bound ions (highlighted in red in the Fig. 5). The remaining ions are considered free (highlighted in blue). The trapezoidal rule is applied to calculate the integral of the atomic density profile and quantify the proportion of adsorbed and free species (as illustrated on the right in Fig. 5).

<figure><img src="image_18.png" /><figcaption>Figure 5: The atomic density profile of chloride ions is examined within a 30 Å pore of a C-S-H model with a Ca/Si ratio of 1.2 and a concentration of 0.5 M. On the right, the cumulative integral is provided to assess the percentage composition of bound and free chlorides.</figcaption></figure>

## S8. Influence of the pore size on the proportion of adsorbed ions

<table><caption>Table 2: The proportion of adsorbed Na⁺ and Cl⁻ ions for different pore sizes for a Ca/Si ratio of 1.8.</caption><thead><tr><th rowspan="2">Pore size (Å)</th><th colspan="2">0.5 M</th><th colspan="2">2 M</th></tr><tr><th>Na</th><th>Cl</th><th>Na</th><th>Cl</th></tr></thead><tbody><tr><td>5</td><td>–</td><td>–</td><td>100 %</td><td>94.1%</td></tr><tr><td>10</td><td>99.4%</td><td>86.8%</td><td>97.0%</td><td>88.0%</td></tr><tr><td>20</td><td>96.8%</td><td>67.1%</td><td>80.4%</td><td>75.1%</td></tr><tr><td>30</td><td>88.3%</td><td>65.1%</td><td>66.0%</td><td>61.2%</td></tr><tr><td>50</td><td>52.8%</td><td>52.5%</td><td>41.8%</td><td>48.3%</td></tr></tbody></table>
## S9. Error calculation for diffusion profiles

In this study, error estimates on the diffusion coefficient profiles are based on the weighted linear regression, assuming that the MSD vary linearly with time (Fickian diffusion regime). The weight $w_i$ is calculated as follows:

$$w_i = \frac{1}{\sigma_{\text{MSD},i}^2} \qquad (2)$$

where $\sigma_{\text{MSD},i}$ represents the uncertainty (or standard deviation) associated with each MSD data point $i$. The $\sigma_{\text{MSD},i}$ are calculated over all single-particle MSD recorded in each spatial block. Weighted linear regression aims to fit a line of the form:

$$y = mx + b \qquad (3)$$

where $m$ is the slope and $b$ is the y-intercept. The equations of weighted linear regression are based on minimizing the sum of the weighted squared errors. This gives:

$$m = \frac{n \sum w_i x_i y_i - \sum w_i x_i \sum w_i y_i}{n \sum w_i x_i^2 - (\sum w_i x_i)^2} \qquad (4)$$

$$b = \frac{\sum w_i y_i \sum w_i x_i^2 - \sum w_i x_i \sum w_i x_i y_i}{n \sum w_i x_i^2 - (\sum w_i x_i)^2} \qquad (5)$$

where $n$ is the total number of data points, and the sums are performed over all the data $(x_i, y_i)$.

Finally, the standard error is calculated as the root square of the variance of the slope $m$.

## S10. Effects of EDL overlap and confinement on the diffusion coefficient

<figure><img src="image_19.png" /><figcaption>Figure 6: Diffusion coefficient profiles of Cl⁻ and Na⁺ ions in pores of 5 Å and 10 Å for Ca/Si=1.8 and a concentration of 2 M.</figcaption></figure>

## References

[1] L. Yuan-Hui, S. Gregory, Diffusion of ions in sea water and in deep-sea sediments, Geochimica et Cosmochimica Acta 38 (5) (1974) 703–714. doi:10.1016/0016-7037(74)90145-8.

[2] P. Vanýsek, Ionic Conductivity and Diffusion, in: CRC Handbook of Chemistry and Physics, 86th Edition, Taylor and Francis, Boca Raton, FL, 2006, pp. 5–76–5–78.