Aiming to identify the complexing mechanisms of heavy metal cations on edge surfaces of 2:1-type clay minerals, systemic first-principles molecular dynamics (FPMD) simulations were conducted and the microscopic structures and complex free energies were obtained. Taking Cd(II) as a model cation, the structures on both (010) and (110) edges of the complexes were derived for the three possible binding sites (≡SiO, ≡Al(OH)2/≡AlOH≡AlSiO, and vacant sites). The stable complexes adsorbed on the three binding sites on both terminations had similar structures. The free energies of the complexes on (010) edges were calculated by using the constrained FPMD method. The free energies of complexes on the ≡SiO and ≡Al(OH)2 sites were similar and they were both significantly lower than the free energy of the complex on the octahedral vacant site. In association with the concept of high energy site (HES) and low energy site (LES) in the 2 Site Protolysis Non Electrostatic Surface Complexation and Cation Exchange (2SPNE SC/CE) sorption model, the vacant site was assigned to HES and the other two sites to LES, respectively.

Investigation of the organization of interlayer water and cations in smectite is a permanent topic in clay science for environmental science, civil engineering, materials science, and industrial applications. Experimental X-ray (or neutron) diffraction methods and molecular simulations are key techniques to probe the organization of the smectite structure at a similar molecular length scale. The combination of both of these experimental and numerical methods represents a complementary approach to reveal the structural heterogeneity of real samples, design and model a wide range of smectite structures, and validate the simulation results through comparison with experimental data.

This paper first revisits establishment of the original interlayer model as developed in the 1930s for the organization of water and ions in the smectite structure using X-ray diffraction (XRD) techniques. Then, based on a simplified approach, key theoretical tools are provided to calculate XRD pattern 00l reflections for a periodic smectite structure with a wide range of interlayer compositions and organizations using conventional spreadsheet software. In addition to educational purposes, this theoretical description is used to describe the principal parameters governing the positions and intensities of experimental XRD 00l reflections. This calculation toolbox is also used to determine better the layer-to-layer distances considered in molecular simulations and to validate these simulations through a detailed collation procedure using experimental data.

Recent examples of the application of such a procedure to collate experimental diffraction data and molecular simulations are presented for the specific case of deciphering the molecular organization of interlayer water and cations in the different smectite hydrates (mono-, bi-, and tri-hydrated layers). The extension of this approach to the interlayer refinement of organo-clays is also detailed, and perspectives regarding the characterization of other lamellar compounds are discussed.

The aqueous chemistry of water films confined between clay mineral surfaces remains an important unknown in predictions of radioelement migration from radioactive waste repositories. This issue is particularly important in the case of long-lived anionic radioisotopes (129I-, 99TcO4-, 36Cl-) which interact with clay minerals primarily by anion exclusion. For example, models of ion migration in clayey media do not agree as to whether anions are completely or partially excluded from clay interlayer nanopores. In the present study, this key issue was addressed for Cl- using MD simulations for a range of nanopore widths (6 to 15 Å) overlapping the range of average pore widths that exists in engineered clay barriers. The MD simulation results were compared with the predictions of a thermodynamic model (Donnan Equilibrium model) and two pore-scale models based on the Poisson-Boltzmann equation under the assumption that interlayer water behaves as bulk liquid water. The simulations confirmed that anion exclusion from clay interlayers is greater than predicted by the pore-scale models, particularly at the smallest pore size examined. This greater anion exclusion stems from Cl- being more weakly solvated in nano-confined water than it is in bulk liquid water. Anion exclusion predictions based on the Poisson-Boltzmann equation were consistent with the MD simulation results, however, if the predictions included an ion closest approach distance to the clay mineral surface on the order of 2.0 ± 0.8 Å. These findings suggest that clay interlayers approach a state of complete anion exclusion (hence, ideal semi-permeable membrane properties) at a pore width of 4.2 ± 1.5 Å.

