Water and Ion Dynamics in Confined Media: A Multi-Scale Study of the Clay/Water Interface

Abstract

This review details a large panel of experimental studies (Inelastic Neutron Scattering, Quasi-Elastic Neutron Scattering, Nuclear Magnetic Resonance relaxometry, Pulsed-Gradient Spin-Echo attenuation, Nuclear Magnetic Resonance Imaging, macroscopic diffusion experiments) used recently to probe, over a large distribution of characteristic times (from pico-second up to days), the dynamical properties of water molecules and neutralizing cations diffusing within clay/water interfacial media. The purpose of this review is not to describe these various experimental methods in detail but, rather, to investigate the specific dynamical information obtained by each of them concerning these clay/water interfacial media. In addition, this review also illustrates the various numerical methods (quantum Density Functional Theory, classical Molecular Dynamics, Brownian Dynamics, macroscopic differential equations) used to interpret these various experimental data by analyzing the corresponding multi-scale dynamical processes. The purpose of this multi-scale study is to perform a bottom-up analysis of the dynamical properties of confined ions and water molecules, by using complementary experimental and numerical studies covering a broad range of diffusion times (between pico-seconds up to days) and corresponding diffusion lengths (between Angstroms and centimeters). In the context of such a bottom-up approach, the numerical modeling of the dynamical properties of the diffusing probes is based on experimental or numerical investigations performed on a smaller scale, thus avoiding the use of empirical or fitted parameters.

1. Introduction

Among colloids, clays are ubiquitous materials exhibiting a large variety of physico-chemical properties, such as adsorption of water and polar molecules [1,2,3], high specific surface area [3], electric surface charge [3], swelling [4], gelling [5], thixotropy, and surface acidity exploited in numerous industrial processes, including waste storing [6,7], cosmetic and paint industries, drilling, cracking, and heterogeneous catalysis [8]. Optimizing these various applications requires a detailed analysis of the impact of the clay composition on the dynamical processes implied in these applications. For that purpose, numerous dynamical studies have been performed to probe, over a broad dynamical range, the mobility of ions [9,10] and molecules [11,12,13] physisorbed at the clay/water interface. As illustrated in Figure 1, we selected to use a complementary set of experiments, including Inelastic Neutron Scattering (INS) [14,15,16,17,18,19], Quasi-Elastic Neutron Scattering (QENS) [12,20,21,22,23,24,25,26,27,28,29,30], Nuclear Magnetic Resonance (NMR) relaxometry [31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49], Pulsed-Gradient Spin-Echo (PGSE) attenuation [11,50,51,52], Magnetic Resonance Imaging (MRI) [51], and macroscopic diffusion experiments as Out-Diffusion (OD) [53,54] or Through-Diffusion (TD) [55] experiments, because they probe complementary time-scales (see Figure 1). These experiments allow investigating a broad dynamical range (from few pico-seconds to few days) and quantifying the impact of the clay/water interface on the dynamics of the interfacial molecular and ionic probes.

Furthermore, as illustrated in Figure 2, the use of swelling clays, such as Na⁺-saturated montmorillonite, allows probing a broad range of interlayer separations [4] varying between a limited number of hydration layers for layer-to-layer distances smaller than 20 Å, in the so-called crystalline swelling regime [2,56,57], and a large continuum of layer-to-layer distances (between 40 Å and 70 Å), in the so-called osmotic swelling regime. By contrast, Ca²⁺-saturated montmorillonite only exhibits crystalline swelling (see Figure 2) because inter-ionic correlation forces strongly reduce the swelling pressure of these charged interfaces [58,59]. For that reason, these dynamical experiments require the use of purified clay samples to avoid coexistence of swollen and unswollen phases within heterogeneous samples.

For that purpose, some of the dynamical studies presented here were performed with synthetic clays (laponite and saponite) in order to avoid impurities within the clay network. In the case of natural clays (kaolinite, beidellite, and vermiculite), the samples were purified according to standard procedures [60] and size separated by centrifugation. Under such conditions, clay/water interfaces are well characterized systems and ideal to investigate the influence of charged interfaces on the thermodynamical and dynamical properties of confined fluids.

