Introduction
Bentonite and argillaceous clay-bearing rocks are considered in the design concepts of several deep geological repositories for the disposal of high level nuclear waste packages as an engineered barrier (Sellin and Leupin, Reference Sellin and Leupin2013) and as the surrounding host rock (Grambow, Reference Grambow2016), respectively. Varying humidity and temperature conditions due to chemical or geologic processes in the near-field could affect the structural properties of swelling clay phases in repository environments. Specifically, smectite clays can undergo intracrystalline expansion or contraction in response to changes in relative humidity (RH). These hydration/dehydration phenomena can lead to changes in clay swelling pressures exerting undue forces on surrounding barrier materials (e.g., waste canister surface) or induce shrinkage crack formation in bentonite material (Soe et al., Reference Soe, Osada, Takahashi and Sasaki2009). This can result in changes in the bulk porosity and permeability of engineered barriers. Therefore, it is important to fully understand the hydration and dehydration processes of expansive clays at the molecular level.
Smectite clays such as montmorillonite (Mnt) are thought to exist in stable hydration states characterized by the number of water layers in between the aluminosilicate layers. As RH increases, smectites undergo stepwise intracrystalline swelling as the d-spacing abruptly increases with increasing number of water layers from zero (0W), one (1W), two (2W), and three (3W) in the interlayers (Bradley et al., Reference Bradley, Grim and Clark1937). More recent X-ray diffraction (XRD) studies have revealed some intermediate peaks rather than distinct peaks for 1W, 2W, or 3W hydration states (Moore and Hower, Reference Moore and Hower1986; Yamada et al., Reference Yamada, Nakazawa, Hashizume, Shimomura and Watanabe1994; Tamura et al., Reference Tamura, Yamada and Nakazawa2000). These intermediate peaks have been interpreted as a result of the random interstratification of two or more distinct hydration states (Sato et al., Reference Sato, Watanabe and Otsuka1992; Sato et al., Reference Sato, Murakami and Watanabe1996; Tamura et al., Reference Tamura, Yamada and Nakazawa2000) and are included in XRD profile modeling (Ferrage et al., Reference Ferrage, Lanson, Michot and Robert2010; Ferrage, Reference Ferrage2016; Chaaya et al., Reference Chaaya, Gaboreau, Milet, Maubec, Tremosa, Raimbourg and Ferrage2023).
A recent molecular dynamics (MD) simulation and experimental study of smectite clay dehydration revealed a two-stage dehydration process (Ho et al., Reference Ho, Coker, Jove-Colon and Wang2022). Based on calorimetric data, the first stage includes a rapid mass loss and partial reduction in crystallographic d spacing. Accompanying MD simulations indicated that water molecules not coordinated to interlayer cations are lost during this first stage. In the second stage, both the rate of mass loss and reduction in d spacing are reduced as water molecules coordinating interlayer cations are removed. Variations in the strength of ion–water and ion–clay interactions have been linked to the hysteresis observed in smectite hydration–dehydration experiments (Fu et al., Reference Fu, Zhang and Low1990; Laird et al., Reference Laird, Shang and Thompson1995), as an explanation for thermo-mechanical irreversibility in smectite swelling behavior (Boek and Coveney, Reference Boek and Coveney1995; Honorio et al., Reference Honorio, Brochard and Vandamme2017; Brochard, Reference Brochard2021). Irregular water loss during smectite dehydration can also result in localized bending of clay layers (Honorio et al., Reference Honorio, Brochard, Vandamme and Lebee2018), particularly near particle edges (Ho et al., Reference Ho, Criscenti and Greathouse2019).
Recent theoretical studies have included clay layer bending related to osmotic swelling and tactoid formation (Whittaker et al., Reference Whittaker, Shoaib, Lammers, Zhang, Tournassat and Gilbert2023), as well as energetic and mechanical aspects of clay layer bending during transitions between thermodynamically stable hydration states (Honorio et al., Reference Honorio, Brochard and Vandamme2017; Honorio et al., Reference Honorio, Brochard, Vandamme and Lebee2018). For example, the transitions from 0W to 1W and from 1W to 2W of Na-Mnt at a confining pressure of 10 MPa and a temperature of 300 K encounter energy barriers of 7 kT nm–2 and 2 kT nm–2, respectively, and the energy barriers decrease with increasing temperature (Honorio et al., Reference Honorio, Brochard and Vandamme2017). These values are consistent with separate simulation studies at room temperature that revealed energy barriers between hydration states in Na-Mnt of ~8−17 kT nm–2 for the 0W–1W transition, decreasing to ~2 kT nm–2 for the 1W–2W transition and <1 kT nm–2 for the 2W–3W transition (Brochard, Reference Brochard2021; Ho et al., Reference Ho, Criscenti and Greathouse2019; Shen and Bourg, Reference Shen and Bourg2021). In the limit of large lateral dimensions for clay particles, Honorio et al. (Reference Honorio, Brochard, Vandamme and Lebee2018) combined plate theory with molecular simulations to show that such hydration energy barriers are overcome by bending of clay layers (i.e. atomic displacement in c) over a relatively small lateral dimension (a or b) of ~1 nm. Thus it is expected that for clay particles with large aspect ratios, the hydration transition involves bending of clay layers in localized areas of the particle, rather than an instantaneous translational displacement of an entire layer.
The vast majority of molecular simulation studies of smectite clays use planar (flat) models of clay layers based on time-averaged structures obtained from diffraction studies of related minerals (Catti et al., Reference Catti, Ferraris, Hull and Pavese1994; Drits et al., Reference Drits, Guggenheim, Zviagina and Kogure2012). Recently, a comprehensive grand canonical Monte Carlo (GCMC) simulation study of water adsorption in Na-Mnt predicted stable and metastable hydration states in good agreement with experiments (Brochard, Reference Brochard2021). However, only standard models consisting of planar (flat) clay layers in stable hydration states (1W–3W) were used. To our knowledge, the hydration behavior of interlayer regions associated with bent clay layers has not been reported previously. Therefore, in this work we compare simulated water adsorption isotherms in Mnt models composed of bent clay layers in mixed hydration states with those consisting of planar clay layers and interstratified hydration states (i.e. 0W–1W). The effect of interlayer cations is included by comparing results for Na- and Ca-exchanged Mnt. Results are also shown for a large-scale MD simulation of a more realistic clay stack that includes dehydrated and hydrated regions in the same interlayer.
Materials and methods
Water adsorption isotherms at 300 K were obtained from GCMC simulations using the Towhee code (Martin, Reference Martin2013) in the grand canonical ensemble. The Clayff parameter set (Cygan et al., Reference Cygan, Liang and Kalinichev2004) together with the rigid extended simple point charge (SPC/E) water model (Berendsen et al., Reference Berendsen, Grigera and Straatsma1987) was used for all Coulombic and van der Waals interactions. Clayff has been used extensively for molecular simulation of clay minerals (Cygan et al., Reference Cygan, Greathouse and Kalinichev2021), including water adsorption in Mnt (Brochard, Reference Brochard2021; Teich-McGoldrick et al., Reference Teich-McGoldrick, Greathouse, Jové-Colón and Cygan2015; Honorio et al., Reference Honorio, Brochard and Vandamme2017). Modified interaction parameters between interlayer Na+ ions and surface O atoms (ob atom type) were used as these modifications resulted in a more realistic description of energy barriers between hydration states (Ho et al., Reference Ho, Criscenti and Greathouse2019). Three-dimensional periodic boundary conditions (3D PBCs) were used with a short-range cut-off distance of 10.0 Å, and long-range electrostatic interactions were evaluated using Ewald summation with a precision of 1×10−4.
The Mnt supercell (Fig. 1) consisted of two clay layers each containing 12 unit cells (4×3×2 repeats in a, b, and c) generated from the pyrophyllite structure (Lee and Guggenheim, Reference Lee and Guggenheim1981). After orthogonalization, the lateral (xy) dimensions of the supercell were 20.8 Å×26.9 Å. A layer charge of −0.83 e/u.c. (elementary charge per unit cell) was introduced in each octahedral sheet by replacing 10 randomly chosen Al atoms with Mg atoms while avoiding Mg-O-Mg close contacts. This layer charge is slightly larger than typically used for Mnt simulations (−0.75 e/u.c.) due to the slightly larger supercell (3b). The layer spacing was expanded to specified values (discussed below), and counterions (Na+ or Ca2+) were randomly inserted in each midplane region for charge balance. The Mnt u.c. formula for Na+ was Na0.833[Si8](Mg0.833Al3.167)O20(OH)4.
Initial 0W–1W models (yz plane) for GCMC simulations of Na-Mnt in the flat (A) and bent (B) configurations. Dimensions are 20.8×26.9×22.0 Å3. Atoms are colored as follows: red (O), white (H), yellow (Si), magenta (Al), green (Mg), and purple (Na). Periodic boundaries are shown in y as black lines (simulations were performed using 3D PBCs). Layer spacings are indicated with arrows.
Adjustments were made to the cell z-parameter and z-coordinate of clay layers so that models corresponded to mixed hydration states (0W–1W, 1W–2W, 2W–3W) with the clay layers in either planar (flat) or bent configurations. Except for the 0W layer spacings, which were determined from separate simulations to be 9.4 Å and 10.0 Å for Ca-Mnt and Na-Mnt, respectively, uniform layer spacing values of 12.0 Å (1W), 15.0 Å (2W), and 18.0 Å (3W) were used to represent hydration states regardless of interlayer cation. The z-dimension of each model system is the sum of the layer spacings of each single hydration state.
Models with flat clay layers were generated by shifting the layers in z with layer spacings corresponding to each stable hydration state (e.g. 1W and 2W). Models with bent clay layers were generated as follows: starting from a layer spacing midway between the stable hydration states (e.g. 13.5 Å for 1W–2W), the z-coordinate of all layer atoms were shifted up or down in z as a function of the atom’s y-coordinate so that the layer spacing varied continuously between the stable hydration states. Using the 0W–1W hydration state as an example (Fig. 1), the layer spacing of the bent configuration varied from 9.4 Å at y=0 to 12.0 Å at the midpoint (y=13.45 Å). The other layer was bent in the opposite direction. In this study bending was arbitrarily chosen to occur only in the y direction, although bending in x is plausible based on similar elastic moduli in lateral dimensions of clay particles (Teich-McGoldrick et al., Reference Teich-McGoldrick, Greathouse and Cygan2012; Zartman et al., Reference Zartman, Liu, Akdim, Pachter and Heinz2010). This type of bent geometry is consistent with a recent theoretical analysis of clay layer bending (Honorio et al., Reference Honorio, Brochard, Vandamme and Lebee2018) which showed that localized bending between hydration states can occur over bending lengths as small as 1 nm. Different degrees of bending would be expected in larger models (see below), but this bending scheme representing two hydration states over a nanometer-scale lateral distance is consistent with conceptual models in the literature (Ferrage et al., Reference Ferrage, Kirk, Cressey and Cuadros2007). Layer spacings and porosities of all models are given in Table 1.
Layer spacing and porosity for GCMC models
a For bent models, values represent the average as layer spacing varies with y (as described in the text). For flat models, layer spacings for smaller and larger interlayers are shown.
b Ratio of free volume to total volume using the Connolly surface method (Connolly, Reference Connolly1983) in BIOVIA Materials Studio (Dassault Systèmes SE) with a probe diameter (3.2 Å) corresponding to a water molecule.
GCMC moves and probabilities included configurational-bias (Siepmann and Frenkel, Reference Siepmann and Frenkel1992) water insertion (20%) and deletion (20%), configurational-bias intrabox water transfer (15%), configurational-bias water molecule regrowth (15%), translation of interlayer water molecules and cations (20%), and rotation of water molecules (10%). No moves involved atoms in the clay lattice as these atoms were held fixed throughout the simulations. Each system was simulated for 2.5×108 steps, with statistics for water loading computed over the final 1.0×108 steps. Values of the water chemical potential (μ) corresponding to RH from 0.5 to 93% were based on previous GCMC simulations of pure SPC/E water (Brochard, Reference Brochard2021).
Large-scale MD simulations at 300 K were performed on Na-Mnt using 3D PBCs in both flat and bent configurations. The simulation was performed with the LAMMPS code (Thompson et al., Reference Thompson, Aktulga, Berger, Bolintineanu, Brown, Crozier, ’tVeld, Kohlmeyer, Moore, Nguyen, Shan, Stevens, Tranchida, Trott and Plimpton2022) using Clayff parameters with the flexible SPC water model (Teleman et al., Reference Teleman, Jonsson and Engstrom1987). Temperature and pressure were controlled with a Nose–Hoover thermostat and barostat with relaxation times of 100 and 1000 fs, respectively. A timestep of 1 fs was used with a 10.0 Å cut-off distance for short-range interactions and a particle-mesh Ewald algorithm (Plimpton et al., Reference Plimpton, Pollock and Stevens1997) with a precision of 1×10−4 for long-range Coulomb interactions. To construct the simulation cell, an orthogonalized Mnt layer with xy dimension 4.16×3.58 nm2 was created in a 4.16×3.58×1.8 nm3 simulation box (8×4×1 u.c. repeats in a, b, and c) as shown in Fig. 2A. The layer charge of the Mnt layer was −0.75 e/u.c. This structure was then replicated 10×2×20 times to create a large clay structure containing 521,600 atoms (Fig. 2B). Water molecules were then randomly inserted in each interlayer, resulting in a water content of 2.34 H2O/u.c. (1W). After minimization, the system was equilibrated for 50 ps at constant volume followed by 1 ns at a constant pressure (Pzz) of 1 atm to equilibrate the layer spacing. Here, Pzz refers to pressure (stress) applied normal to the basal surfaces, so that volume changes only involved changes in the supercell c parameter. The final supercell with dimensions 41.6×71.6×23.8 nm3 is shown in Fig. 2C. The average layer spacing was ~11.9 Å, confirming the monolayer (1W) hydration state. Water uniformly distributes in each interlayer and each clay layer is flat.
Model preparation for large-scale MD simulation of Na-Mnt with heterogeneous interlayers. (A) Initial simulation box with one clay layer (gray) and Na+ ions (yellow). (B) Large clay structure by replicating clay structure in A. (C) Homogeneous 1W configuration (water O red) with dimensions ~ 41.6×7.1×23.3 nm3.
To create a bent clay configuration, the simulation box in Fig. 2A was filled randomly with water molecules (i.e. 2.36 H2O/u.c.) and five Lennard-Jones (LJ) particles (
$ \unicode{x025B} $
=0.294076 kcal mol–1,
$ \unicode{x03C3} $
=8.1656 Å). The large
$ \unicode{x03C3} $
value is required to locally expand the interlayer to accommodate the 2W state. This system was then replicated 10×2×20 times to create a large structure similar to that shown in Fig. 2C. The large structure was then equilibrated for 10 ps at constant volume, followed by 500 ps at constant pressure (Pzz=100 atm), and finally 2 ns at constant pressure (Pzz=1 atm). As seen in Fig. 3, the LJ particles aggregated in the interlayer, acting as a spacer that caused localized bending of the clay layers. Note that interactions between water molecules and LJ particles were turned off, so water molecules can coexist with the LJ particle aggregates. The LJ particles were then removed, and the system was equilibrated during a 100 ps constant pressure simulation (Pzz=100 atm) followed by 3 ns Pzz=1 atm. The final box dimensions are 41.6×71.6×23.4 nm3.
Equilibrium snapshot from a large-scale MD simulation of Na-Mnt showing interlayers containing spacer particles (pink spheres) to create the bent structure.
XRD patterns were computed from the MD-equilibrated large-scale models in both flat and bent configurations. The positions of atoms obtained from MD simulation were used to compute the XRD intensity I (Coleman et al., Reference Coleman, Spearot and Capolungo2013):
$$ I= Lp\left(\unicode{x03B8} \right)\frac{F\ast F}{N}, $$
where
$ Lp\left(\unicode{x03B8} \right) $
is the Lorentz-polarization factor,
$ \unicode{x03B8} $
is the scattering angle of diffraction, N is the number of atoms. The structural factor F was calculated as:
$$ F(k)=\sum \limits_{\mathrm{j}=1}^N{f}_{\mathrm{j}}\left(\theta \right)\exp \left(2\pi ik{r}_{\mathrm{j}}\right) $$
where k is the location of the reciprocal lattice node, f
j is the atomic scattering factor (Peng et al., Reference Peng, Ren, Dudarev and Whelan1996), and r
j is the position of each atom obtained from the MD simulation. The XRD reflection was calculated with the spacing of the reciprocal lattice node of 0.1 Å in the z direction and the wavelength is 1.541838 Å. The XRD profiles were recorded from 3 to 33° 2
$ \unicode{x03B8} $
with a resolution of 0.1° 2
$ \unicode{x03B8} $
. Note that the XRD calculation was computationally expensive.
Results and Discussion
The effect of clay-layer bending on interlayer hydration was evaluated from GCMC simulations. Water adsorption isotherms for Na-Mnt and Ca-Mnt in mixed hydration states are shown in Fig. 4. As expected, water loading increases with accessible volume (i.e. 0W–1W<1W–2W<2W–3W; see Table 1 for porosity). Curvature of the isotherm profiles as a function of RH is consistent with previous GCMC simulations of water adsorption in stable hydration states (i.e. 1W, 2W, 3W) (Hensen et al., Reference Hensen, Tambach, Bliek and Smit2001; Smith et al., Reference Smith, Wang, Chaturvedi and Whitley2006; Ferrage et al., Reference Ferrage, Sakharov, Michot, Delville, Bauer, Lanson, Grangeon, Frapper, Jiménez-Ruiz and Cuello2011; Dazas et al., Reference Dazas, Lanson, Delville, Robert, Komarneni, Michot and Ferrage2015; Teich-McGoldrick et al., Reference Teich-McGoldrick, Greathouse, Jové-Colón and Cygan2015). Specifically, the RH required to achieve maximum water loading increases with accessible interlayer volume. Models with the smallest accessible volume (0W–1W) achieve maximum loading at very low RH, while models with the largest accessible volume (2W–3W) do not achieve maximum water loading until ~60% RH. However, isotherms for the Ca-Mnt models for these states show a more gradual increase in water loading at low RH, because there are fewer interlayer cations that serve as initial adsorption sites.
Water adsorption isotherms from GCMC simulations at 300 K for mixed hydration states of (A) Na-Mnt and (B) Ca-Mnt in bent and flat configurations.
Models with bent clay layers (representative of transition hydration states) show similar or identical water adsorption behavior compared with the standard models (flat clay layers with single hydration states). However, bent models with the smallest layer spacing (0W–1W) show slightly lower water loading than the corresponding flat model at all RH values. As seen in Table 1, the 0W–1W bent models have significantly lower porosities than the flat models. For the 1W–2W and 2W–3W systems, the bent and flat models have similar porosities and water loadings at all RH.
GCMC snapshots for the 0W–1W models (Fig. 5) illustrate the effect of reduced porosity on water adsorption. A large proportion of the bent interlayer region has little or no adsorbed water, indicating that there is not enough accessible volume for a water layer to form. Additional snapshots from all GCMC simulations are shown in the Supplementary material (Figs S1–S3). All of the interlayer regions are accessible to water in the expanded states (1W–2W and 2W–3W), so clay-layer bending has little effect on water loading at higher RH when these hydration states are thermodynamically stable (Brochard, Reference Brochard2021).
Equilibrium 0W–1W snapshots (yz plane) of Na-Mnt in the flat (A) and bent (B) configurations. Atoms are colored as follows: red (O), white (H), yellow (Si), magenta (Al), green (Mg), and purple (Na).
As further evidence of the minimal impact of clay-layer bending on montmorillonite structure, a comparison of lattice bond lengths and angles (see Supplementary material, Figs. S4 and S5) shows nearly identical distributions for bent and flat models. The small degree of layer bending in this study results in a small vertical (z) displacement of 1.0–1.5 Å over a 1.3 nm lateral dimension, which has no discernable effect on the bond or angle distributions. Of course, more significant bending of the layers would cause more measurable increases in bonds and angles, resulting in lattice strain. In addition, short 50 ps MD simulations with identical water loadings in each bent/flat pair (see Supplementary material, Fig. S6) show that layer bending has no effect on potential energy in the 1W–2W and 2W–3W models and only increases slightly in the potential energy in the 0W–1W models. Note that these potential energy comparisons are limited to the interlayer species as the clay layers are held fixed.
This depiction of variable layer spacing within a single interlayer is consistent with a recent study of bending moduli for clay layers (Honorio et al., Reference Honorio, Brochard, Vandamme and Lebee2018), which showed that isolated clay layers are quite flexible even during stirring in an exfoliated state, with bending moduli in the lateral dimensions ~2×10−17 N m. For stacks of clay particles with large aspect ratios, the authors combined plate theory with molecular simulations to show that hydration energy barriers (<10 kT nm–2) are overcome by bending of clay layers over a relatively small bending length of ~1 nm (Honorio et al., Reference Honorio, Brochard, Vandamme and Lebee2018).
The bent GCMC models reformed their original flat configurations during subsequent constant-Pzz MD simulations (1 atm), resulting in homogeneous interlayers with layer spacings intermediate between stable hydration states. However, the large bent model (41.6×71.6×23.8 nm3) maintained its configuration throughout the final 3 ns simulation (Fig. 6). Removal of the LJ spacer particles in Fig. 3 resulted in some local rearrangement of clay layers and interlayer water, but dry interlayer regions are clearly visible. Interlayer regions corresponding to 0W, 1W, and even 2W can exist within a single interlayer of the structure. Note that the average water content in each interlayer was the same (i.e. 2.36 H2O/u.c.). Results from the large-scale MD simulation indicate that (1) the persistence of bent clay layers depends on simulation system size, and (2) localized layer bending in larger clay stacks appears to be long lived. Similar layer bending was also seen during large-scale simulations of organically modified smectites with interlayer ethanol (Metz et al., Reference Metz, Anderson, Geatches, Suter, Lines and Greenwell2015). Additionally, thermal undulations in layered double hydroxides seen in large-scale MD simulations enabled the calculation of bending moduli (Thyveetil et al., Reference Thyveetil, Coveney, Suter and Greenwell2007).
Equilibrium snapshot from a large-scale MD simulation of Na-Mnt showing heterogeneous 0W–1W–2W interlayers after removal of spacer particles.
In Fig. 7 we compare the XRD profile for the bent (heterogeneous) model (Fig. 6) with that of the flat homogeneous model (Fig. 2C; 2.34 H2O/u.c.). The results indicate that heterogeneity in the same interlayer at the length scale studied here does not result in any significant difference in the (001) XRD reflection, but such heterogeneity does result in a broadening of the peaks at larger values of 2θ. This broadening has been attributed to crystal strain due to variability in atomic z coordinates in XRD samples (i.e. layer bending), and can be quantified with a variance parameter sz (Ferrage et al., Reference Ferrage, Lanson, Sakharov and Drits2005). Importantly, larger values of sz are associated with incomplete hydration transitions, which is accommodated by layer flexibility (Ferrage et al., Reference Ferrage, Lanson, Sakharov and Drits2005).
Comparison of XRD profiles from large-scale MD simulations of Na-Mnt with heterogeneous (black) and homogeneous (red) interlayer hydration.
It is possible that water distribution in the same interlayer is not uniform (Wang and Xu, Reference Wang and Xu2006). A conceptual model of smectite dehydration based on kinetic analysis of XRD patterns at different temperatures includes localized layer bending as water molecules diffuse out of the interlayers (Ferrage et al., Reference Ferrage, Kirk, Cressey and Cuadros2007). Similarly, analysis of X-ray scattering and cryogenic microscopy indicates a two-step dehydration process involving diffusional mixing and localized layer bending (Whittaker et al., Reference Whittaker, Lammers, Carrero, Gilbert and Banfield2019). As discussed above, XRD patterns of heterogeneous hydration structures are characterized by irrational, asymmetrical, and broad reflections (Ferrage et al., Reference Ferrage, Lanson, Sakharov and Drits2005). The issue of heterogeneous hydration states in smectites becomes more relevant when considering continuous hydration processes (i.e. not stepwise swelling mechanism) of natural smectite due to heterogeneous charge distribution on the clay layers (Chipera et al., Reference Chipera, Carey, Bish, Gilfrich, Noyan, Jenkins, Huang, Snyder, Smith, Zaitz and Predecki1995; Tamura et al., Reference Tamura, Yamada and Nakazawa2000; Morodome and Kawamura, Reference Morodome and Kawamura2011; Schleicher et al., Reference Schleicher, Hofmann and van der Pluijm2013).
Conclusions
We used molecular simulations to compare water adsorption isotherms for Mnt models corresponding to transitions between standard hydration states. Hydration transitions may cause bending and undue forces in the near field of clay phases, so it is important to quantify such effects in relevant scenarios such as underground repository settings.
GCMC simulations on small model systems (nanometer dimensions) were performed to compare the effect of layer bending on water adsorption. For models with expanded interlayers (1W–2W, 2W–3W), porosities and simulated water adsorption isotherms were nearly identical between models with bent clay layers representing mixed hydration states and those with flat clay layers with (interstratified) standard hydration states. Bent models involving dry states (0W–1W) have reduced porosities and water adsorption compared with the corresponding flat models. Similar results were obtained for models with monovalent (Na+) and divalent (Ca2+) interlayer cations, indicating that ion hydration properties have a minor influence layer bending, within the models studied. However, an analysis of (001) peak irrationality suggests that more fluctuations in layer thickness in the presence of divalent cations may be due to the lower cation density in the interlayer (Ferrage et al., Reference Ferrage, Lanson, Sakharov and Drits2005).
Bond lengths and angles in the bent models with a small vertical displacement (up to 1.5 Å) are essentially unchanged compared with the standard flat models. Although the small models used in the GCMC simulations do not retain the bent configuration during unconstrained MD simulation, large-scale simulations (with dimensions similar to actual clay particles) indicate that models with localized bending of clay layers are stable over MD timescales. The tendency for bent layers to flatten clearly depends on system size. Consistent with a previous analysis of mechanical properties of clay-layer bending (Honorio et al., Reference Honorio, Brochard, Vandamme and Lebee2018), our results indicate that stacks of Mnt layers are able to exist in mixed hydration states (and localized bending) with minimal effect on macroscopic clay structure or water adsorption properties. Additional research is needed to confirm the coexistence of these configurations under hydrating conditions. Finally, this study has established techniques for creating model systems at various length scales which can be used to explore the effect of layer bending on smectite properties.
Supplementary material
The supplementary material for this article can be found at http://doi.org/10.1017/cmn.2024.22.
Data availability statement
The authors confirm that the data supporting the findings of this study are available within the article and its Supplementary material, or may be made available through appropriate requests.
Acknowledgements
None
Author contribution
J.A.G.: Conceptualization, Methodology, Simulations, Writing, Visualization. T.A.H.: Conceptualization, Methodology, Simulations, Writing, Visualization. C.F.J–C.: Conceptualization, Resources, Writing–original draft.
Financial support
This work was supported by the U.S. Department of Energy Office of Nuclear Energy, through the Office of Spent Fuel and Waste Science and Technology (SFWST) Research and Development Campaign. This work was produced by National Technology & Engineering Solutions of Sandia, LLC (NTESS) for the U.S. Department of Energy under contract DE-NA0003525. This work is authored and owned by an employee of NTESS and is responsible for its contents. Any subjective views or opinions expressed in the written work do not necessarily represent the views of the U.S. Government. The publisher acknowledges the U.S. Government license to provide public access under the DOE Public Access Plan.
Competing interests
The authors declare that they have no competing interests.