A better understanding of the thermodynamics of radioactive cesium uptake at the surfaces of phyllosilicate minerals is needed to understand the mechanisms of selective adsorption and help guide the development of practical and inexpensive decontamination techniques. In this work, molecular dynamics simulations were carried out to determine the thermodynamics of Cs+ adsorption at the basal surface of six 2:1 phyllosilicate minerals, namely pyrophyllite, illite, muscovite, phlogopite, celadonite, and margarite. These minerals were selected to isolate the effects of the magnitude of the permanent layer charge (⩽2), its location (tetrahedral vs. octahedral sheet), and the octahedral sheet structure (dioctahedral vs. trioctahedral). Good agreement was obtained with the experiments in terms of the hydration free energy of Cs+ and the structure and thermodynamics of Cs+ adsorption at the muscovite basal surface, for which published data were available for comparison. With the exception of pyrophyllite, which did not exhibit an inner-sphere free energy minimum, all phyllosilicate minerals showed similar behavior with respect to Cs+ adsorption; notably, Cs+ adsorption was predominantly inner-sphere, whereas outer-sphere adsorption was very weak with the simulations predicting the formation of an extended outer-sphere complex. For a given location of the layer charge, the free energy of adsorption as an inner-sphere complex varied linearly with the magnitude of the layer charge. For a given layer charge location and magnitude, adsorption at phlogopite (trioctahedral sheet structure) was much less favorable than at muscovite (dioctahedral sheet structure) due to electrostatic repulsion between adsorbed Cs+ and the H atom of the OH- ion directly below the six-membered siloxane ring cavity. For a given layer charge magnitude and octahedral sheet structure, adsorption to celadonite (octahedral sheet layer charge) was favored over adsorption to muscovite (tetrahedral sheet layer charge) due to the increased distance to the surface K+ ions and the decreased distance to the O atom of the OH- ion directly below the surface cavity.

A detailed structural characterization of organo-clays is a key in understanding their properties. In this work, mono-, di-, tri-, and tetra-butylammonium (nBA; n = 1–4) cations intercalated in the layered clay mineral montmorillonite (Mnt) have been studied for the first time by combining a theoretical approach based on density functional theory (DFT) and infrared spectroscopy. The DFT calculations revealed the detailed structure and position of nBA cations in the interlayer space. A relation between the basal spacing (d001 parameter) and the cation size and structure was found, and explained with respect to the structure, composition, and size of the organic cations. Hydrogen bonds between -NHx/-CH3/-CH2 groups of the nBA cations and oxygen atoms of the basal planes of the Mnt layers were found to be an important factor for the arrangement and energetic stabilization of cations in the interlayer space. The N-H-O hydrogen bonds are stronger than C-H-O hydrogen bonds and the stabilization decreases with decreased number of bands. Analysis of DFT-calculated vibrational modes helped in understanding a problematic region of the experimental infrared spectra (4000–3000 cm-1), in which assignment of all vibrational modes unambiguously was not possible because of a significant overlap of broad bands.

Some aspects of the crystal structure of illite are not understood properly yet, in spite of its abundance and significance as a component of soils, sediments, and low-grade metamorphic rocks. The present study aimed to explore the role of hydronium cations in the interlayer space of illite in a theoreticalexperimental approach in order to clarify previous controversial reports. The infrared spectroscopy of this mineral has been studied experimentally and by means of atomistic calculations at the quantum mechanical level. The tetrahedral charge is critical for the stability of the hydronium cations, the presence of which has probably been underestimated in previous studies. In the present study, computational studies have shown that the hydronium cations in aqueous solutions are likely to be intercalated in the interlayer space of illite, exchanging for K cations. During the drying process these cations are stabilized by hydrogen bonds in the interlayer space of illite.

