ARTICLE pubs.acs.org/est
Water Dynamics in Hectorite Clays: Infuence of Temperature Studied by Coupling Neutron Spin Echo and Molecular Dynamics Virginie Marry,†,‡ Emmanuelle Dubois,*,†,‡ Natalie Malikova,§ Serge Durand-Vidal,†,‡ Stephane Longeville,§ and Josef Breu|| †
UPMC Univ Paris 06, UMR 7195, PECSA, F-75005 Paris, France CNRS, UMR 7195, PECSA, F-75005 Paris, France § LLB (CEA-CNRS), CEA Saclay, 91191 Gif-sur-Yvette, France Lehrstuhl f€ur Anorganische Chemie, Universit€at Bayreuth, Universit€atsstr. 30, 95440 Bayreuth, Germany
)
‡
bS Supporting Information ABSTRACT: Within the wider context of water behavior in soils, and with a particular emphasis on clays surrounding underground radioactive waste packages, we present here the translational dynamics of water in clays in low hydrated states as studied by coupling molecular dynamics (MD) simulations and quasielastic neutron scattering experiments by neutron spin echo (NSE). A natural montmorillonite clay of interest is modeled by a synthetic clay which allows us to understand the determining parameters from MD simulations by comparison with the experimental values. We focus on temperatures between 300 and 350 K, i.e., the range relevant to the highlighted application. The activation energy Ea experimentally determined is 6.6 kJ/mol higher than that for bulk water. Simulations are in good agreement with experiments for the relevant set of conditions, and they give more insight into the origin of the observed dynamics.
’ INTRODUCTION Investigating water dynamics in anisotropic and confined media such as clays is an essential key for the understanding of transport properties in such materials, but, more generally, in all divided solids, used for example in catalysis, or porous media, which are of interest for example as soil models in environmental sciences.1 For clay materials, the role of water is especially important in the frame of the storage of radioactive wastes, the clay being around the waste and surrounded by water and ions that can diffuse inside the clay, phenomena which need to be predicted on large time scales. This implies the knowledge of the behavior of water on different scales, among which the nanoscale is the smallest one, and also for a range of temperatures between 300 and 350 K, a range which has been seldom considered.2-4 In such systems, the quasielastic neutron scattering techniques (QENS) are particularly appropriate to study water dynamics on the nanoscale (several angstroms and several picoseconds), scales which can also be probed by molecular dynamics simulations. Although many studies were performed with QENS (see in ref 2, which is a review on the topic), very few studies associated simulations and experimental measurements.5,6 This coupled approach is nevertheless very rich as it can improve the interpretation of experiments and the understanding of the determining parameters in the clay materials. We use here this coupled approach on a synthetic hectorite clay that we previously compared at room temperature to the natural montmorillonite reference clay (MX80) used in radioactive waste storage. From these results, we concluded that the two clays behave similarly. Also, the simulated water dynamics gave a reasonable qualitative agreement with the experimental water dynamics, however overestimated by a factor of 2-3.6 r 2011 American Chemical Society
Therefore, we propose here to focus on the thorough comparison of experiments and simulations. For this purpose, the synthetic hectorite is especially appropriate thanks to its characteristics similar to the natural montmorillonite of interest and because it is a very well-defined material with no interstratification (monolayer and bilayer states are obtained pure).7 This highly improves the ability to compare the measured quantities with the ones calculated from simulation, which can now be performed with new types of clay force fields giving dynamics slower than previous calculations.8 Water dynamics is determined here for several temperatures on the hectorite clay hydrated with one water layer. We use the technique of neutron spin echo (NSE), seldom applied to such clays and more time-consuming than the usually used neutron time of flight (TOF) technique. We focus on spatial scales for which the measured signal can be interpreted as due to water translation only (rotation can be neglected); therefore, NSE is more appropriate than TOF for this purpose. Only one parameter is extracted: the translational diffusion coefficient of water. It is obtained from the whole set of data thus giving a reliable determination, the data being fitted using a 2D model adapted to the case of an unoriented powder as the one we study here.9 In parallel, these diffusion coefficients can be extracted from molecular dynamics simulations. From the analysis of the water trajectories on long time scales, the translational dynamics is isolated and can be compared to the 2D experimental diffusion. Received: September 19, 2010 Accepted: February 15, 2011 Revised: February 7, 2011 Published: March 07, 2011 2850
dx.doi.org/10.1021/es1031932 | Environ. Sci. Technol. 2011, 45, 2850–2855
Environmental Science & Technology Simulation allows us to quantify the changes in water diffusion induced by parameters such as the relative horizontal shift of adjacent clay layers. It is impossible to control such parameters systematically in the experimental clay system. Our simulations are performed with two currently widely used clay force fields, in order to test their influence on the physical properties obtained.
’ EXPERIMENTAL SECTION Hectorite Clays. The studied clay is a synthetic hectorite with a half unit cell composed of [Na0.4]inter[Mg2.6Li0.4]oct[Si4]tetO10F2. This composition corresponds to a charge per half unit cell of 0.4 e, chosen to be similar to the charge of the natural montmorillonite MX80, which is the reference material for the studies of radioactive waste storage in France.10 The negative charge on clay layers is compensated by Naþ cations in the interlayer. The material is fluorated, meaning that the structural OH- groups in the clay are replaced by F- atoms. This hectorite clay is synthesized from a melt at high temperature, which produces a homogeneous material in terms of composition and charge density especially.7 Individual clay platelets have typical dimensions around 12 μm and a thickness of 6.3 Å. The water adsorption and desorption isotherms show two clear plateaux corresponding to pure monolayer and pure bilayer of water, which implies no water adsorption in micropores. We focus here on the monolayer state, obtained by hydration of the dry clay with RH = 43% (relative humidity), which corresponds to an uptake of 4 water molecules per Naþ cation from the dry state and to an interlayer distance of 12.4 Å.11 In the dry state, the water amount is less than 1 molecule per cation according to ref 12. Also, from an XRD study on the present hectorite following the method used in ref 13, the best fit of the positions and intensities of the (00l) peaks of the complete XRD pattern is obtained with 4H2O/ Naþ.14 Therefore, the total number of H2O per cation lies between 4 and 5. Neutron Spin Echo. The water dynamics is determined by measurement of the incoherent scattering of hydrogen atoms with the neutron spin echo (NSE) technique, which gives access to the intermediate scattering function. The sample is placed in a flat aluminum cell of thickness 5 mm in order to obtain sufficient signal (equivalent water thickness 0.4 mm). The hectorite powder inside the aluminum cell is not oriented, which has been checked by neutron diffraction measurements on similar samples: the intensity of the main diffraction peak is independent of the orientation of the cell with respect to the direction of the neutron beam (G4.1., LLB, CEA Saclay). Neutron spin echo (NSE)15 is performed on MUSES spectrometer (LLB, Saclay, France).16 By combination of both the NSE and NRSE (neutron resonance spin echo) technique,17 measurements for correlation times up to approximately 1000 ps are achieved. Measurements have been performed at 300, 323, and 347 K, and it has been checked that the monolayer state is kept whatever the temperature. The wavevectors Q are chosen in the range 0.3 < Q < 1 Å -1, far enough from the Bragg peaks so that the scattered intensity is mainly incoherent. Therefore, NSE gives here the incoherent intermediate scattering function, Sinc(Q, t). Simulations. Interaction Potentials. Molecular simulations imply the evaluation of the interaction potential between all the atoms of the system, which depends on the force field chosen to describe the system. Thus, the influence of the force field is the first point that we will examine (see details in Supporting
ARTICLE
Information). We chose here the SPC/E model (simple point charge/extended) for water molecules,18 known to describe quite well the dynamics of bulk water. SPC/E adapted ion parameters have been employed.19 Several force fields have been proposed to describe clay materials,20-27 and some of them have been already adjusted or tested with SPC/E water for specific clays. They thus need to be adapted to our synthetic clay. In order to examine the influence of the potential on the calculated dynamics, we used two different force fields: the SSFF (Skipper and Smith force field), adapted to SPC/E water by Smith,24,28 which we have extensively used on montmorillonite29-32 and also compared with experimental data,6,33 and the clayFF force field, created by Cygan et al.26 and recently used with SPC/E by several authors.8,34-36 These two potentials give structural results and swelling behaviors in clays in fair agreement with experiments.26,37-39 However, the former potential seems to overestimate the diffusion coefficient of water in montmorillonite,6 and the latter gives diffusion coefficients closer to experiments.8 In the case of SSFF force field, we applied the parameters already used to simulate natural hectorites.40,41 In the case of clayFF, we took the same parameters as those used for talc,42 which is the uncharged version of hectorite, and modified them to take into account the substitution of Mg2þ by Liþ and the presence of fluorine (see Supporting Information where the parameters of the potentials used are given in Table S1 for both force fields). Simulation Details. The simulation boxes contain 2 clay layers in the xy plane, with water and ions in the interlayer. The crystallographic positions of atoms in the clay unit cell are well-defined (taken from ref 43). However, several parameters are not perfectly determined in the real clay: the number N of H2O units per cation, the interlayer spacing, and the preferential positions of adjacent clay surfaces. Therefore, we explore their influence on dynamic properties using simulations. In the present study, N is fixed either to 4H2O/Naþ or to 5H2O/ Naþ (with clayFF only), which corresponds to the range of possible experimental values. The pressure is fixed to 1 bar. For the chosen water content, for a fixed temperature and pressure perpendicular to the clay layer, allowing for vertical and horizontal translation of clay layers, the equilibrium interlayer spacing and clay surface arrangements can be calculated. Simulations at T = 298 K with N = 4 led to equilibrium interlayer spacings of L = 12.2 and 12.05 Å for clayFF and SSFF, respectively. It is lower than our experimental value 12.4 Å,11 especially with SSFF. We chose to perform molecular dynamics simulations with L = 12.2 Å with both force fields in order to allow a more direct comparison between them. Simulations with L = 12.3 Å were performed in order to approach the experimental interlayer spacing. The simulations with N = 5 gave L = 12.5 Å with clayFF, which is a little higher than the experimental value. As the temperature has been seen to have a very small influence on the value of interlayer spacings, especially in the monohydrated case,30 we kept the above spacings at all the simulated temperatures: 298, 323, and 347 K. The last parameter explored is the preferential position of adjacent clay surfaces. It is found from Monte Carlo that the preferred position is configuration B (Figure 1 and Figure S1 in Supporting Information) where a silicon atom of one of the surfaces is facing the center of the hexagonal cavity of the other surface. Most simulations were performed with configuration B; however, configuration D was probed. Figure 1 summarizes the combinations of the 4 parameters studied in molecular dynamics simulations for clays. For comparison, bulk water simulations were also performed. 2851
dx.doi.org/10.1021/es1031932 |Environ. Sci. Technol. 2011, 45, 2850–2855
Environmental Science & Technology
ARTICLE
Figure 1. Conditions of simulations: the columns correspond to the force field used, the number N of H2O per cation Naþ, the interlayer spacing L, and CLP, which is the clay arrangement used. The positions of the clay layers in B and D configurations are plotted on the right: hexagons are the hexagonal cavities, the vertices of which are the silicon atoms of the clay surface.
Neutron Spectrometry. The intermediate scattering functions S(Q, t) measured by NSE are plotted in Figure 2 for Q = 0.4 Å-1 for the three temperatures studied (see the other data in Figures S3, S4, and S5 in Supporting Information). As this incoherent intermediate scattering function corresponds to the self-correlation function of the positions of hydrogen atoms in the sample, the vibrations of the atoms as well as the rotational and translational dynamics contribute to S(Q, t). However, given the probed time scale, the vibrations are too fast and can only introduce a coefficient A lower than 1 in front of S(Q, t) (see eq 1). We focus here on Q values smaller than 1 Å-1 in order to isolate the translational diffusion of water molecules in the interlayer space, which can be obtained accurately in this Q range.9 Rotational diffusion can be extracted for larger Q0 s; however, there is no model for rotation in clays, and the uncertainty of this determination is very large.44 To decrease the number of fitting parameters, some constraints can be introduced, for example in the form of an oriented sample with individual measurements done in different directions with respect to clay layers.5 Concerning the translation, as already underlined, the TOF technique is often used and preferred to NSE in order to determine the water dynamics in clays, among which is the translational diffusion. We compared the two techniques on the studied monohydrated hectorite. TOF was performed at several resolutions (fwhm of 28 and 15 μeV); however, the broadening of S(Q, ω) was very small even using the best resolution.11 S(Q, t) can then be modeled considering or not an elastic contribution, which can correspond to immobile H atoms or H2O molecules too slow for the time window of the technique. As we cannot discriminate between these hypotheses by the quality
)
)
’ RESULTS AND DISCUSSION
where θ is the angle between Q and the direction perpendicular to the clay layer. This model is the best suited here as the diffusion is strongly anisotropic in such a system, which is nevertheless randomly oriented according to our measurements.45 Note that the use of the standard 3D model to interpret the data results in an underestimation of D by 25%.11 For each temperature, the whole set of S(Q, t) is fitted in order to determine the value of D , reported in the first column of Table 1. The activation energy Ea can be obtained by a linear fit of the logarithm of D as a function of 1/T according to the following: )
by taking the slope of the mean-squared displacement (MSD) in the xy plane as a function of time between 40 and 1000 ps, which corresponds to the diffusive regime of water molecules8 (an example of MSD versus time is given in Figure S2 of Supporting Information).
lnðD Þ ¼
-Ea þb kB T
In Figure 3, the plot of ln(D /Dbulk) directly gives ΔEa = Ea(clay) - Ea(bulk). Here, the experimental values are Ea(clay) = 23.3 kJ/mol and Ea(bulk) = 16.7 kJ/mol; therefore, ΔEa = 6.6 kJ/mol. Only a few studies in the literature report on measurements on the water dynamics at room T or higher T in clays with a monolayer of water. Our values at room T can nevertheless be compared with the results for a synthetic saponite in the monolayer state.5 Their values (10-10 to 4 10-10 m2 s-1), obtained for 2 different monolayer states and different types of measurements in oriented samples, are on average smaller than ours, which may be due to the higher charge in their saponite (0.7e per half unit cell although ours has 0.4e). Also, our )
t sf¥
hx 2 þ y2 i 4t
of the fit, we obtain a range for D, the extreme values differing by a factor of 3 (see Figure S6 of the Supporting Information). Therefore, we need the NSE, which gives access to longer times, and shows a nonzero level in S(Q, t) for large t if immobile H atoms exist or atoms slower than the instrumental resolution, the height of the plateau corresponding to the proportion of these slow species. Such a plateau does not appear in our data meaning that all H2O species are mobile on the timescale probed. The NSE curves are fitted by the powder averaged 2D model of diffusion, in which the amplitude A, which takes into account short time motions, and the diffusion coefficient D in the plane parallel to the clay layers (D^ = 0 perpendicular to the clay layers) are the fitting parameters Z π 1 2 2 ð1Þ e-D Q ½sin ðθÞt sin θ dθ ÆS2d ðQ , tÞæθ ¼ A 2 0
)
)
D ¼ lim
Figure 2. Experimental S(Q, t) for Q = 0.4 Å-1 and the fits according to a powder average model (see text for details).
)
After a phase of equilibration, the dynamic quantities were calculated based on 5 ns molecular dynamics trajectories. Horizontal diffusion coefficients were calculated from the Einstein relation
2852
dx.doi.org/10.1021/es1031932 |Environ. Sci. Technol. 2011, 45, 2850–2855
Environmental Science & Technology
ARTICLE )
)
)
Table 1. D (109 m2 s-1) for the Clay from NSE Experiments, D Calculated at Several T Values and Simulation Conditions (Given in Figure 1), and D for Bulk Water from Simulations and from Experiments Dexp clay
Dexp bulka
298
0.28 ( 0.03
2.4 ( 0.025
323 347
0.53 ( 0.05 3.95 ( 0.04 1.0 ( 0.1 5.97 ( 0.06
T (K)
ΔEa (kJ/mol):
6.7 ( 2.0
III
IV
0.41 ( 0.03 0.49 ( 0.03 0.62 ( 0.03 0.71 ( 0.04 0.85 ( 0.04 1.10 ( 0.05 1.14 ( 0.06 3.0 ( 1.9
0
II
0.99 ( 0.06 1.33 ( 0.06
0.5 ( 1.8
-0.9 ( 1.6
V
0.271 ( 0.020 0.36 ( 0.02 0.51 ( 0.04 0.84 ( 0.04
0.60 ( 0.03 0.89 ( 0.04
5.5 ( 1.9
1.7 ( 1.6
VI
VII
Dsim bulk
0.45 ( 0.03 0.55 ( 0.03 2.60 ( 0.06 0.71 ( 0.04 0.86 ( 0.04 4.19 ( 0.06 1.04 ( 0.05 1.25 ( 0.05 5.88 ( 0.08 0.5 ( 1.5
0.1 ( 1.5
0
From ref 46.
obtained in bulk conditions is specific for the chosen water model; thus, it will also be present in clay. Consequently, we compare the experimental and simulated D /Dbulk and their corresponding values for ΔEa = Ea(clay) - Ea(bulk), which are more representative of diffusion processes characteristic of confinement in clays. The results are given in Table 1 and plotted in Figure 3. Concerning the diffusion coefficients D, the first conclusion is that their values for water are lower in the clay than in the bulk. Second, if D values were expected to depend on the force field, they appear also very sensitive to the conditions of simulations for the same force field. The interlayer spacing has a large influence: when increased by 0.1 Å (only parameter modified between simulations II and III or V and VI), it induces an increase by 15-35% depending on T and the force field (see left graphs in Figure 3). Also, the arrangement of clay layers is decisive: at 298 K, water diffuses more slowly by 15-25% when the clay layers are in position D than in position B. Finally, in the same conditions of simulation, clayFF always gives a diffusion coefficient lower than SSFF. From the variation of D versus T, the Ea values deduced (right graphs in Figure 3) are in all but one case higher than in the bulk. As D, they depend on the conditions of the simulations: for the same configuration; clayFF always gives Ea slightly higher than that given by SSFF. Moreover, Ea is strongly influenced by the clay layers' arrangement: ΔEa for position D is about 3 kJ/mol higher than ΔEa for position B. Finally, an increase of the interlayer spacing by 0.1 Å decreases Ea. In the case of SSFF in particular, ΔEa is slightly negative when the interlayer spacing is 12.3 Å. This latter result is in agreement with our previous simulations on a montmorillonite with the same force field, where the monohydrated state was defined by 6 water molecules per cation and an interlayer spacing equal to about 12.4 Å:30 the approximate activation energy was found to be in the range 1215 kJ/mol, thus, close to the bulk SPC/E value. The analysis of the influence of the conditions on the calculated dynamics allows us to propose explanations of the mechanism of water diffusion in clays and the reasons for a higher activation energy than in the bulk. Let us look first at the mechanism of the bulk water motion before discussing the water dynamics in the clay. Already in the bulk, the motion of water is complex with continual formation and breaking of H-bonds between adjacent H2O molecules. A recent model of Laage et al.,47 on the basis of simulations of bulk SPC/E water, shows that a water molecule can break an H-bond with an old partner when another partner approaches via a fast jump mechanism leading to its reorientation. The limiting step is shown to be the concomitant translational motion of the departing partner and of the incoming one. Moreover, the dependence of the characteristic reorientation time determined by both NMR and simulations with the temperature gives Ea = 14.6 kJ/mol, which is close )
a
I
)
)
)
Figure 3. D /Dbulk versus T and ln(D /Dbulk) = -(ΔEa)/(RT) þ cst versus 1000/T for experiments and simulations: blue crosses, NSE; red symbols, D arrangement; black symbols, B arrangement; circles, SSFF; triangles, clayFF; empty symbols, L = 12.2 Å; filled symbols, L = 12.3 Å. For clarity, the plots giving ln(D /Dbulk) were vertically shifted with respect to each other. The linear regression was done on the same temperature range as the experiment.
)
measurements can be compared with the results obtained on a natural montmorillonite with a similar amount of water:3 at 300 K, D = 5.0 ( 0.6 10-10 m2 s-1, and at T = 318 K, D = 7.4 ( 0.6 10-10 m2 s-1, which leads to Ea = 17 kJ/mol. These D values are on the contrary larger than ours, and Ea is close to bulk water and thus smaller than our value. It could be due to two reasons as we observed previously in our measurements on natural montmorillonite:6 first, a mixture of monolayers and bilayers can exist in the natural material; second, as shown above, even the highest TOF resolution (fwhm = 13 μeV) cannot see the slowest water molecules. Thus, the characteristic times are underestimated, especially at low T, which decreases Ea. However, our Ea values, higher than bulk water, are in agreement with the macroscopic measurements in ref 3 using tracers, which give values around 20 kJ/mol for montmorillonite with 2 water layers. Our results will now be compared with the same quantities calculated from MD simulations. Simulations. The calculated diffusion coefficients D are collected in Table 1 for all the conditions of simulation (see Figure 1), together with diffusion coefficients of bulk water at the same temperatures. For bulk water, although the SPC/E model gives diffusion coefficients in reasonable agreement with experiment, their values can differ from experiment by about 10% depending on the temperature. Therefore, considering only the 3 temperatures 298, 323, and 347 K measured by NSE, Easim = 14.8 kJ/mol, which is slightly lower than Eaexp = 16.7 kJ/mol determined by ref 46 over the same T range. This difference
2853
dx.doi.org/10.1021/es1031932 |Environ. Sci. Technol. 2011, 45, 2850–2855
Environmental Science & Technology
)
)
influence of parameters that cannot be tuned experimentally. The confinement between the layers and the strength of interactions with clay surface can explain that D < Dbulk and that Ea > Ebulk a . It shows that when compared to a reliable experiment, simulations give information about the diffusion processes involved in clays. The next step will be to improve the force fields which still give L = 12.2 Å (Monte Carlo with clayFF) lower than the experimental value L = 12.4 Å. Therefore, the influence of the very small variation of L with temperature and the role of cation specificity should be taken into account in future studies. The present study and its future extensions, for example to the case of the bilayer of water in the clay, open the possibility to use simulations to study detailed processes occurring in clays which cannot be seen experimentally, to analyze data mixing several dynamics (rotation and librations), and to model accurately more and more complex systems approaching the real natural clays. )
to the experimental and simulated Ea obtained from the diffusion coefficient in the bulk. From these results, we can make the hypothesis that the translational and reorientation processes are closely related. In the clay, we observed previously32 that water can form H-bonds with oxygen atoms of clay surfaces and that their formation and breaking have different kinetics from the bulk. Laage et al. showed that the slowdown in H-bond exchange dynamics at the vicinity of a hydrophilic H-bond acceptor surface can be due to two main reasons: an excluded volume effect and a stronger H-bond formed with the acceptor.48 The first effect is due to the presence of the acceptor surface which hinders the approach of a new water H-bond acceptor. It is an entropic effect which hardly affects the energy barrier for the jump reorientation of the molecule but is rather responsible for a reduction of the number of accessible transition state configurations, therefore leading to a reorientational slowdown. However, the second effect has a strong enthalpic contribution because the cost to elongate an H-bond depends on its strength and thus has a strong influence on activation energies. A similar mechanism could be considered in clays. More generally, the force field created by the clay surfaces, which is responsible for attractive or repulsive interactions of mobile species with clay atoms, is expected to play a significant role in this enthalpic contribution, thus influencing strongly both diffusion coefficients and activation energies. In the case of SSFF, the charge of surface oxygens is only slightly less negative than SPC/ E oxygen charge, thus, a priori slightly less attractive toward water hydrogens. As a result, the slowdown of D compared to Dbulk is expected to come mainly from (i) the geometrical excluded volume effect due to the presence of surfaces, and (ii) the cations, as water is involved in their hydration shell. Thus, Ea = or e Ebulk a is expected. This is actually the case with L = 12.3 Å (simulation III). On the contrary, for clayFF, the charge of the surface oxygen is more negative than the bulk by -0.202e, leading to stronger H-bonds with the surface. This explains the Ea values always being larger for clayFF than for SSFF and the bulk. However, it with SSFF for the cannot explain Ea values larger than Ebulk a lowest L = 12.2 Å (simulations I and II). Here it can be put forward the proximity of a second surface, which can strongly modify the interaction potential created by the clay (thus the H-bond strength) compared with the case of only one surface. It is proven by the strong influence of clay arrangements on the xy distributions of oxygen and hydrogen atoms between the clay layers (see Figure S7 in Supporting Information). This effect, which increases activation energies, is shown to rapidly decrease when L increases (simulations II/III and V/VI). In conclusion, when directly compared with experiment, the best simulated Ea is obtained in simulation IV, corresponding to L = 12.2 Å, CLP D with clayFF, but simulation V, with CLP B, gives diffusion coefficients closer to experiment even if Ea is smaller than the experimental value. However, considering the strong influence of L on Ea, taking into account even very small variations of L with temperature would lead to an increase of simulated Ea. Therefore, the Ea values in Table 1 have to be considered as lower bounds. This study has shown that D for water is very sensitive to the detailed parameters of MD simulations, in particular to the chosen force field, although structural properties such as the clay swelling behavior are much less sensitive and lead always to similar values of interlayer spacings L for the same thermodynamic conditions.26,37 Also, L and the arrangement of clay layers highly influence D and Ea. Therefore, simulations enlighten the
ARTICLE
’ ASSOCIATED CONTENT
bS
Supporting Information. Additional details on simulations and experimental data not shown in the text. This material is available free of charge via the Internet at http://pubs.acs.org/.
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
’ ACKNOWLEDGMENT The authors thank J. M. Zanotti for his help during TOF (Mibemol, LLB, CEA Saclay), G. Andre for his help during the diffraction experiments (G4.1, LLB, CEA Saclay), A. Botan for the realization of Figure S7 in Supporting Information, and A. G. Kalinichev for fruitful discussions. ’ REFERENCES (1) Handbook of Clay Science; Bergaya, F., Theng, B. K. G., Lagaly, G., Eds.; Elsevier: New York, 2006. (2) Bihannic, I.; Delville, A.; Deme, B.; Plazanet, M.; Villieras, F.; Michot, L. Clay swelling: new insights from neutron-based techniques. In Neutron Applications in Earth, Energy and Environmental Sciences; Liang, L., Rinaldi, R., Schober, H., Eds.; Springer: New York, 2009; Chapter 18, p 521. (3) Sanchez, F. G.; Gimmi, T.; Juranyi, F.; van Loon, L.; Diamond, L. Linking the diffusion of water in compacted clays a two different time scales: tracer through-diffusion and quasielastic neutron scattering. Environ. Sci. Technol. 2009, 43, 3487–3493. (4) Swenson, J.; Bergman, R.; Longeville, S. A neutron spin-echo study of confined water. J. Chem. Phys. 2001, 115, 11299–11305. (5) Michot, L.; Delville, A.; Humbert, B.; Plazanet, M.; Levitz, P. Diffusion of water in a synthetic clay with tetrahedral charges by combined neutron time-of-flight measurements and molecular dynamics simulations. J. Phys. Chem. C 2007, 111, 9818–9831. (6) Malikova, N.; Cadene, A.; Marry, V.; Dubois, E.; Turq, P. Diffusion of water in clays on the microscopic scale: modeling and experiment. J. Phys. Chem. B 2006, 110, 3206–3214. (7) Breu, J.; Seidl, W.; Stoll, A.; Lange, K.; Probst, T. Charge homogeneity in synthetic fuorohectorite. Chem. Mater. 2001, 13, 4213–4220. (8) Bourg, I.; Sposito, G. Connecting the molecular scale to the continuum scale for diffusion processes in smectite-rich porous media. Environ. Sci. Technol. 2010, 44, 2085–2091. 2854
dx.doi.org/10.1021/es1031932 |Environ. Sci. Technol. 2011, 45, 2850–2855
Environmental Science & Technology (9) Bee, M. Quasi-Elastic Neutron Scattering; Principles and Applications in Solid State Chemistry, Biology and Material Science; Adam Hilger: Bristol, England, 1988. (10) ANDRA, Referentiel Materiaux, Tome 1: Contexte et objet; Rapport C.RP.AMAT.01.060, 2001. (11) Malikova, N.; Cadene, A.; Dubois, E.; Marry, V.; Durand-Vidal, S.; Turq, P.; Breu, J.; Longeville, S.; Zanotti, J.-M. Water diffusion in a synthetic hectorite clay studied by quasielastic neutron scattering. J. Phys. Chem. C 2007, 111, 17603–17611. (12) Michot, L.; Bihannic, I.; Pelletier, M.; Rinnert, E.; Robert, J. Hydration and swelling of synthetic Na-saponites: Infuence of layer charge. Am. Mineral. 2005, 90, 166–172. (13) Ferrage, E.; Lanson, B.; Michot, L. J.; Robert, J.-L. Hydration properties and interlayer organization of water and ions in synthetic Nasmectite with tetrahedral layer charge. Part 1. Results from X-ray diffraction profle modeling. J. Phys. Chem. C 2010, 114, 4515–4526. (14) Lanson, B. Private communication. (15) Mezei, F. Lecture Notes in Physics; Springer Verlag: New York, 1980; Vol. 128. (16) Lechner, R. E.; Longeville, S. In Quasielastic Neutron Scattering in Biology Part I: Methods in “Neutron Scattering in Biology, Techniques and Applications”; Fitter, J., Gutberlet, T., Katsaras., J., Eds.; Springer: New York, 2006; pp 309-354. (17) G€ahler, R.; Golub, R. J. Phys. (Paris) 1988, 49, 1195–1202. (18) Berendsen, H.; Grigera, J.; Straatsma, T. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269–6271. (19) Koneshan, S.; Rasaiah, C.; Lynden-Bell, R.; Lee, S. Solvent structure, dynamics, and ion mobility in aqueous solutions at 25 C. J. Phys. Chem. B 1998, 102, 4193–4204. (20) Skipper, N.; Refson, K.; McConnell, J. Computer calculation of water-clay interactions using atomic pair potentials. Clay Miner. 1989, 24, 411–425. (21) Delville, A. Structure of liquids at a solid interface: An application to the swelling of clay by water. Langmuir 1992, 8, 1796–1805. (22) Teppen, B.; Rasmussen, K.; Bertsch, P.; Miller, D.; Schafer, L. Molecular dynamics modeling of clay minerals. 1. Gibbsite, kaolinite, pyrophyllite, and beidellite. J. Phys. Chem. B 1997, 101, 1579–1587. (23) Bougeard, D.; Smirnov, K.; Geidel, E. Vibrational spectra and structure of kaolinite: A computer simulation study. J. Phys. Chem. B 2000, 104, 9210–9217. (24) Smith, D. Molecular computer simulations of the swelling properties and interlayer structure of cesium montmorillonite. Langmuir 1998, 14, 5959–5967. (25) Sainz-Diaz, C.; Hernandez-Laguna, A.; Dove, M. Modeling of dioctahedral 2:1 phyllosilicates by means of transferable empirical potentials. Phys. Chem. Miner. 2001, 28, 130–141. (26) Cygan, R.; Liang, J.-J.; Kalinichev, A. Molecular models of hydroxide, oxyhydroxide, and clay phases and the development of a general force field. J. Phys. Chem. B 2004, 108, 1255. (27) Heinz, H.; Koerner, H.; Anderson, K. L.; Vaia, R. A.; Farmer, B. Force field for mica-type silicates and dynamics of octadecylammonium chains grafted to montmorillonite. Chem. Mater. 2005, 17, 5658–5669. (28) Skipper, N.; Chang, F.-R.; Sposito, G. Monte Carlo simulation of interlayer molecular structure in swelling clay minerals. 1. Methodology. Clays Clay Miner. 1995, 43, 285–293. (29) Marry, V.; Turq, P. Microscopic simulations of interlayer structure and dynamics in bihydrated heteroionic montmorillonites. J. Phys. Chem. B 2003, 107, 1832–1839. (30) Malikova, N.; Marry, V.; Dufr^eche, J.-F.; Simon, C.; Turq, P.; Giffaut, E. Temperature effect in a montmorillonite clay at low hydration - microscopic simulation. Mol. Phys. 2004, 102, 1965–1977. (31) Rotenberg, B.; Marry, V.; Vuilleumier, R.; Malikova, N.; Simon, C.; Turq, P. Water and ions in clays: Unraveling the interlayer/ micropore exchange using molecular dynamics. Geochim. Cosmochim. Acta 2007, 71, 5089. (32) Marry, V.; Rotenberg, B.; Turq, P. Structure and dynamics of water at a clay surface from molecular dynamics simulation. Phys. Chem. Chem. Phys. 2008, 10, 4802.
ARTICLE
(33) Rotenberg, B.; Morel, J.-P.; Marry, V.; Turq, P.; Morel-Desrosiers, N. On the driving force of cation exchange in clays: Insights from combined microcalorimetry experiments and molecular simulation. Geochim. Cosmochim. Acta 2009, 73, 4034–4044. (34) Ateeque Malani, K. G. A. Adsorption isotherms of water on mica: Redistribution and film growth. J. Phys. Chem. B 2009, 113, 1058–1067. (35) Tournassat, C.; Chapron, Y.; Leroy, P.; Bizi, M.; Boulahya, F. Comparison of molecular dynamics simulations with triple layer and modifed Gouy-Chapman models in a 0.1 M NaCl montmorillonite system. J. Colloid Interface Sci. 2009, 339, 533–541. (36) Croteau, T.; Bertram, A. K.; Patey, G. N. Simulation of water adsorption on kaolinite under atmospheric conditions. J. Phys. Chem. A 2009, 113, 7826–7833. (37) Tambach, T.; Hensen, E.; Smit, B. Molecular simulations of swelling clay minerals. J. Phys. Chem. B 2004, 108, 7586–7596. (38) Smith, D.; Wang, Y.; Chaturvedi, A.; Whitley, H. Molecular simulations of the pressure, temperature, and chemical potential dependencies of clay swelling. J. Phys. Chem. B 2006, 110, 20046–20054. (39) Wang, J.; Kalinichev, A.; Kirkpatrick, R.; Cygan, R. Structure, energetics, and dynamics of water adsorbed on the muscovite(001) surface: A molecular dynamics simulation. J. Phys. Chem. B 2005, 109, 15893–15905. (40) Skipper, N. T.; Lock, P. A.; Titiloye, J. O.; Swenson, J.; Mirza, Z. A.; Howells, W. S.; Fernandez-Alonso, F. The structure and dynamics of 2-dimensional fluids in swelling clays. Chem. Geol. 2006, 230, 182–196. (41) Skipper, N. T.; Chang, F.-R. C.; Sposito, G. Monte Carlo simulation of interlayer molecular structure in swelling clay minerals. 1. Methodology. Clays Clay Miner. 1995, 43, 285. (42) Wang, J.; Kalinichev, A.; Kirkpatrick, R. Effects of substrate structure and composition on the structure, dynamics, and energetics of water at mineral surfaces: A molecular dynamics modeling study. Geochim. Cosmochim. Acta 2006, 70, 562–582. (43) Breu, J.; Seidl, W.; Stoll, A. Disorder in smectites in dependence of the interlayer cation. Z. Anorg. Allg. Chem. 2003, 629, 503–515. (44) Laage, D. Reinterpretation of the liquid water quasi-elastic neutron scattering spectra based on a nondiffusive jump reorientation mechanism. J. Phys. Chem. B 2009, 113, 2684–2687. (45) Malikova, N.; Longeville, S.; Zanotti, J.-M.; Dubois, E.; Marry, V.; Turq, P. Signature of low-dimensional diffusion in complex systems. Phys. Rev. Lett. 2008, 101, 265901. (46) Holz, M.; Heil, S.; Sacco, A. Temperature-dependent selfdiffusion coeffcients of water and six selected molecular liquids for calibration in accurate 1H NMR PFG measurements. Phys. Chem. Chem. Phys. 2000, 2, 4740–4742. (47) Laage, D.; Hynes, J. T. On the molecular mechanism of water reorientation. J. Phys. Chem. B 2008, 112, 14230–14242. (48) Sterpone, F.; Stirnemann, G.; Hynes, J. T.; Laage, D. Water hydrogen-bond dynamics around amino acids: The key role of hydrophilic hydrogen-bond acceptor groups. J. Phys. Chem. B 2010, 114, 2083–2089.
2855
dx.doi.org/10.1021/es1031932 |Environ. Sci. Technol. 2011, 45, 2850–2855