Finally, numerical modeling is required to relate the organization of these clay/water interfaces and the dynamical information extracted from the various experimental data obtained by INS, QENS, PGSE, MRI, and macroscopic diffusion experiments (see Figure 1). In that context, an atomic description of the clay/water interface is required for modeling the structural and dynamical properties of these interfacial systems probed in the crystalline swelling regime. For that purpose, quantum calculations are used to directly determine the dynamical properties of confined ions and water molecules by using ab initio DFT quantum methods [61,62]. Quantum calculations are also useful to evaluate atomic properties of the clay network, such as atomic partial charges, by using semi-empirical [63] or ab initio [64] quantum calculations. These parameters are further exploited in classical statistical treatments of the clay/water interfaces by using classical Monte Carlo (MC) simulations and Molecular Dynamics (MD) to investigate the structural and dynamical properties of the confined ions and water molecules within the crystalline swelling regime [2,12,17,28,40,41,46,63,64,65,66]. By contrast, in the so-called osmotic regime, one may easily neglect the atomic structure of the clay network and the solvent molecules by using classical Brownian Dynamics (BD) in order to determine the effective mobility [67] of the neutralizing counterions in the presence of structureless charged surfaces [68]. In the same context, numerical integration of differential equations can also be used to study the mobility of diffusing probes at the macroscopic scale (Refs. [54,69], and references therein).

In the present review, these techniques are described following a path of increasing probed time-scales, together with the associated numerical methods, in order to highlight the methodology used to get information regarding the multi-scale properties of clay/water interfaces. The purpose of that study is to exploit such a multi-scale analysis of the dynamical properties of ions and water molecules diffusing within clay sediments in order to develop a bottom-up approach of such complex interfacial dynamical processes. In such analysis, we avoid to use fitted parameters derived from empirical models by performing realistic numerical modeling exploiting structural or dynamical properties, such as average residence time or local mobility, extracted from preliminary studies performed on a smaller scale.

2. Multi-Scale Approaches for Studies of Water and Ions Dynamics in Clay Systems

2.1. Inelastic Neutron Scattering

The intensity of the INS spectra is directly proportional to the so-called Generalized Density of States, labeled $GDOS\left(\omega \right)$, that is the density of states of a system taking into account the scattering cross section of the various atoms:

$$GDOS\left(\omega \right) = \sum _j\frac{{\sigma }_j}{M_j}{g}_j\left(\omega \right),$$

(1)

where $M_j$ and ${\sigma }_j$ are, respectively, the atomic mass and the cross section of the atom labeled j. An important point is that the scattering cross section of Hydrogen is almost 2 order of magnitude that of the other atoms of the clay (Si, Al, O), so the $GDOS\left(\omega \right)$ gives mainly information on the vibrational modes of the H atoms. For the hydrated systems, the simulated spectra were calculated with Density Functional Theory Molecular Dynamics simulations since the confined water molecules remain liquids and, thus, disordered compared to monocrystalline clay network. After solving classical Newton equations, the $GDOS\left(\omega \right)$ was approximated by evaluating the velocity auto-correlation function of the atoms:

$$g_j\left(\omega \right) = \sqrt{\frac{1}{2\pi }}{\int }_{-\infty }^{+\infty }{e}^{-i\omega t}\sum _j\frac{〈{\overrightarrow{v}}_j\left(0\right)·{\overrightarrow{v}}_j\left(t\right)〉}{〈{\overrightarrow{v}}_j\left(0\right)·{\overrightarrow{v}}_j\left(0\right)〉}dt,$$

(2)

where ${\overrightarrow{v}}_j\left(t\right)$ is the instantaneous velocity of atom j. Before comparing the predictions obtained from the ab initio MD simulations with experimental data, the contributions from the dry clay network is subtracted [18]. The resulting spectra are displayed in Figure 3 for two hydration states corresponding, respectively, to the formation of one and two hydration layers [70]. In both longitudinal and perpendicular directions, the agreements between predicted $g_{H-water}\left(\omega \right)$ and experimental spectra are semi-quantitative with a fair matching of the characteristic frequencies (see Figure 3), partially validating the use of the above-mentioned approximation in the treatment of such local motions of the confined water molecules.

2.2. Quasi-Elastic Neutron Scattering

Quasi-Elastic Neutron Scattering (QENS) experiments are performed [12,28] to probe the diffusive motions of the confined water molecules and their orientational dependence regarding the clay surface. These measurements are performed with the same saponite synthetic clay at various hydration states and neutralized by sodium or calcium counterions. Since the dynamical range probed by these QENS measurement extend up to a few hundred pico-seconds (see Figure 1), the molecular modeling of the clay/water interface is performed by using the classical CLAYYFF force field [64], with a minor correction of its Lennard/Jones parameters introduced to improve the agreement between the structure of the clay/water interface obtained by X-Ray and neutron scattering experiments and simulated data [2,57,70]. Grand Canonical Monte Carlo (GCMC) simulations [71] are first performed to determine the number of confined water molecules as a function of the experimentally fixed water chemical potential [12,28] and the layer-to-layer distance previously determined by X-Ray Scattering (XRS) [70]. The mobility of the confined water molecules and neutralizing sodium or calcium counterions is then modeled by performing numerical simulations of Molecular Dynamics (MD) [72] (see Figure 4), starting from an equilibrium configuration previously determined by GCMC simulations.

