Introduction
In many nuclear waste repository concepts, metallic canisters (e.g. steel) are buried in deep tunnels with bentonite used as backfill material. Upon repository closure, elevated temperatures and variable water-saturation conditions promote canister corrosion, releasing iron species into the bentonite. Numerous studies (Kaufhold et al., Reference Kaufhold, Hassel, Sanders and Dohrmann2015; Latta et al., Reference Latta, Neumann, Premaratne and Scherer2017), including the FEBEX in situ heater test at the Grimsel Test Site in Switzerland (Hadi et al., Reference Hadi, Wersin, Serneels and Grenèche2019) have documented the resulting alteration of bentonite in close proximity to steel material. Diffusion of iron species into the bentonite has resulted in the formation of large (>140 mm wide) red, orange, and blue regions due to the formation of iron minerals such as goethite. Excess Fe2+ ions may also diffuse far away from the container interface and adsorb onto the clay surface.
In contrast to the formation of goethite in the work by Kaufhold et al. (Reference Kaufhold, Hassel, Sanders and Dohrmann2015), Milodowski et al. (Reference Milodowski, Cave, Kemp, Green, Williams, Shaw, Gowing and Eatherington2007) concluded that the bentonite immediately adjacent to the steel canister was significantly enriched in iron, not as iron oxide or iron oxyhydroxide but through the substitution of iron into the octahedral layers of bentonite, forming minerals such as nontronite. The cation exchange capacity (CEC) of altered bentonite is slightly reduced, with depletions in exchangeable Ca²⁺ and Na⁺ and an increase in exchangeable Fe2+. In general, bentonite interaction with iron has been shown to reduce shrinkage significantly, possibly due to its conversion to nontronite and chlorite (Jones et al., Reference Jones, Murphy, Waite and Collins2017).
Laboratory studies further showed that sorbed Fe2+ on clays can undergo redox transformations. Under anoxic conditions, Fe2+ sorbed onto low-iron clay surfaces can be partly oxidized to Fe3+, forming Fe(III)-bearing minerals (Géhin et al., Reference Géhin, Grenèche, Tournassat, Brendlé, Rancourt and Charlet2007; Van Groeningen et al., Reference Van Groeningen, Thomas Arrigo, Byrne, Kappler, Christi and Kretzschmar2020). When exposed to air, this oxidation proceeds rapidly and completely, demonstrating that clay minerals were unable to protect ferrous iron from oxidation. Nonetheless, Fe2+ can persist in clay interlayers under reducing conditions as the conversion to Fe(III) in anoxic conditions is only partial. Furthermore, Manjanna (Reference Manjanna2008) prepared Fe2+-Mnt successfully by reducing Fe3+-Mnt and suggested that, for long-time preservation, Fe2+-montmorillonite must be stored in deoxygenated and/or low relative humidity (<40%) conditions at room temperature or below. A further report suggested that iron adsorption at 2<pH<7.5 occurs both via cation exchange on interlayer surfaces and as sorption onto highly reactive surface sites (e.g. edge and defect) (Géhin et al., Reference Géhin, Grenèche, Tournassat, Brendlé, Rancourt and Charlet2007). The conversion from Fe2+ to Fe3+ increases with pH. This reaction is reversible (e.g. reducing the pH can reduce Fe3+ to Fe2+ in solution).
The adsorption of Fe2+ onto clay depends significantly on structural Fe. On a low-Fe content clay or dithionite-reduced clays, Fe2+ uptake generally follows the sorption behavior of other divalent metals (Soltermann et al., Reference Soltermann, Marques Fernandes, Baeyens, Dähn, Joshi, Scheinost and Gorski2014a), primarily through cation exchange and surface complexation with clay edges or surface defects, which are highly pH dependent. When Fe2+ adsorbs onto Fe(III)-rich clays, the amount of Fe adsorbed is unexpectedly large (Schaefer et al., Reference Schaefer, Gorski and Scherer2011; Soltermann et al., Reference Soltermann, Marques Fernandes, Baeyens, Dähn, Joshi, Scheinost and Gorski2014a), possibly due to electron transfer between sorbed Fe2+ and structural Fe(III). In low-Fe clays, sorbed Fe2+ can be partially oxidized to mineral-bound Fe(III) complexes, whereas in high Fe-clays, sorbed Fe2+ is fully oxidized.
Despite this progress, the reaction pathways remain unclear. Aqueous Fe2+ may: (1) adsorb to basal surfaces by cation exchange; (2) form adsorption complexes at edge sites; (3) precipitate as secondary minerals; or (4) undergo interfacial redox reactions with structural Fe(III) (O’Reilly et al., Reference O’Reilly, Watkins and Furukawa2005; Merola et al., Reference Merola, Fournier and McGuire2007; Jaisi et al., Reference Jaisi, Liu, Dong, Blake and Fein2008; Schaefer et al., Reference Schaefer, Gorski and Scherer2011; Soltermann et al., Reference Soltermann, Marques Fernandes, Baeyens, Dähn, Joshi, Scheinost and Gorski2014a; Soltermann et al., Reference Soltermann, Marques Fernandes, Baeyens, Miehé-Brendlé and Dähn2014b; Van Groeningen et al., Reference Van Groeningen, Thomas Arrigo, Byrne, Kappler, Christi and Kretzschmar2020). The latter pathway suggests Fe2+ intercalation into clay interlayers with electron transfer to structural Fe(III), a process influenced by Fe2+ concentration and Fe(III) accessibility. Such redox-driven intercalation could alter surface charge, CEC, and swelling/shrinkage behavior, while simultaneously promoting secondary mineral formation such as that of chlorite. This background points to a key knowledge gap: the behavior and hydration of Fe2+ in clay interlayers under repository-relevant conditions. Accordingly, while Fe-montmorillonite interactions may involve surface sorption, structural substitution, redox transformations, and secondary mineral precipitation, the present study focused specifically on Fe2+ interlayer exchange and hydration, and does not address long-term structural incorporation, Fe-bearing clay phase stabilization, or secondary mineral formation. The hypothesis is that Fe2+ generated from canister corrosion can exchange with interlayer Na+ in montmorillonite, modifying its hydration and structural properties. Evidence has been found for Fe2+/Na+ exchange in smectite interlayers under reducing conditions (Kamei et al., Reference Kamei, Oda, Mitsui, Shibata and Shinozaki1999), reinforcing this possibility. In the present study, molecular dynamics simulation was used to investigate the hydration structural properties and energetics of iron ion intercalation into the clay interlayer space.
Materials and methods
To address the knowledge gap about the hydration and structure of iron-exchanged Mnt, molecular dynamics simulations were used to investigate Fe-Mnt hydration properties. The simulation system including two clay layers with the x–y surface area of 4.1 nm×3.5 nm is presented in Fig. 1a. The clay structure was prepared by replicating a unit cell by 8×4×2 in x, y, and z-dimensions. The montmorillonite structure adopts a dioctahedral cis-vacant configuration, which is prevalent in engineered clay barrier materials (Drits and Zviagina, Reference Drits and Zviagina2009; Subramanian et al., Reference Subramanian, Whittaker, Ophus and Lammers2020). A previous study (Gao et al., Reference Gao, Zhang, Liu and Lu2025) has shown that cis- and trans-vacant montmorillonite exhibit similar interlayer hydration behaviors; nevertheless, the vacancy configuration is provided here for completeness. For all systems, isomorphic substitution occurred in octahedral sheets. Different numbers of water molecules were placed in the interlayer spaces to create hydrated structures, and the amount of water in the interlayers was expressed as the number of water molecules per unit cell of clay, i.e. w u.c.–1. Interlayer ions, Na+, Fe2+, or Fe3+, were placed in the interlayer space also to balance the structural charge. There were 24 Na+ ions in each interlayer, and when Fe2+ or Fe3+ replaced Na+, their number was adjusted according to the charge.
(a) Simulation snapshot demonstrating the periodic clay system including two clay layers and two interlayers with water and ions to study hydration properties of Na+, Fe2+, and Fe3+. (b) Clay system with edge to study the exchange reaction of Fe2+ (pink) with Na+ (yellow). Cl– (blue) balances the charge of Fe2+.
In the simulations in Fig. 1a, Na+ ions were replaced by Fe ions in the interlayer spaces without studying the exchange reaction explicitly. The simulation in Fig. 1b focused on ion exchange in clay because it is one of the fundamental processes in the bentonite buffer, e.g. radionuclide exchange for Na+ in bentonite; Ca2+, Mg2+ exchange for Na+ ions in the interlayers leading to a change in the concentration of the ions (Kaufhold et al., Reference Kaufhold, Dohrmann, Ufer and Kober2018). The system in Fig. 1b is different from that in Fig. 1a as the clay particles have edges along the x-direction which expose the interlayers to aqueous solutions while the system in Fig. 1a is periodic (e.g. infinitely long) in the x-direction. The system in Fig. 1b allows the exchange of water and ions in the interlayer spaces with those in bulk aqueous solution. The method for constructing clay particles with edges can be found in previous studies by the present authors (Ho et al., Reference Ho, Criscenti and Greathouse2019; Ho et al., Reference Ho, Wang, Jové Colón and Coker2020; Ho et al., Reference Ho, Coker, Jové-Colón and Wang2022; Ho et al., Reference Ho, Yang, Rempe, Dasgupta, Johnston, Xu, Zwier and Mills2023b). Exchange reactions via equilibrium simulations and free-energy calculation were studied here. In the equilibrium simulations, iron chloride was added to the bulk regions and the amount of iron ions moving into the interlayer was monitored. In the free-energy calculation, a metadynamics technique (Barducci et al., Reference Barducci, Bussi and Parrinello2008) was applied using the Colvars module (Fiorin et al., Reference Fiorin, Klein and Hénin2013) to calculate the free energy of Fe2+ or Fe3+ ions moving from bulk water into the interlayer space. Such free-energy landscapes provide a comprehensive evaluation about energy barrier and thermodynamics of iron intercalations. This technique has been applied to calculate gas intercalation into the clay interlayer space (Ho et al., Reference Ho, Yang, Rempe, Dasgupta, Johnston, Xu, Zwier and Mills2023b; Choudhary and Ho, Reference Choudhary and Ho2024); this includes hydrogen generated from canister corrosion being intercalated into bentonite which is relevant to the nuclear repository (Ho et al., Reference Ho, Jove-Colon and Wang2023a). A summary of simulated systems is reported in Table 1.
Systems simulated and the corresponding number of ions
The ClayFF force field (Cygan et al., Reference Cygan, Liang and Kalinichev2004) was used to simulate clay with additional Si-O-H and Al-O-H angle terms describing edge hydroxyl groups (Pouvreau et al., Reference Pouvreau, Greathouse, Cygan and Kalinichev2017; Pouvreau et al., Reference Pouvreau, Greathouse, Cygan and Kalinichev2019). Water molecules were simulated using a flexible SPC water model (Teleman et al., Reference Teleman, Jonsson and Engstrom1987), and the Na+ ion was modeled using the force field proposed by Smith and Dang (Reference Smith and Dang1994). Iron-water force-field parameters were taken from Babu and Lim (Reference Babu and Lim2006). Fe2+ and Fe3+ have the same Lennard-Jones parameters and have charges of +2 and +3, respectively. Fe2+ and Fe3+ were modeled using non-reactive force fields that treat Fe as fully solvated cations interacting electrostatically with water molecules and clay surfaces. These force fields were selected because they have been applied widely in classical molecular dynamics studies of hydrated metal cations and ion-clay systems, and they provide a computationally efficient framework for comparing interlayer hydration and swelling behavior systematically across various exchangeable cations. The present model does not account for hydrolysis reactions or for the formation of Fe-hydroxyl coordination complexes in aqueous solution. Instead, the simulations were designed to represent Fe cations under conditions where aqueous speciation and redox reactions are suppressed, allowing the effects of ionic charge and hydration strength on interlayer structure to be isolated. The Lennard-Jones (LJ) potential was expressed as:
$ {V}_{\mathrm{LJ}}=4\unicode{x025B} \left[{\left(\frac{\unicode{x03C3}}{r}\right)}^{12}-{\left(\frac{\unicode{x03C3}}{r}\right)}^6\right], $
where r is the distance between two atoms, and ε and σ are the depth of the potential energy well and the distance at which the LJ potential is zero, respectively. LJ interactions among atoms were calculated using the Lorentz-Berthelot mixing rules:
$ {\unicode{x025B}}_{ij}=\sqrt{\unicode{x025B}_{ii}{\unicode{x025B}}_{jj}} $
and
$ {\unicode{x03C3}}_{ij}=\left({\unicode{x03C3}}_{ii}+{\unicode{x03C3}}_{jj}\right)/2 $
. Short-range interactions were calculated using a cut-off distance of 10 Å. Long-range electrostatic interactions were calculated using the PPPM (particle-particle-particle-mesh) solver (Hockney and Eastwood, Reference Hockney and Eastwood1988). All simulations were performed using the LAMMPS package (Plimpton, Reference Plimpton1995) at a temperature of 300K, pressure of 1 atm using the Nose-Hoover thermostat (Nose, Reference Nose1984) and barostat, and a timestep of 1 fs.
Results and Discussion
The d-spacing of Mnt systems containing 100% Na+, Fe2+, or Fe3+ ions in the interlayers is shown in Fig. 2a. These results were obtained from the simulation system shown in Fig. 1a. The results indicate that Fe ions influence the d-spacing primarily at low water contents (e.g. <7 water molecules per unit cell, or w u.c.–1). In other words, when Na+ ions were completely replaced by either Fe2+ or Fe3+, changes in the d-spacing of Mnt were observed only when the interlayer contains less than a bilayer of water (d-spacing <14 Å). At higher water contents, beyond the bilayer region (d-spacing >14 Å), the type of interlayer cation (Na+, Fe2+, or Fe3+) has little effect on the swelling behavior of the clay. When the water content is less than a bilayer, the d-spacing of Fe-Mnt is ~1–2 Å larger than that of Na-Mnt at equivalent hydration levels. This difference arises from the greater hydration energy of Fe ions compared with Na+, which leads to a larger hydration shell surrounding the Fe ions. Simulation snapshots in Fig. 2c, d show that at a water content of 3.7 w u.c.–1, water molecules are more evenly distributed in the Fe2+ system than in the Fe3+ system. Because there are fewer interlayer ions in the Fe3+ case and Fe3+ attracts water more strongly due to its higher charge, water molecules are more tightly coordinated around the ions. This leads to a smaller d-spacing in Fe3+ compared with Fe2+.
(a) d-spacing as a function of the number of water molecules per clay unit cells (or water mass per clay mass without considering mass of ions) for Fe2+, Fe3+, and Na+ montmorillonite. (b) Comparison of d-spacing in montmorillonite when 50% of the Na+ ions in the interlayer space are replaced by Fe2+ ions, along with those of pure Na+- and Fe2+-montmorillonite. (c, d) Simulation snapshots for Fe2+ and Fe3+ systems at 3.71 w u.c.–1, respectively.
The expansion of the clay varies with increasing amount of Fe2+ in the interlayer space compared with Na+ ions. For example, the d-spacing of Fe2+-Mnt and Na-Mnt is compared with that of Mnt containing both Fe2+ and Na+ in the interlayer (i.e. Fe–Na-Mnt) in Fig. 2b. At low water contents, the d-spacing increases in the order Na-Mnt<Fe–Na-Mnt<Fe-Mnt. Therefore, depending on the extent of Fe ion exchange for Na+ in the interlayer space, the clay can exhibit different swelling behaviors. The present d-spacing results for Fe-Mnt and Na-Mnt are consistent with the experimental findings of Manjanna (Reference Manjanna2008).
The hydration energies for Na+-, Fe2+-, and Fe3+-montmorillonite obtained from the simulation system (Fig. 1a) are shown in Fig. 3. The hydration energy was defined as the change in potential energy, U, associated with adding n water molecules to a dry (n=0) clay system:
$$ \varDelta U(n)=\frac{U(n)-U\left(n=0\right)}{n}. $$
Comparison of hydration energies for Na+-, Fe2+-, and Fe3+-montmorillonite.
The results for Na+ agree well with those of Teich-McGoldrick et al. (Reference Teich-McGoldrick, Greathouse, Jové-Colón and Cygan2015). As the number of water molecules increased, the hydration energy approached the value characteristic of bulk water (~10 kcal mol–1). The hydration energy exhibited a minimum at low water content (<1 monolayer of water). The profile for Fe2+ is similar to that of other divalent cations (Teich-McGoldrick et al., Reference Teich-McGoldrick, Greathouse, Jové-Colón and Cygan2015), showing a sharper minimum compared with the broader minimum observed for monovalent ions. The hydration energy of Fe3+ is generally lower (i.e. more negative) than that of Na+ but greater than that of Fe2+. Note that the hydration energy differs from hydration free energy, which includes entropic contributions. Such entropic effects have been examined in previous studies (Whitley and Smith, Reference Whitley and Smith2004; Zhou et al., Reference Zhou, Chen, Zhu, Lu, Zhu and He2020; Li et al., Reference Li, Lu, Chen, Zhang, Cheng, Liu and Yin2023).
The diffusion coefficients of water and ions are reported in Fig. 4 for the interlayer space of the system shown in Fig. 1a, to examine how the exchange of Fe ions for Na+ ions affects the transport properties of the ions and water. The results indicate that both water (Fig. 4a) and ions (Fig. 4b) diffuse more slowly in Fe-Mnt than in Na-Mnt. While the difference in water diffusion between Na-Mnt and Fe-Mnt is not substantial, Fe2+ ions diffuse much more slowly, by at least an order of magnitude, compared with Na+ ions. The results in Fig. 4 also suggest that the diffusion of both water and ions depends heavily on the amount of water present in the interlayers.
(a) Comparison of water diffusion coefficients in Na+- and Fe2+-Mnt. (b) Diffusion coefficients of Na+ and Fe2+ ions in the interlayers.
As mentioned above, in the simulations shown in Fig. 1a, the Na+ ions were replaced with Fe ions in the interlayer spaces without explicitly modeling the exchange reaction. To investigate the intercalation of iron ions from the bulk solution into the interlayer space initially occupied by Na+ ions, the system shown in Fig. 1b was used. For example, the positions of five Fe3+ ions that were initially placed in the bulk solution are shown in Fig. 5a. The presence of Fe³⁺ ions in the interlayer indicates that, within a few nanoseconds of the simulation, Fe3+ ions had already migrated into the interlayer. The free energy of ion exchange between the bulk and the interlayer is reported in Fig. 5b. The data indicate that both Fe2+ and Fe3+ are thermodynamically more stable within the interlayer space than in the bulk solution, with Fe3+ being more favorable. In the edge-exposed simulation set-up used to calculate the potential of mean force (PMF), both ions and water molecules were free to exchange dynamically between the interlayer and the bulk solution. As Fe2+ or Fe3+ migrates from the bulk across the clay edge into the interlayer space, this process is inherently coupled to water redistribution, including partial dehydration in the bulk clay, reorganization of hydration shells, and rehydration under interlayer confinement. These solvent rearrangements were captured explicitly in the PMF and contribute to both energetic and entropic components of the free-energy landscape. Consequently, the calculated PMF profiles reflect not only ion–clay electrostatic interactions, but also the cost of water displacement and confinement effects during ion intercalation, explaining the distinct barriers and stability observed for Fe2+ and Fe3+.
(a) Positions of five Fe3+ ions throughout the simulation of the system shown in Fig. 1b. (b) Potential of mean force for the exchange of an iron ion between the bulk solution and the interlayer space.
Conclusions
The simulations shown here demonstrated that the exchange of Na+ by Fe2+ or Fe3+ in montmorillonite interlayers influences the clay’s structural and transport properties significantly, particularly at low water content. Fe ions increase the d-spacing by 1–2 Å compared with Na-Mnt when the interlayer contains less than a bilayer of water, reflecting their stronger hydration. Hydration-energy calculations reveal that Fe2+ and Fe3+ form more stable hydration shells than Na+, with the order Fe2+<Fe3+<Na+. The diffusion of both water and ions decreases in Fe-Mnt, with Fe2+ ions exhibiting much slower mobility than Na+, suggesting that Fe enrichment can impede interlayer transport. Free-energy profiles indicate that both Fe2+ and Fe3+ intercalate preferentially into the clay interlayers, highlighting the thermodynamic favorability of Fe–Na exchange under repository-relevant conditions. These findings imply that Fe-enriched smectites can form near corroded canisters, altering bentonite swelling (e.g. at the same water content as in the interlayer, Fe-exchanged Mnt swells more than Na-exchanged Mnt), CEC, and potentially its long-term performance as a buffer material in nuclear waste repositories.
Author contributions
TAH and CFJ designed the research and wrote the paper. TAH conducted simulations and analyzed data.
Acknowledgments
The authors thank Dr Jeffery Greathouse of Sandia National Laboratories for valuable discussions.
Data availability statement
The authors confirm that the data supporting the findings of this study are available within the article or can be made available through appropriate requests.
Financial support
This work was supported by the DOE Spent Fuel and Waste Science & Technology (SFWST) Program. This article was co-written by an employee of National Technology & Engineering Solutions of Sandia, LLC under Contract No. DE-NA0003525 with the U.S. Department of Energy (DOE). The employee owns all right, title and interest in and to the article and is solely responsible for its contents. The United States Government retains, and the publisher, by accepting the article for publication, acknowledges that the United States Government retains a non-exclusive, paid-up, irrevocable, world-wide licence to publish or reproduce the published form of this article or allow others to do so, for United States Government purposes. The DOE will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan: https://www.energy.gov/downloads/doe-public-access-plan.
Competing interests
The authors declare that they have no competing interests.