The results from mesoscale simulations of the formation and evolution of microstructure for assemblies of Na-smectite particles based on assumed size distributions of individual clay platelets are presented here. The analyses predicted particle arrangements and aggregation (i.e. platelets linked in face—face configurations) and are used to link geometric properties of the microstructure and mechanical properties of the particle assemblies. Interactions between individual ellipsoidal clay platelets are represented using the Gay-Berne potential based on atomistic simulations of the free energy between two Na-smectite clay-platelets in liquid water, following a novel coarse-graining method developed previously. The current study describes the geometric (aggregate thickness, orientation, and porosity) and elastic properties in the ‘jammed states’ from the mesoscale simulations for selected ranges of clay particle sizes and confining pressures. The thickness of clay aggregates for monodisperse assemblies increases (with average stack thickness consisting of n = 3–8 platelets) with the diameter of the individualclay platelets and with the level of confining pressure. Aggregates break down at high confining pressures (50–300 atm) due to slippage between the platelets. Polydisperse simulations generate smaller aggregates (n = 2) and show much smaller effects of confining pressure. All assemblies show increased order with confining pressure, implying more anisotropic microstructure. The mesoscale simulations are also in good agreement with macroscopic compression behavior measured in conventional 1-D laboratory compression tests. The mesoscale assemblies exhibit cubic symmetry in elastic properties. The results for larger platelets (D = 1000 Å) are in good agreement with nano-indentation measurements on natural clays and shale samples.

Sorption of U(VI) on clay and related minerals has been inspected experimentally and computationally because of its central role in safety considerations of geological repositories for highly radioactive waste. Np(V), which also has long half-life isotopes, has received considerably less attention. The purpose of the present study was to investigate computationally the adsorption of Np(V) on a clay-mineral surface and to compare it to adsorption of U(VI). As a sample case study, adsorption of Np(V) at the (110) edge surface of the common clay mineral montmorillonite was modeled. The density functional approach and periodic supercell models were applied. Mono- and bidentate adsorption complexes with coordination numbers 4 and 5 were inspected and compared to corresponding U(VI) species. While U(VI) prefers bidentate adsorption complexes with varying coordination numbers, Np(V) is more stable when monodentate-coordinated with a coordination number of four. In line with its smaller hydrolysis constant in aqueous solution, Np(V) shows a lower tendency to form monohydroxides on the mineral surface compared to U(VI). As no experimental geometry parameters are available for Np(V) adsorbed on montmorillonite, the results were compared tentatively to EXAFS data for adsorption at kaolinite and good agreement for the geometry changes due to adsorption was found for the more preferred adsorbed species.

Experimental studies have shown that a sharp, high-frequency IR band at ~3615 cm-1 (in H2O form) and at ~2685 cm-1 (in D2O form) is a common feature for all smectites, and its position correlates with layer charge. In order to explain the molecular origin of this band in terms of total layer charge, charge localization, as well as nature of interlayer cations influencing the position and intensity of this peak, a series of classical molecular dynamics (MD) simulations was performed for several smectite models. The smectite layers were described using a modified CLAYFF force field, where the intramolecular vibrations of H2O were described more accurately by the Toukan-Rahman potential. The power spectra of molecular vibrations of hydrogens were calculated for selected sub-sets of interlayer H2O to analyze quantitatively their contribution to the observed spectral features. The statistics of hydrogen bonds in the smectite interlayers were also analyzed to support the spectral calculations.

The simulation results demonstrated clearly that only the H2O molecules in close proximity to the smectite surface are responsible for the sharp vibrational band observed. Other hypotheses for the possible origins of this band were considered carefully and eventually rejected. Two orientations of H2O molecules donating one or two H bonds to the basal oxygens of the smectite surface (monodentate and bidentate orientations, respectively) were observed. In both orientations, these H bonds are quite weak, pointing to a generally hydrophobic character of the smectite surface. Both orientations contributed to the high-frequency band, but the monodentate orientation provided the predominant contribution because surface H2O molecules in this orientation were much more abundant. In good agreement with experiment, only a small difference in the peak position was observed between smectites with different charge localization. The effect of the total layer charge, i.e. the red-shift for higher-charge smectites, was also confirmed. This shift arose from the decrease in the H-bonding distances of H2O in monodentate and bidentate orientation.