The QENS spectra displayed in Figure 5 are directly evaluated from the Fourier Transform of the Intermediate Scattering Function (ISF):

$$F_s({\overrightarrow{Q}}_{\alpha },t) = \frac{1}{2N_w}\sum _{j=1}^{N_w}\sum _{nh=1}^2\exp \left[i{\overrightarrow{Q}}_{\alpha }·\left({\overrightarrow{r}}_{i,nh}\left(0\right) - {\overrightarrow{r}}_{i,nh}\left(t\right)\right)\right],$$

(3)

convoluted by the spectral resolution of INS experimental set-up [28]. As displayed in Figure 5, an excellent agreement is obtained in all cases between the experimental spectra and the simulated numerical data. As a consequence, one may be quite confident in the validity of the water mobility determined by these numerical simulations on a time-scale slightly shorter than nano-second. The water self-diffusion coefficients displayed in

Table 1. Water dynamical parameters, $D_{//}$ and $D_{\perp }$, respectively, the self-diffusion coefficient defined in the direction parallel and perpendicular to the clay surface, deduced from MD simulations for all samples investigated at 300 K. The data are from Table 3 in the Ref. [28] with permission from Ref. [28]. Copyright (2012) American Chemical Society.

Sample	${\mathit{D}}_{//}$ (10${}^{-10}$ m${}^2·$s${}^{-1}$ )	${\mathit{D}}_{\mathit{\perp }}$ (10${}^{-10}$ m${}^2·$s${}^{-1}$ )
Na+-saponite 0.8 e/uc-1W	2.4 ± 0.20
Na+-saponite 0.8 e/uc-2W	7.2 ± 0.50	3.0 ± 1.0
Na+-saponite 1.4 e/uc-1W	0.6 ± 0.05
Na+-saponite 1.4 e/uc-2W	3.1 ± 0.10	2.0 ± 0.5
Ca2+-saponite 0.8 e/uc-1W	1.6 ± 0.20
Ca2+-saponite 0.8 e/uc-2W	6.5 ± 0.30	3.0 ± 1.0
Ca2+-saponite 1.4 e/uc-1W	0.3 ± 0.05
Ca2+-saponite 1.4 e/uc-2W	2.4 ± 0.10	2.0 ± 0.5

are obtained by integrating the velocity auto-correlation function [73] of the water molecules,

$$D_{\alpha } = \underset{\tau \to +\infty }{lim}{\int }_0^{\tau }〈{\overrightarrow{v}}_{\alpha }\left(0\right)·{\overrightarrow{v}}_{\alpha }\left(\tau \right)〉d\tau ,$$

(4)

allowing to determine the local water mobility in the direction parallel and perpendicular to the clay surface [74]. In Figure 5, the wave numbers q probed by these QENS experiments vary between 0.3 and 1.8 Å${}^{-1}$, thus restricting the corresponding diffusion length between 21 and 3.5 Å. Furthermore, the same QENS measurements probe a limited dynamical range characterized by a spectral window of 0.25 meV with a resolution of 0.0025 meV (see Figure 5), corresponding to maximum diffusion times of 132 ps. By taking into account the order of magnitude of the apparent mobility of the confined water displayed in Table 1, this maximum diffusion period restricts the order of magnitude of the probed diffusion length to typically a few Angstroms. Furthermore, as displayed in Figure 4, the time-scale probed by the QENS measurements totally matches the dynamical range probed by MD numerical simulations. As a consequence, the agreement between these QENS experimental data and the numerical results obtained by MD simulations depends not only on the quality of the numerical model used in that study but also on the matching between the time-scales and diffusion-ranges probed by both experimental and numerical methods. As a consequence, the self-diffusion coefficients reported here in the direction perpendicular to the clay surface quantifies only the local mobility of the confined water molecules and does not correspond to macroscopic displacements totally impeached by the confining clay platelets (see Figure 6a). By contrast, the water mobility of the water molecules quantified by the same QENS measurements in the direction parallel to the clay surface can be used to determine the order of magnitude of the residence time of the same water molecules (labeled ${\tau }_C$ in Figure 6a) confined in the interlamellar space between two clay platelets. As displayed in Table 1, confinement under low electrostatic coupling (i.e., for 0.8 electron per unit cell) simply reduces the water mobility by a factor ten by comparison with bulk water. By contrast, at higher electrostatic coupling (i.e., for 1.4 electron per unit cell), the water mobility may be reduced by two orders of magnitude. These numerical results illustrate how geometric and energetic constraints can significantly modify the dynamics of confined water molecules.

2.3. Nuclear Magnetic Resonance Relaxometry

In the last decade, Nuclear Magnetic Resonance (NMR) relaxometry appeared as a powerful tool [75,76] to investigate dynamics of fluids under confinement by measuring the frequency variation of the NMR relaxation rates. Usual NMR spectrometers allow to measure easily the relaxation of the nuclear magnetization in the direction parallel and perpendicular to the static magnetic field, leading, respectively, to the longitudinal and transverse relaxation rates. For bulk liquids, these two independent measurements lead to the same relaxation rates. By contrast, the transverse relaxation rates of confined fluids largely exceed their longitudinal counterpart due to their large frequency variations. The purpose of NMR relaxometry is to probe that frequency variation of the relaxation rates in order to extract dynamical information on the molecular motions (rotation and diffusion) inducing the fluctuations of the magnetic coupling monitoring the NMR relaxation processes. For I = 1/2 spin nuclei [35,75,76], field-cycling NMR spectrometers are perfectly appropriate to perform such investigations since the order of magnitude of the probed relaxation rates is generally compatible with the time required to switch the magnetic field of these spectrometers (typically, milli-seconds). As an example, field cycling NMR relaxometry measurements were only performed with dilute dispersions of synthetic clays [37]. By contrast, the same field cycling spectrometers become inefficient for quadrupolar nuclei (spin I > 1/2) because of the large enhancement of their relaxation rates under confinement. In that context, spin-locking relaxation measurement [77,78,79,80] is a powerful tool since it may be performed with usual spectrometer and is simply limited by the short dead time required to generate a sequence of NMR pulses (typically, micro-seconds).

Furthermore, a complete understanding of the NMR relaxation processes requires a quantum analysis of the spin. For I = 1/2 spin nuclei, this is easily performed by using the four Pauli matrices [75] to describe the three independent states of the magnetization (M${}_x$, M${}_y$ and M${}_z$) completed by the unit operator, in order to obtain a complete basis set. For quadrupolar nuclei (spin I > 1/2), this basis set becomes incomplete and new operators must be added [81] (five for I = 1 spin nuclei like ²H, twelve for I = 3/2 spin nuclei like ⁷Li or ²³Na and up to sixty for I = 7/2 spin nuclei like ¹³³Cs). These supplementary operators describe new quantum states, also called coherences, accessible to the nuclear magnetization by using specific pulse sequences [82,83,84]. By measuring the corresponding multi-quanta relaxation rates [85], it becomes possible to separately identify the contributions from the different relaxation mechanisms of confined quadrupolar nuclei [10,41,43,44,51,80,86]. Finally, by performing selective multi-quanta spin-locking relaxation measurements [85], one can extract dynamical information on confined fluids from the frequency variations of the relaxation rates of these additional coherences.

In the case of heavy water, the electrostatic field gradient felt by the deuterium atom is oriented along the $\overrightarrow{OD}$ director. As a consequence, the NMR relaxation rate of these deuterium atoms may be determined by the correlation time describing the reorientation of that $\overrightarrow{OD}$ director by reference with the static magnetic field. In bulk water, that correlation time is about 3.6 pico-seconds, leading to a NMR relaxation rate of 2.4 s${}^{-1}$. By applying the same analysis to the NMR relaxation rate of confined water molecules, the measured relaxation rate (typically, 10${}^4$ s${}^{-1}$; see Figure 7) requires reorientation correlation time of the $\overrightarrow{OD}$ director with an order of magnitude of 15 nano-seconds, i.e., three order of magnitude larger than the value measured by QENS (see Michot et al. [12]). That interpretation is, thus, invalid since it neglects the influence of confinement on the specific orientation of the confined water molecules (see Figure 6d), leading to the splitting of the deuterium resonance line (see Figure 6b,c). As a consequence, the decorrelation of the electric field gradient felt by the confined water molecules is incomplete until they remain confined inside the interlamellar space of the same clay sediment (see Figure 6a). By contrast, after desorption from one oriented micro-domain, the residual quadrupolar coupling vanishes and the long-time decrease of the correlation function monitoring the water quadrupolar coupling is then driven by the water residence time ${\tau }_C$, as displayed in Figure 6a.

In that context, a first set of ²H multi-quanta spin-locking relaxometry experiments were performed to determine the average residence time of water molecules confined within the interlamellar space of clay platelets forming aggregates (Figure 6a) of dense sediments of Na⁺-beidellite clay [80]. As illustrated in Figure 6b, the apparent splitting ${\nu }_Q$ of the ²H resonance line varies as a function of the sediment orientation ${\beta }^{LF}$ into the static magnetic field, according to the expected relationship (Figure 6c):

$${\nu }_Q\left({\beta }^{LF}\right) = \frac{3{cos}^2{\beta }^{LF} - 1}{2}{\nu }_Q\left(0\right),$$

(5)

where ${\beta }^{LF}$ is the angle between the direction of the static magnetic field and the normal $\overrightarrow{n}$ to the surface of the clay film. This residual splitting of the ²H NMR resonance line results from the specific orientation of the confined water molecules as illustrated by an equilibrium configuration of the confined ions and water molecules obtained by Grand Canonical Monte Carlo simulations (Figure 6d).

By measuring the specific relaxation rates of the so-called $T_{10}\left(s\right)$, $T_{11}(a,s)$, $T_{20}$, and $T_{22}(a,s)$ coherences [82] of this I = 1 spin nucleus as a function of the sediment orientation into the static magnetic field (Figure 7), it becomes possible to quantify the relative contributions of the quadrupolar and heteronuclear dipolar NMR relaxation mechanisms [80]. Three specific angles (0°, 30°, and 90°) are then selected to determine the time evolution of four independent coherences (labeled $T_{11}\left(s\right)$, $T_{21}\left(a\right)$, $T_{21}\left(s\right)$, and $T_{22}\left(a\right)$ [82,83,84] under spin-locking conditions by applying an irradiation field of calibrated irradiation power ${\omega }_1$. For I = 1/2 spin nuclei, such spin-locking relaxometry measurements probe only a single angular velocity ($2{\omega }_1$) for each irradiation power [77,78,79]. By contrast, in the presence of a residual quadrupolar coupling (${\omega }_Q$), I = 1 spin nuclei probe three different angular velocities:

$$k_1 = \sqrt{{\omega }_Q^2 + 4{\omega }_1^2},k_2 = {\displaystyle \frac{{\omega }_Q + \sqrt{{\omega }_Q^2 + 4{\omega }_1^2}}{2}}\mathrm{and}{k}_3 = \frac{-{\omega }_Q + \sqrt{{\omega }_Q^2 + 4{\omega }_1^2}}{2},$$

(6)

for a single irradiation power ${\omega }_1$ [82].

Some raw data are displayed in Figure 8 in addition to their Fourier transforms illustrating the specific angular velocities ($k_1$, $k_2$, and $k_3$) probed by these spin-locking relaxometry measurements. As also illustrated by Equation (6), the three apparent velocities ($k_1$, $k_2$ and $k_3$) vary as a function of the sediment orientation which modulates the residual quadrupolar splitting (${\omega }_Q$); see Equation (5). As a consequence, it becomes possible to investigate a broad dynamical range by using a limited number of irradiation powers ${\omega }_1$ (see Figure 9). Figure 9 exhibits a transition between a plateau, at low angular velocities, and a power-law decrease, at high angular velocities [80]. By using the water self-diffusion coefficient of the confined water molecules previously determined by QENS measurements under equivalent conditions (see Section 2.2), numerical simulations of Brownian Dynamics are performed to determine the average residence time of water molecules confined within the interlamellar space between two beidellite platelets. That residence time fully corresponds [43] to the critical time ${\tau }_C$ evaluated from the measured angular velocity ${\omega }_C$ (${\tau }_C = 1/{\omega }_C$; see Figure 9). This result illustrates how multi-quanta spin-locking NMR relaxometry of quadrupolar nuclei largely extends up to $\sim {10}^{-4}$ s, the time-scale previously probed by QENS measurements.

The same spin-locking relaxometry measurements were performed by ⁷Li NMR spectroscopy to study the mobility of lithium counterions within aqueous dispersions of synthetic laponite clays [10]. Dense dispersions are prepared by uni-axial oedometric compression; more details are given in the Supplementary Materials. Above the Isotropic/Nematic phase transition [5], the ⁷Li NMR spectra exhibit a splitting of the resonance line (see Figure 10) that perfectly cancels at the so-called magic angle; see Equation (5). This behavior results from the occurrence of a single nematic phase, with a director parallel to the compression axis. Since ⁷Li is a I = 3/2 spin nucleus, six different angular velocities may now be probed by using ⁷Li NMR relaxometry with a single irradiation power ${\omega }_1$ [84] (see Figure 11):

$$\begin{array}{cc}\hfill {\lambda }_1 = \sqrt{{\omega }_Q^2 + 2{\omega }_1{\omega }_Q + 4{\omega }_1^2}& {\lambda }_2 = \sqrt{{\omega }_Q^2 - 2{\omega }_1{\omega }_Q + 4{\omega }_1^2}\hfill \\ \hfill {\lambda }_3 = {\omega }_1 + \sqrt{\frac{{\omega }_Q^2 + 4{\omega }_1^2 + {\lambda }_1{\lambda }_2}{2}}& {\lambda }_4 = {\omega }_1 + \sqrt{\frac{{\omega }_Q^2 + 4{\omega }_1^2 - {\lambda }_1{\lambda }_2}{2}}\hfill \\ \hfill {\lambda }_5 = {\omega }_1 - \sqrt{\frac{{\omega }_Q^2 + 4{\omega }_1^2 + {\lambda }_1{\lambda }_2}{2}}& {\lambda }_6 = {\omega }_1 - \sqrt{\frac{{\omega }_Q^2 + 4{\omega }_1^2 - {\lambda }_1{\lambda }_2}{2}}\hfill \end{array}.$$

(7)

As illustrated in Figure 12, a broad range of angular velocities is then investigated by using a limited number of angular velocities. Figure 12 also exhibits the dispersion curve of the NMR relaxation rates. Here, we also detect a transition between a low-frequency plateau and a high-frequency continuous decrease.

Numerical simulations of Brownian Dynamics (BD) are then required to interpret these results. In order to reproduce ionic condensation on the clay surface, a mean force potential [73] is evaluated from preliminary Monte Carlo simulations of the distribution (see Figure 13) of thousands Lithium counterions neutralizing the negative charge of an isolated laponite platelet. By using the weak overlap approximation [68], BD simulations are further performed for modeling lithium diffusion within a collection of partially oriented laponite platelets. As illustrated in Figure 14, the decorrelation of the memory function, induced by the fluctuations of the quadrupolar relaxation mechanism [41], evolves according to a $t^{-1}$ power law for diffusion time longer than ${10}^2$ ns. This behavior results from the strong electrostatic coupling between the highly charged clay platelet and the lithium counterions trapped in the so-called condensation layer (see Figure 14). This phenomenon is responsible for the low frequency logarithmic decrease of the ⁷Li NMR relaxation rates displayed in Figure 12. As a consequence, two different confinements may be detected by NMR relaxometry within such solid/liquid interfacial systems, induced either by geometrical constraints, as in the case of water molecules adsorbed within the interlamellar space of dense sediments, or by energetic coupling, as in the case of ionic condensation within dilute clay dispersions.

2.4. Pulsed-Gradient Spin-Echo Attenuation

Attenuation of the NMR signal measured after an echo pulse-sequence encoded by using gradients of the static magnetic field leads to dynamical information on the mobility of the NMR probes in the direction parallel to the applied field gradient [87]. These attenuation measurements directly probe the Intermediate Scattering Function; see Equation (1) in the Supplementary Materials for wave-vectors q defined by the relationship:

$$q = \frac{\gamma \delta G}{\pi },$$

(8)

where $\gamma$ is the gyromagnetic ratio of the NMR nucleus ($2.6752 \times {10}^8$ rad/s for ¹H and 1.037 $\times {10}^8$ rad/s for ⁷Li), and $\delta$ is the duration of the applied field gradient and G its strength; see Figure S2 in the Supplementary Materials for more details.

Because of the maximum strength of the applied field gradient (typically, $G_{max}$ = 1.6 T/m) and its duration (typically, $\delta$ = 500 $\mathsf{\mu }$s), the maximum resolution defined as resolution${}_{max} = 1/q_{max}$ is not better than 15 $\mathsf{\mu }$m, thus corresponding to macroscopic mobility measurements. Since the magnetic field gradient may be applied along any direction within the sample, ¹H PGSE NMR measurements were first performed to detect the Isotropic/Nematic phase transition of laponite aqueous dispersions [11,50,88] by measuring the anisotropy of the water self-diffusion tensor. These measurements were possible thanks to the great anisotropy of the laponite platelets (diameter∼300 Å, thickness∼10 Å). As illustrated in Figure 15, the Isotropic/Nematic transition of these laponite suspensions does not vary as a function of the neutralizing monovalent counterions and is well reproduced by numerical simulations of Brownian Dynamics mimicking water diffusion within a collection of partially oriented hard disks [50].

The same PGSE measurements were performed by ⁷Li NMR spectroscopy within aqueous dispersions of laponite platelets neutralized by lithium counterions. By contrast with the continuous decrease of the measured water mobility normalized by the bulk water mobility, Figure 16 exhibits a strong initial decrease of the lithium mobility by reference to its value measured within simple LiCl aqueous solutions. This phenomenon is the fingerprint of the ionic condensation of the lithium counterions at the surface of the highly charged laponite platelets. Numerical simulations of Brownian Dynamics may reproduce that phenomenon by including an effective mean force potential as described previously (see Figure 13). While these BD simulations agree perfectly with experimental data at low clay concentrations (see Figure 16), large differences appear at higher concentrations because of the limitation of the weak overlap approximation [68] implied in these BD simulations [50].

Because of the time scale required to build these magnetic field gradients (typically, 60 $\mathsf{\mu }$s), such PGSE measurements become inefficient to measure the mobility of water molecules or counterions confined within the interlamellar space of clay aggregates due to the corresponding enhancement of their transverse relaxation rates. A last illustration of the use of ¹H PGSE concerns the analysis of anisotropy of water diffusion in compacted clay samples made of micrometer-sized kaolinite particles [51,52]; see the Supplementary Materials for more details. Since kaolinite is unswelling clay, these measurements are appropriate to detect the macroscopic mobility of the water molecules diffusing around the clay particles, leading to the so-called tortuosity factor [89,90]. In addition, this mobility is also strongly impacted by the preferential orientation of lamellar particles in the system (see Figures 18 and 19) [52]. These properties can be quantified based on X-Ray Scattering technique on bi-dimensional detector plate (2D-XRS) (see Figure 17) and the analysis of intensity distribution along the diffraction ring of the 001 reflection [52]. This measurement allows extracting the order parameter $〈P_2〉$ through the second order Legendre polynomial function [91], a parameter ranging from 0 for isotropic organization of particle to 1 if particles are perfectly oriented. The influence of preferred orientation for a constant porosity value can reveal a reduction of water mobility by a factor 2 when increasing anisotropy in particle orientation (see Figure 18). Such anisotropy in water mobility is well highlighted using ¹H PGSE-NMR technique and can be quantitatively reproduced by DB simulations using virtual porous media mimicking the size, shape, and orientation of the particles in the porous media (see Figure 18). By accounting for the self-diffusion of water probes in the interparticle volume and elastic collisions with particle interfaces, the simulated data are found to well reproduce experimental measurements (see Figure 18 and Figure 19), thus allowing to decouple the contribution of porosity and particle orientation on the overall water diffusion process.

2.5. Magnetic Resonance Imaging

Magnetic Resonance Imaging (MRI) experiments were also performed by using ¹H NMR spectroscopy to obtain long time-scale information on the macroscopic mobility of water molecules within clay porous media made of kaolinite particles [51]; see the Supplementary Materials for more details. The experimental conditions are selected to probe a total length of 25 mm. The porosity of the dry kaolinite sample (${\varphi }_{dry}$ = 0.52) is evaluated by helium pycnometry [51]. The clay sample is initially saturated with bulk water, leading to a total length (noted L${}_1$) of 12.3 mm (see Figure 20a). The MRI experiment starts just after the addition of a limited amount of heavy water (length L${}_2$ = 7 mm) on the top of the hydrated clay sediment, leading to a total length (L = 19.3 mm) fitting the experimental window. As illustrated in Figure 20a, the ¹H water concentration profiles exhibit some attenuation of the ¹H NMR signal at the bottom of the clay sample because a significant section of the clay sample (5.3 mm) is located outside the optimal detection area of the detection coil. We also detect alterations of the water concentration profile near the clay/heavy water interfacial area, further limiting the useful length of the hydrated sample to 5 mm. By integrating the ¹H NMR signal over the corresponding limited area (see Figure 21), we obtain information of the long time-scale, $\tau \approx (40 \pm 2) \times {10}^3$ s, implied in this exchange between bulk and heavy water. Numerical simulation are performed [51] to fully validate that analysis (see Figure 20b and Figure 21) by using the porosity of the dry clay sample and the water mobility previously measured by ¹H PGSE attenuation measurements (see Section 2.4 above). That result clearly illustrates how MRI experiments may be used as complementary tool to probe the mobility of confined molecules over a broad macroscopic time-scale (from minutes up to day).

2.6. Macroscopic Diffusion Experiments

Out-Diffusion experiments were performed to determine the time-scale characterizing the exchange between Ca²⁺, initially located in the interlayer of vermiculite disks and neutralizing its charge, and Na⁺ or Sr²⁺, present in aqueous solutions [53,54] (see Figure 22 and the Supplementary Materials for more details). Firstly, the resulting exchange rates between these counterions were well predicted by BD simulations using self-diffusion coefficients for interlayer and aqueous bulk cations, calculated from preliminary MD simulations and the geometry of both aqueous reservoir and vermiculite disks (see Figure 22a) in the case of the Ca²⁺/Na⁺ exchange [53]. The results of this study showed that the time-evolution of the concentration gradient of interlayer Ca²⁺ was the main driving force interpreting exchange rates measured at all aqueous reservoir salinity investigated, while limitation at the solid/liquid interface, due to diffusion in the aqueous reservoir, can acts as a limiting factor to interpret exchange rates measured at low salinity conditions. As displayed in Figure 22b, the same experimental data were also well reproduced by a classical finite-volume model based on continuous differential equations, constrained by a self-diffusion coefficient of the slowest interlayer cation. That approach uses also selectivity coefficients obtained independently in batch and describing the ionic exchange at thermodynamic equilibrium (see Ref. [54] for more details). Using the two numerical approaches (BD and classical finite-volume model based on continuous differential equations), data were always based on self-diffusion coefficients obtained from preliminary numerical simulations of Molecular Dynamics (MD) without the use of any fitted parameter [53,54]. As displayed in Figure 22, the agreement between experimental and numerical data validates the use of self-diffusion coefficients obtained from MD simulations to predict the dynamic of the exchange between interlayer cations and aqueous ones over several days (up to 20 days) and then the multi-scale approaches proposed.

As a second illustration of long-time diffusion process, Through-Diffusion experiments are commonly used to measure the flux of a tracer, traversing a compacted sample and to extract an effective diffusion coefficient. Such effective diffusion coefficient is directly linked to the pore diffusion parameter commonly derived from NMR measurements probing only the porosity accessible for the investigated probe. A large number of experiments were performed in literature to investigate the influence of porosity on the diffusion of water and ions in porous media made of clayey particles [69,92,93,94]. Water diffusion coefficients obtained using Through-Diffusion experiments and PGSE-NMR analyses performed with compacted kaolinite samples, characterized by the same porosity and preferential particle orientation, revealed similar results (see Figure 23), despite the different time scale probed by the two methods (i.e., ms for PGSE-NMR versus several days for Through-Diffusion [52]). This finding provides interesting information regarding representative time and length scales for the investigation of water dynamics in clay porous media. Moreover, this consistency allows taking advantage of both techniques for the study of clay/water interfaces. On the one hand PGSE-NMR cannot be used when sample contains paramagnetic elements, whereas Through-Diffusion is insensitive to the chemical composition of the sample. On the other hand, PGSE-NMR provides the full tensor of water dynamics whereas Through-Diffusion experiments requires to be repeated when changing the direction of diffusion probed. For the analysis of clay porous media with both interparticle and interlayer porosity (dual porosity media), Through-Diffusion was successfully applied on vermiculite to investigate the role played by the different types of porosities and preferential orientation of particle on water dynamics [55,94]. Among the large set of different numerical methods to interpret the experimental data [69,95,96], Asaad et al. [55] used BD simulations on representative virtual porous media to mimic the diffusion of water probes in both interparticle and interlayer volume, using self-diffusion coefficients extracted from MD simulations. The good agreement between experiments and simulations (see Figure 23) allowed highlighting the marginal role of interlayer volume on the flux of water tracer traversing samples, and the similar role of preferential orientation of particle when dealing with single or dual porosity media.

3. Conclusions

This review focuses on dynamical properties of the water molecules and neutralizing counterions confined within clay particles. By exploiting a large panel of complementary experimental studies, it becomes possible to investigate a large time scale varying between pico-second and days, implying local motions, as well as macroscopic displacements, within the clay aggregates. In addition to these experimental studies, we also performed multi-scale modeling of these clay/water interfacial systems in order to better extract dynamical information from the raw experimental data and further validate the conclusions of our analysis. This allowed deriving a complete bottom-up analysis of the dynamical properties of ions and water molecules confined within dense clay sediments. In that context, the interpretation of such complex dynamical phenomena is based on numerical modeling exploiting realistic models and dynamical information (residence time, local mobility) extracted from preliminary measurements performed on a smaller scale. Such complementary approaches should be useful for the analysis of other solid/liquid interfacial systems, such as cementeous and zeolitic materials, porous silicates, membranes and micelles, and aqueous dispersions of macromolecules and biopolymers.