Many important properties of clay minerals are defined by the species of charge-balancing cation. Phenomena such as clay swelling and cation exchange depend on the cation species present, and understanding how the cations bind with the mineral surface at a fundamental level is important. In the present study the binding affinities of several different charge-balancing cations with the basal surface of the smectite mineral, montmorillonite, have been calculated using molecular dynamics in conjunction with the well-tempered metadynamics algorithm. The results follow a Hofmeister series of preferred ion adsorption to the smectite basal surfaces of the form:

K+ > Na+ > Ca2+ > Cs+ > Ba2+

The results also revealed the energetically favorable position of the ions above the clay basal surfaces. Key features of the free-energy profiles are illustrated by Boltzmann population inversions and analyses of the water structures surrounding the ion and clay surface. The results show that weakly hydrated cations (K+ and Cs+) preferentially form inner-sphere surface complexes (ISSC) above the ditrigonal siloxane cavities of the clay, while the more strongly hydrated cations (Na+) are able to form ISSCs above the basal O atoms of the clay surface. The strongly hydrated cations (Na+, Ca2+, and Ba2+), however, preferentially form outer-sphere surface complexes. The results provide insight into the adsorption mechanisms of several ionic species on montmorillonite and are relevant to many phenomena thought to be affected by cation exchange, such as nuclear waste disposal, herbicide/pesticide-soil interactions, and enhanced oil recovery.

Organo-clays represent a special challenge for molecular simulations because they require accurate representation of the clay and the organic/aqueous sections of the model system and accurate representation of the interactions between them. Due to the broad range of force-field models available, an important question to ask is which sets of parameters will best suit the molecular modeling of the organo-intercalated smectites? To answer this question, the structure of the ethylene glycol (EG)-smectite complex is used here as a testing model because the intercalation of EG in smectites provides a stable interlayer complex with relatively constant basal spacing.

Three smectite samples with substantially different layer charge and charge localization were selected for X-ray diffraction (XRD) measurements. Their molecular models were built and molecular-dynamics simulations performed using various combinations of the organic force fields (CGenFF, GAFF, CVFF, and OPLS-aa) with ClayFF and INTERFACE force fields used to describe smectites. The simulations covered a range of different EG and water contents. For every structure, the density distribution of interlayer species along the direction perpendicular to the layer plane was calculated and then used to optimize the XRD patterns for these simulated models.

A comparison of these results with experimental XRD patterns shows very large discrepancies in the structures and basal spacings obtained for different layer charges as well as for different force fields and their combinations. The most significant factor affecting the accuracy of the calculated XRD patterns was the selection of the clay-mineral force-field parameters. The second important conclusion is that a slight modification of the basal oxygen parameters for non-electrostatic interactions (increase of their effective atomic diameters) may be a simple and straightforward way to improve significantly the agreement between the modeled XRD patterns with experiments, especially for high-charge smectites. Generally, among organic force fields, the least accurate results were obtained with CGenFF. For unmodified ClayFF, its combination with GAFF gave the best results, while the two other sets (OPLS-aa and CVFF) gave the best results in combination with ClayFFmod. The INTERFACE and INTERFACEmod produced much better results for low-charge montmorillonite than for high-charge smectites.

Interlayer swelling of hydrated montmorillonite is an important issue in clay mineralogy. Although the swelling behavior of montmorillonite under ambient conditions has been investigated comprehensively, the effects of basin conditions on the hydration and swelling behaviors of montmorillonite have not been characterized thoroughly. In the present study, molecular dynamics simulations were employed to reveal the swelling behavior and changes in the interlayer structure of Na-montmorillonite under the high temperatures and pressures of basin conditions. According to the calculation of the immersion energy, the monolayer hydrate becomes more stable than the bilayer hydrate at a burial depth of 7 km (at a temperature of 518 K and a lithostatic pressure of 1.04 kbar). With increasing burial depth, the basal spacings of the monolayer and bilayer hydrates change to varying degrees. The density-distribution profiles of interlayer species exhibit variation in the hydrate structures due to temperature and pressure change, especially in the structures of bilayer hydrate. With increasing depth, more Na+ ions prefer to distribute closer to the clay layers. The mobility of interlayer water and ions increases with increasing temperature, while increasing pressure caused the mobility of these ions to decrease.