A Molecular Dynamics Study of Structure and Short-time Dynamics of

Angle brackets denote averaging over molecules α and time origins t0. Note that 1/3T(0) gives the self-diffusion coefficient of molecules, according ...
0 downloads 0 Views 181KB Size
5266

J. Phys. Chem. B 1999, 103, 5266-5273

A Molecular Dynamics Study of Structure and Short-time Dynamics of Water in Kaolinite Konstantin S. Smirnov† and Daniel Bougeard* LASIR-CNRS, Baˆ t. C5, Centre d’Etudes et de Recherches Lasers et Applications, UniVersite´ de Sciences et Technologies de Lille, 59655 VilleneuVe d’Ascq Cedex, France ReceiVed: January 4, 1999; In Final Form: March 17, 1999

The method of molecular dynamics was used to investigate the structure and short-time dynamics of the interlayer water molecules in kaolinite structures with 8.5 and 10 Å clay layer spacings. Two types of adsorbed water molecules were identified, that are characterized by different orientations with respect to the surface of the clay sheets. The first type of the molecules is adsorbed on the surface of the silicate clay sheet, and oriented with its molecular HH vector parallel to the surface and with molecular dipole inclined by 30° to the surface normal. For the second type of molecules the HH vectors and the molecular dipoles are perpendicular to the surface and the surface normal, respectively. The increase of the interlayer spacing and of the number of water molecules in the interlamellar space leads to the formation of an intermediate layer. Molecules of this layer form weak hydrogen bonds with hydroxyls of the octahedral clay sheet and strongly interact with the molecules adsorbed on the surface of the silicate sheet. The dynamics of the water molecules in the interlamellar space are analyzed in terms of the molecular center-of-mass motion, the orientational relaxation times, and the diffusion coefficient of the molecules. The results show significant decrease of the diffusion coefficient and increase of the relaxation times relative to those in bulk water.

1. Introduction Interest in studies of water-clay systems arises from the role that interlayer water plays as a solvent for intercalated species. Understanding this role is closely related to the knowledge of the molecular structure and the dynamical behavior of the water at clay-water interfaces and the extent to which the clay surface affects the water properties. Modeling techniques can significantly contribute to obtain this information at the microscopic level. This allows a more detailed interpretation of the experimental data and often a description that can hardly be derived by experimental methods. Several studies on water-clay systems have been conducted with the use of different computational techniques.1-9 The Monte Carlo (MC) method was used by Skipper et al.1-5 to study water confined between the layers of 2:1 clay minerals with model systems based on Li, Na, and Mg smectites. Some of the water molecules were found to form a hydration shell around the charge balancing cations.1 Due to their large positive charges, the Mg2+ cations were shown to affect the structure of interlayer water to a larger extent than Na+ ions. The water molecules not involved in the hydration of the cations were calculated to be hydrogen bonded to the clay surface; they reside in hexagonal rings formed by SiO4 tetrahedra of the clay structure.1 Later on, using the methodology developed in ref 2, the authors showed that the structure of the interlayer water is influenced by both the charge of the silicate layers and the location of the charge.3 In a recent paper,4 this approach provided an insight into the role of K+ ions as swelling inhibitors of clays. A Monte-Carlo study of water molecules at an uncharged surface of talc5 revealed that the water molecules † Permanent address: Institute of Physics, St. Petersburg State University, St. Petersburg 198904, Russian Federation. * Corresponding author: LASIR-CNRS, Baˆt. C5, Universite´ de Sciences et Technologies de Lille, 59655 Villeneuve d’Ascq Cedex, France. Fax: +33 (0)3 20.43.67.55. E-mail: Daniel.Bougeard@univ_lille1.fr.

are adsorbed in the center of rings formed by SiO4 tetrahedra of the surface and form a hydrogen bond with the clay’s hydroxyl groups. The calculation also showed that the structural modification of hydrogen-bond network is limited by a second layer of adsorbed water molecules. Delville6 used grand canonical MC simulations to investigate water molecules confined between two potassium mica particles. The structure of the interlayer water was calculated to be better organized between neutral mica particles (with intercalated K+ ions) than between charged particles with potassium cations solvated in the interlamellar space. In agreement with experimental data, the calculations revealed an unswelling behavior of the mica-potassium-water system. In a recent work,7 the author studied wetting of clay surfaces and showed that their wettability is related to the surface charge and the chemical composition. The method of molecular dynamics was used by Aicken and co-workers8 in their study of anionic species intercalated in layered double hydroxides. The variation of the interlayer spacing with the water content was in a good agreement with experimental data for all systems and the authors could determine the arrangement and the orientation of the guest species inside the layered double hydroxides. Ab initio total energy calculations were performed by Bridgeman et al.9 to find out the position of a water molecule on an uncharged surface of clay minerals. For talc structure the calculations showed that water molecule is situated between the clay sheets with its plane parallel to the clay surface. The oxonium ion H3O+ was found to correspond to an unstable state. Semiempirical and density-functional studies carried out by Chatterjee et al.10 on a montmorillonite cluster interacting with a water molecule revealed the lowest energy for a configuration with two hydrogen atoms pointing toward the clay surface. Many of the simulation studies mentioned above have been carried out for structures that contain charge balancing interlayer

10.1021/jp9900281 CCC: $18.00 © 1999 American Chemical Society Published on Web 06/09/1999

Water in Kaolinite

Figure 1. . Structure of a kaolinite sheet formed by 3 × 2 × 1 unit cells.

ions. The latter significantly affect the characteristics of the interlayer water. In smectites, for example, the water behaves nearly like an aqueous solution.11 Thus, it is difficult to determine the extent to which the surfaces of clay sheets influence the adsorbed water, although some studies indicate that this influence is rather limited.5 Among the clay structures, mineral kaolinite provides an example of a solid with a 1:1 uncharged layer (one sheet of tetrahedrally coordinated atoms: one sheet of octahedrally coordinated atoms) and thus, without charge-balancing cations in the interlayer void. The structure of the layer consists of a sheet formed by six-membered silicate rings (6R) linked to a sheet of octahedrally coordinated Al atoms (Figure 1). One basal surface of the layer is therefore composed of the oxygen atoms of the 6R rings whereas the opposite surface of the layer is covered by OH groups bound to the hexacoordinated Al atoms of the octahedral sheet. Thus, modeling water molecules confined between the kaolinite layers gives access to the influence of the silicate surface on the interlayer water. The data obtained for the water-kaolinite system can therefore serve as a base for comparison with other water-clay systems. This work presents results of an MD simulation study of the water-kaolinite system that is aimed at determining the structure of the interlayer water and the dependence of the structure on the interlayer spacing, as well as getting insight into the shorttime dynamics of water molecules confined in the clay mineral. The results of the simulations are compared with available experimental data and, whenever possible, with other modeling studies of water-clay systems. Although kaolinite and its polymorphs with additional water molecules between the layers are commonly named halloysites,12 we shall keep the name “kaolinite” for the specific structure used in the calculations throughout this paper. 2. Experimental Information

J. Phys. Chem. B, Vol. 103, No. 25, 1999 5267 The dynamical structure of water in halloysites was studied by the NMR technique.14,15 For endellite, Cruz et al.14 found that the motion of the interlayer water can be separated into a rapid 180° flipping motion around the C2 axis of the molecule and a slow tumbling motion about an axis perpendicular to the C2 axis. The activation energies for both motions were found to be 4.6 and 3.6 kcal mol-1 for the flipping and tumbling, respectively. In addition, their IR measurements allowed to conclude that there exists only a weak hydrogen bond between the water molecules and the surface hydroxyls, whereas the molecules are strongly hydrogen bonded to each other. A similar study performed by Lipsicas and co-workers for kaolinite15 showed that the flipping motion of associated water molecules has an activation energy of 4.5 kcal mol-l, whereas for the hole water it is equal to 3.1 kcal mol-1. An activation energy of 2.4 kcal mol-1 was found for the tumbling motion. The results of the study revealed that the translational motion of both types of water is slower than in bulk water. The results of both Cruz et al.14 and Lipsicas et al.15 showed that the correlation time for the flipping motion is smaller than that for the tumbling motion at least by one order of magnitude. 3. Model and Computational Procedure The kaolinite structure used in the present work is based on data of Young and Hewat16 and it has a unit cell of P1 symmetry with lattice parameters a ) 5.149 Å, b ) 8.934 Å, c ) 7.384 Å, and R ) 91.93°, β ) 105.04°, and γ ) 89.79°. The unit cell has a chemical composition Si4Al4O10(OH)8 and the simulation box used in the calculations consisted of 3 × 2 × 2 unit cells (in total 408 atoms). During the calculations, the atoms were kept fixed at their crystallographic positions. The potential model employed for water molecules is a harmonic version of modified SPC model of Berendsen et al.17 that was used by Toukan and Rahman18 in their MD study of the atomic dynamics in liquid water. The geometry of water molecules in the SPC model is described by the O-H equilibrium distance ROH ) 1 Å and the H-O-H equilibrium angle RHOH ) 109.5°. The interaction energy UAB between molecules A and B is the sum of an electrostatic Coulombic potential between charges located on the oxygen and hydrogen atoms of different molecules and a Lennard-Jones (12-6) potential arising from interaction between the oxygen atoms

UAB )

qiqj

∑ ∑ iA jB r

ij

Extensive experimental investigation of the structural characteristics of interlayer water in kaolinite has been carried out by Costanzo and co-workers13 (see also references cited therein). Their study led to the identification of two types of interlayer water. For kaolinite structures with a clay layer spacing d(001) ) 8.4-8.6 Å the water molecules were found to form hydrogen bonds with the basal oxygens of the surface of the tetrahedral silicate sheet. It was suggested that the C2 axes of the molecules are nearly perpendicular to the surface normal. This type of molecules was called hole water. Further, when the hydroxyl groups of the OH covered clay surface were substituted by fluorine atoms, the IR measurements indicated that another type of hole water with a different orientation might be present in the interlayer space. For structures with 10 Å spacing, in addition to the hole water, a layer of associated water is present in the structure that interacts with both the hole water and the hydroxyls of the octahedral sheet. This associated water was suggested to form a strong hydrogen bond with other water molecules and to only weakly interact with the clay surface.

+ 4OO

[( ) ( ) ] σOO rOO

12

-

σOO rOO

6

(1)

where qi and qj are the effective atomic charges of atoms i and j, respectively, and OO and σOO are Lennard-Jones potential parameters. All the parameters used for the intramolecular and intermolecular interactions corresponded to the flexible SPC model.17,18 The interaction of atoms of a water molecule with atoms of the clay was described with the sum of the electrostatic potential between the H and O atoms of molecules and atoms of the clay, and of the Lennard-Jones (12-6) potential between the oxygen of the molecule and atoms of the clay structure (eq 1). Effective atomic charges employed for the clay atoms were taken to be equal to the half of the Mulliken 6-21G charges computed in a periodic Hartree-Fock study of kaolinite.19 The used values of charges are of the same order as those employed for lattice atoms in previous modeling studies of water-clay systems.1-3 The small charge of the hexacoordinated Al atom as compared to the charge of tetrahedral Si is accounted for an effective screening of the aluminium by six surrounding oxygens.

5268 J. Phys. Chem. B, Vol. 103, No. 25, 1999

Smirnov and Bougeard

TABLE 1: Parameters of Equation 1 Used in Calculations atoma

charge, |e|

, K

σ, Å

Ow Hw Si Al Oc Hc

-0.82 +0.41 +1.25 +1.05 -0.85 +0.28

78.195

3.166

47.803 32.707 78.195

3.951 4.112 3.166

a Subscripts “w” and “c” denote atoms of water molecule and clay structure, respectively.

The SPC model was used to describe interactions between the atoms of the water molecule and the O and H atoms of the clay structure. The similar values of the charges used for both water and clay oxygens justify this approximation. Such a description is often invoked to model water-clay systems,1-3 when the interaction between the OH groups of solid and water molecules is represented via intermolecular water potential. The parameters of the Lennard-Jones potential between the Si and Al atoms and the oxygen atom of water were calculated from the standard combination rules using atomic Lennard-Jones parameters of the Dreiding20 and SPC17 force field models for Si, Al atoms, and for water oxygens, respectively. All potential parameters are listed in Table 1. The molecular dynamics calculations were carried out as follows. At the beginning of a calculation the water molecules were uniformly distributed in the interlayer space and initial velocities were chosen from a Maxwell-Boltzmann distribution at temperature 750 K. During subsequent 40000 time steps (20 ps) the system was cooled down to 300 K by rescaling the atomic velocities. The calculation was then continued in the NVE statistical ensemble for 45960 time steps (22.98 ps) and information about coordinates and velocities of the atoms in the system was stored each 4th step for the last 40960 time steps. Relative root mean square fluctuations of the total energy during the NVE stage of the calculation were less than 10-4. The classical equations of atomic motion were integrated with the velocity form of the Verlet algorithm21 with a time step of 0.5 fs and the periodical boundary conditions were applied to the simulation box in three dimensions. A cut off radius of 7.7 Å was used for the Lennard-Jones interactions and a shiftedforce modification22 of the potential was employed to smooth out both the potential and the force at the cut-off distance. The electrostatic interactions were treated by the Ewald sum technique.22 The structure of interlayer water was characterized by (i) density profiles of the oxygen and hydrogen atoms of water molecules along the direction perpendicular to the surface of the clay sheet. This direction (p) was defined as vector product of the a and b crystallographic axes vectors p ) [a b]; (ii) the orientations of the molecular dipole (θ angle) and the HH interatomic vector (φ angle) with respect to the p axis (Figure 2) as a function of a molecular center-of-mass (CoM) position along the p direction:

θ ) arccos eµ ‚ ep and φ ) arccos |eHH ‚ ep|

(2)

where eµ, eHH, and ep are unit vectors along the dipole µ, HH vector, and p directions, respectively; (iii) the pair radial distribution functions (RDF) gRβ(r) between atoms of type R and β

gRβ(r) )

V Nβ(r) Nβ4πr2∆r

(3)

Figure 2. Definition of the orientation of a water molecule with respect to the p axis. See text for detail.

where Nβ(r) is the number of atoms of type β in a spherical layer of thickness ∆r at distance r from an atom of type R chosen as a center of a frame, Nβ is the number of atoms of type β in the system, and V is the volume of the system. The RDF can be used to calculate the running coordination number nRβ(r) of atoms R:

nRβ(r) )

∫0rgRβ(F)F2 dF

4π Nβ V

(4)

The dynamical behavior of water molecules was investigated by calculating (i) the power spectrum of molecular CoM

T(ω) )

∫dt eiωt 〈VR(t0) ‚ VR(t0 + t)〉,

and VR )

1 M

∑mi vRi

(5)

where vRi is the velocity of the atom i with a mass mi, and M signifies the molecular mass of the water molecule. Angle brackets denote averaging over molecules R and time origins t0. Note that 1/3T(0) gives the self-diffusion coefficient of molecules, according to the time-correlation function formalism.23 Due to the relatively short duration of the MD calculations and the low mobility of the interlayer water, the coefficient cannot be calculated by the Einstein relationship; (ii) the power spectrum of the molecules in the molecular frame:

V R(ω) )

∫dt eiωt 〈uRi (t0) ‚ uRi (t0 + t)〉,

where uRi ) vRi - VR (6)

Indexes R and i run over molecules and atoms in a molecule, respectively. (iii) the orientational correlation functions

Cx2 ) 〈P2[ex(t0) ‚ ex(t0 + t)]〉

(7)

where ex is a unit vector directed along the x axis (x being either the molecular dipole µ or the HH interatomic vector) and P2 signifies the second-order Legendre polynomial. The zerofrequency component of the Fourier transform of Cx2(t) yields the relaxation time τx of the corresponding molecular motion. When Cx2(t) is fitted by an exponential function A exp(-t/tx), the relaxation time can be calculated as τx ) A ‚ tx. If x corresponds to the axes directed along molecular dipoles (rotation about the HH vector), the relaxation time can then be compared with that measured in NMR experiments.24 Two structures with the clay layer distances of 8.5 and 10 Å were considered that correspond to the structures studied experimentally.13,15 For the former system, 12 water molecules per aluminosilicate layer (in total 24 molecules) were calculated in the simulations that corresponds to one molecule per 6R ring, a loading estimated by Costanzo et al.13 For the 10 Å system,

Water in Kaolinite

J. Phys. Chem. B, Vol. 103, No. 25, 1999 5269

Figure 4. Orientation of the two types of adsorbed water molecules with respect to a six-membered silicate ring of kaolinite structure. Dotted lines show H-bonds between the molecules and the ring.

Figure 3. Density profiles of the atoms of the water molecules and of the molecular CoM along the p-coordinate for kaolinite structure d(001) ) 8.5 Å. (b) Dependence of the angles θ and φ on the p coordinate. Arrow on the left side of the p axis indicates the positions of the basal hydrogen atoms of the OH covered surface of the octahedral clay sheet. That on the right side of the p axis indicates the positions of the basal oxygen atoms of the silicate clay sheet.

40 molecules per layer (in total 80 molecules) were simulated. According to experimental data13 the 10 Å system contains twice as much water as 8 Å kaolinite. Thus, the latter model studied can serve as a fully hydrated form of 10 Å halloysite. Control calculations performed for the structures with different numbers of water molecules in the interlayer space showed qualitatively similar results. An additional calculation of a system of 108 SPC water molecules at T ) 300 K and density of 1 g/cm3 was carried out, to compare the characteristics of confined water molecules with those of liquid water. 4. Results 4.1. Structure. 8.5 Å System. Figure 3a shows the density profiles of the CoM, hydrogen, and oxygen atoms of the water in the water-kaolinite systems with 8.5 Å interlayer spacing. The figure demonstrates that the molecules form two layers in the proximity of the silicate and OH covered clay surfaces, with most of the molecules adsorbed on the former surface (on average 9 molecules out of 12). Note that the CoM density profile for this layer has a bimodal distribution. The difference of ca. 0.1 Å between the maxima in the distribution roughly corresponds to a difference that can be expected when a water molecule changes its orientation from that with the HH vector perpendicular to the surface to an orientation with the dipole perpendicular to the surface. Figure 3b shows the distribution of the θ and φ angles and indeed reveals that, according to the distribution of the CoM profile, there are two types of water molecules in the silicate layer. The first type is characterized by mean values θ ) 30° and φ ) 85°, i.e., the HH vector of the molecules is parallel to the surface while the C2 axis is slightly tilted with respect to the p direction. The second type of molecules has an orientation with θ ) 100° φ ) 30°. The

Figure 5. Density profiles of the atoms of the water molecules and of the molecular CoM along the p coordinate for the kaolinite structure with d(001) ) 10 Å. Dependence of the angles θ and φ on the p coordinate. See legend to Figure 2.

water molecules adsorbed on the OH covered surface are characterized by an orientation with mean values of 90° and 12° for θ and φ angles, respectively, similar to the molecules of the second type. The visualization of the structure of the interlayer water during the MD simulation revealed that molecules of the first type are adsorbed in the six-membered silicate rings with both hydrogen atoms of the molecules involved in H-bonds to the oxygens of the rings. This interaction orients the molecules so that the HH vector is parallel to the surface. Molecules of the second type are bound by one hydrogen atom to an oxygen of a 6R silicate ring while a second hydrogen bond is formed between the oxygen of the molecule and a hydrogen of a hydroxyl group of the OH covered kaolinite surface. This leads to the orientation of the molecules with the HH vector nearly perpendicular to the surfaces. A snapshot illustrating the two types of adsorbed molecules and their orientation with respect to a 6R silicate ring is presented in Figure 4. 10 Å System. Figure 5a shows density profiles for CoM, O, and H atoms of the interlayer water. The figure demonstrates that, in addition to the two above mentioned layers, a third layer is formed in the intersheet space. Molecules in this layer are

5270 J. Phys. Chem. B, Vol. 103, No. 25, 1999

Smirnov and Bougeard

Figure 6. A snaphot of a layer of water molecules in the 10 Å waterkaolinite system.

characterized by orientation (Figure 5b) with mean values of the θ and φ angles close to 50°. Note that the formation of the intermediate layer mostly affects positions of the water molecules adsorbed on the OH covered surface of the octahedral sheet. The peak at -4.5 Å in the hydrogen atoms distribution (Figure 5a) shifts by ca. 0.35 Å with respect to the position of the hydrogen atoms of the OH groups as compared to that in the 8.5 Å system. A still larger shift is observed for the position of the corresponding peak in the O atoms distribution: it shifts by ca. 0.6 Å. A snapshot of the interlayer water during the MD simulation is presented in Figure 6 and it clearly demonstrates both the positions and the orientation of the molecules in the space between the clay sheets. The figure reveals that the orientation of the molecules in the third layer is due to the formation of hydrogen bonds between the H atom of the molecule and the oxygen of molecules adsorbed in 6R silicate rings. The oxygen and the second hydrogen atoms of the molecule in the intermediate layer are able to interact with both atoms of the molecules adsorbed on the OH covered surface of the clay sheet and with atoms of other molecules of the layer. This stronger interaction with the molecules adsorbed on the surface of the octahedral sheet obviously accounts for the changes observed in the distribution of atoms in the interlamellar space. 4.2. Dynamics. 8.5 Å System. Figure 7a shows the power spectrum of the CoM of the interlayer water as compared to the CoM spectrum of molecules in the liquid state. They strongly differ from each other. Firstly, the spectral intensity in the former spectrum is concentrated in a narrower frequency region shifted to higher wavenumbers as compared to the latter. Secondly, a structure is present in the spectra for different reasons. The two bands at 48 and 190 cm-1 in the CoM spectrum of molecules in liquid water belong to the O-O-O and O-O vibrations, respectively,25 whereas the bands at 180 and 260 cm-1 in the spectrum of interlayer water are due to the anisotropy of the motion of confined molecules. The analysis of the components of the spectrum reveals that the high-frequency band belongs to a CoM vibration perpendicular to the surfaces, whereas vibrations parallel to the surfaces appear as a broad band at 180 cm-1. This anisotropy in the CoM motion is obviously due to the anisotropy of the interlayer space in the clay: motions of molecules perpendicularly to the clay sheets is strongly

Figure 7. CoM power spectrum of water molecules in the (a) 8.5 Å and (b) 10 Å structures. Solid line and dashed line signify total power spectrum and spectrum calculated for molecular motion perpendicular to the sheet surface, respectively. The solid circles connected by line in the pannel “a” show the spectrum for molecules in liquid.

restricted. Note that a negligible intensity in the CoM spectrum of the interlayer water indicates the absence of diffusion, i.e., molecules are strongly fixed at their adsorption positions. A wide vibrational band occurs at ca. 500 cm-1 in the spectrum of the atoms in the molecular frame of both interlayer and liquid water (not presented). Similarly to the CoM spectrum of confined molecules, the spatial anisotropy leads to a higher value of the vibrational frequency for motions perpendicular to the clay surfaces. 10 Å System. The CoM spectrum of the molecules is shown in Figure 7b. In comparison to the spectrum of 8.5 Å system, that of the 10 Å system loses its features, and the figure displays that the spectrum of motions perpendicular to the surface broadens and shifts to lower wavenumbers. Although the spectrum is still different from that for liquid water, the increase of interlayer spacing and number of molecules results in more isotropic CoM motions, despite a relatively ordered arrangement of the molecules in the interlayer space. The increase of the interlayer spacing results in a broadening of the vibrational band in the spectrum of the molecules in the molecular frame and the spectra of the different spatial components become practically undistinguishable. For both systems the self-diffusion coefficient D was calculated to be less than 10-10 m2 s-1 as compared to D ) 3 × 10-9 m2 s-1 computed for bulk water molecules. An estimation of the diffusion coefficient is given rather than an exact value, because the latter might be subject to a substantial uncertainty due to the vibrational behavior of the CoM motion. Figure 8 presents the Cx2(t) functions calculated for both systems. The functions have a short-range oscillatory part that is followed by a long-range slow decay. The relaxation time of the corresponding molecular motions were calculated by fitting the curves by the exponential function on interval 0.2-4 ps to avoid a bias due to the fast decrease of the functions during the

Water in Kaolinite

J. Phys. Chem. B, Vol. 103, No. 25, 1999 5271

Figure 9. gHWOC radial distribution function for the kaolinite structures with 8.5 Å (solid line) and 10 Å (dashed line) interlayer spacing. The arrow indicates the distance used for the calculation of the running coordination number given Table 3.

Figure 8. Cµ2 (t) functions (upper panel) and CHH 2 (t) functions (lower panel) calculated for water molecules in kaolinite structures with the interlayer spacing of (a) 8.5 Å, (b) 10 Å, and for molecules in a (c) liquid state.

TABLE 2: Relaxation Times (in Picoseconds) of the molecular dipole (τµ) and the HH vector (τHH) relaxation time

8.5 Å

10 Å

liquid H2O

τµ τHH

29.2 24.2

13.4 10.6

1.30 1.54

initial time. Table 2 lists relaxation times calculated for the systems and compare them with the relaxation times computed for molecular motions in liquid water. 5. Discussion The molecules which are hydrogen bonded to the oxygens of the 6R silicate rings with their C2 axes slightly inclined with respect to the [001] direction correspond to the hole water of Costanzo et al.13 The existence of this type of adsorbed water was also inferred in MC calculations of water-smectites systems by Skipper and co-workers.1 A minimum-energy configuration with two water hydrogens pointing towards the clay surface was calculated by Chatterjee et al.10 The existence of a second type of adsorbed water with its HH vector perpendicular to the surface was also suggested by Costanzo and co-workers13 in 8.4-8.6 Å kaolinite samples with fluorine atoms substituting the OH groups of the clay structure. The present calculations suggest that this type of adsorbed water can also exist in hydrated kaolinites. Note that the OH group is isoelectronic to a fluorine atom and a similar strength of interaction between the hydroxyl water molecule and the surface can therefore be expected in both cases. The formation of the intermediate layer of water at the increase of the interlayer spacing is also in line with experimental data showing the presence of associated water in 10 Å kaolinite.13 In the case of smectites the molecules of this layer form a coordination sphere around the charge balancing cations and behave like an aqueous solution.1,11 The present work shows that in the absence of the charge-balancing cations the molecules of the associated water form hydrogen bonds with the adsorbed molecules, and the structure of the layer is determined by this

Figure 10. gOWHC radial distribution function for the kaolinite structures with 8.5 Å (solid line) and 10 Å (dashed line) interlayer spacing. The arrow indicates the distance used for the calculation of the running coordination number given Table 3.

interaction. It is interesting to note that the layer of associated water lies closer to the OH covered surface of the clay (Figures 5a and 6). This is due to the fact that both the OH covered surface and the water molecules adsorbed on it provide more possibilities of forming hydrogen bonds than the molecules of hole water, since both hydrogen atoms of the latter are already involved in H bonds with clay oxygens. The analysis of the radial distribution function between atoms of the molecules and atoms of the clay calculated for both systems allows us to elucidate the influence of the interlayer spacing on the structure of the adsorbed layers. Thus, the RDF for water (W) hydrogens and clay (C) oxygens (gHWOC) and that for water oxygens and clay hydrogens (gOWHC) demonstrate the influence of the spacing on water molecules adsorbed on the silicate and OH covered surfaces, respectively. Figure 9 shows that the formation of the layer of the associated water mostly causes a decrease of intensity of the first peak in gHWOC, while the position of the peak remains unchanged. In both 8.5 Å and 10 Å systems the peak in the RDF appears at a distance of 1.7 Å that is close to the H-bond length in liquid water and indicates a similar strength of the interaction. The 0.05 Å shift of the peak to lower distances is obviously due to the confinement of the molecules in the interlayer space. More drastic changes occur with the variation of the spacing in the gOWHC RDF (Figure 10). Thus, for the system with 8.5 Å interlayer spacing, bonding between the atoms is characterized by hydrogen bonds with two characteristic bond lengths 1.6 and 1.9 Å. The formation of the intermediate layer in the 10 Å system leads to a weakening of these hydrogen bonds because

5272 J. Phys. Chem. B, Vol. 103, No. 25, 1999

Smirnov and Bougeard TABLE 3: Calculated Running Coordination Numbers of Atoms of Interlayer Water atoma

atom

8.5 Å

10 Å

Hw

Ow Oc Hw Hc

0.12 1.00 0.23 0.72

0.49 0.54 0.98 0.42

Ow

a Subscripts “w” and “c” denote atoms of water molecule and clay structure, respectively.

Figure 11. gHWOW radial distribution function for the kaolinite structures with 8.5 Å (solid line) and 10 Å (dashed line) interlayer spacing. Open circles connected by line show gHO RDF for atoms in liquid water. The arrow indicates the distance used for the calculation of the running coordination number given Table 3.

new bonds are formed between the adsorbed molecules and molecules of the associated water. Although the orientation for the absorbed molecules does not change significantly (Figure 5b), the position of the first peak in the RDF at 2.04 Å (Figure 10) suggests that the molecules become less bound to the OH covered surface, whereas they more strongly interact with molecules of the associated water layer. Therefore, the calculated RDFs indicate that the formation of the intermediate layer mostly influences the molecules adsorbed on the OH covered surface of the clay sheets. In fact, these molecules cannot be considered any longer as being “adsorbed” but they rather form, together with molecules of associated water, a unique layer weakly interacting with the hydroxyls of the octahedral clay sheet. Assuming that the position of the peaks in RDFs reflects the strength of the hydrogen bonds, one can infer that the energy of the hole water is stronger than that for molecules interacting with the OH covered surface. This conclusion is in agreement with experimental observations.13,14 Differences in the structure of the interlayer water for the two systems become evident by analyzing the gHWOw radial distribution function (Figure 11). In the case of the 8.5 Å system, the first peak in the distribution function occurs at a distance of 2.1 Å that is 0.35 Å larger than the 1.75 Å hydrogen bond peak in the RDF of liquid water (Figure 11). This indicates that no real hydrogen bond between water molecules is formed. In this system the formation of hydrogen bonds is driven by the ability of the molecules to migrate on the clay surface; it is very limited since the molecules are strongly bound to the surfaces. The increase of both interlayer spacing and number of interlayer molecules leads to a structure that is characterized by a gHWOW similar to that of liquid water. However, the geometrical confinement of the molecules results in a noticeable increase of intensity of the first peak as compared to its intensity in the RDF of liquid water.26 The confinement leads to the formation of hydrogen bonds with a bond length of 1.68 Å, i.e., 0.07 Å shorter than in the liquid state. Table 3 lists the running coordination numbers of atoms of the interlayer water: 3.19 and 3.46 hydrogen bonds per molecule are formed in 8.5 and 10 Å systems, respectively. Both values are noticeably lower than 3.78 H-bonds per molecule calculated for liquid water. The comparison of the relaxation times listed in Table 2 shows that the confinement of the molecules leads to an increase of the relaxation times at least by one order of magnitude as compared with those for molecules in a liquid state. The calculations also reveal that the tumbling motion, following the

notation of refs 14 and 15, has a longer relaxation time (τµ) than that of the flipping motion (τHH). This is in line with the data of the experimental studies.14,15 However, the simulations fail to reproduce the large difference in these relaxation times found in the experiments, although Table 2 indicates that the difference between the relaxation times for these two motions increases with the decrease of the interlayer spacing. Computational studies of the vibrational spectra of complex aluminosilicate lattices27 show that the dynamical behavior of a system is more sensitive to the potential model than structural characteristics. Calculations using different potential/structural models would be of interest to check the results of the present calculations. So far, most of the modeling studies of waterclays systems were aimed at the analysis of the structural characteristics of these systems. 6. Conclusions The structure and short-time dynamics of the interlayer water molecules in kaolinite structures with 8.5 and 10 Å interlayer spacings were studied by the molecular dynamics method. For small interlayer spacing, the calculations reveal two types of adsorbed water molecules that are characterized by different orientation with respect to the surface of the clay sheets. The first type of molecules is adsorbed on silicate surface and oriented with its molecular HH vector parallel to the surface and with a mean value of the angle between molecular dipole and normal to the surface equal to 30°. The second type of the adsorbed water has its HH vector perpendicular to the surface, whereas the dipoles of the molecules are almost perpendicular to the surface normal. Such molecules exist on surfaces of both the tetrahedral and octahedral sheets. Due to the confinement effect, the dynamics of the molecules behaves strongly anisotropic. The increase of the interlayer spacing and of the number of water molecules in the interlamellar space leads to the formation of a layer of associated water. This layer is characterized by a weak interaction with hydroxyls of the octahedral clay sheet, and by a strong interaction with the molecules adsorbed on the surface of the silicate sheet. The formation of this layer results in a softening of the motion of the molecular center of mass and to a more isotropic dynamics of the interlayer water. The calculations show that the diffusion coefficient of the interlayer water is significantly lower than that of liquid water, whereas the relaxation times of the molecular motions are increased. Acknowledgment. K.S.S. gratefully acknowledges research fellowships from the “Re´gion Nord-Pas de Calais” and “Centre National de la Recherche Scientifique”. References and Notes (1) ) Skipper, N. T.; Refson, K.; McConnell, J. D. C. J. Chem. Phys. 1991, 94, 7434.

Water in Kaolinite (2) Skipper, N. T.; Chang, F.-R. C.; Sposito, G. Clay Clays Miner. 1995, 43, 285. (3) Skipper, N. T.; Sposito, G.; Chang, F.-R. C. Clay Clays Miner. 1995, 43, 294. (4) Boek, E. S.; Coveney, P. V.; Skipper, N. T. J. Am. Chem. Soc. 1995, 117, 12608. (5) Bridgeman, C. H.; Skipper, N. T. J. Phys.: Condens. Matter 1997, 9, 4081. (6) Delville, A. J. Phys. Chem. 1993, 97, 9703. (7) Delville, A. J. Phys. Chem. 1995, 99, 2033. (8) Aicken, A.; Bell, I. S.; Coveney, P. V.; Jones, W. AdV. Mater. 1997, 9, 496. (9) Bridgeman, C. H.; Buckingham, A. D.; Skipper, N. T.; Payne, M. C. Mol. Phys. 1996, 89, 879. (10) ) Chatterjee, A.; Iwasaki, T.; Ebina, T.; Hayashi, H. Appl. Surf. Sci. 1997, 121-122, 167. (11) ) Sposito, G.; Prost, R. Chem. ReV. 1982, 82, 553. (12) 2) Newman, A. C. D.; Brown, G. In Chemistry of Clays and Clay Minerals; Newman, A. C. D., Ed.; John Wiley & Sons: New York, 1987; p 1. (13) 3) Costanzo, P. M.; Giese, R. F., Jr.; Lipsicas, M. Clays Clay Miner. 1984, 32, 419. (14) 4) Cruz, M. I.; Letellier, M.; Fripiat, J. J. J. Chem. Phys. 1978, 69, 2018.

J. Phys. Chem. B, Vol. 103, No. 25, 1999 5273 (15) 5) Lipsicas, M.; Straley, C.; Costanzo, P. M.; Giese, R. F., Jr. J. Colloid Interface Sci. 1985, 107, 221. (16) 6) Young, R. A.; Hewat, A. W. Clays and Clay Miner. 1988, 36, 225. (17) 7) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. In Intermolecular Forces; Pullman, B., Ed.; Reidel: Dordrecht, 1981; p 331. (18) 8) Toukan, K.; Rahman, A. Phys. ReV. B 1985, 31, 2643. (19) 9) Hess, A. C.; Saunders, V. R. J. Phys. Chem. 1992, 96, 4367. (20) ) Mayo, S. L.; Olafson, B. D.; Goddard, W. A., III J. Phys. Chem. 1990, 94, 8897. (21) ) Swope, W. C.; Andersen, H. C.; Berens, H. C.; Wilson, K. R. J. Chem. Phys. 1982, 76, 637. (22) Allen, M. P.; Tildesley, D. Computer Simulation of Liquids; Clarendon Press: Oxford, 1987. (23) McQuarrie, D. A. Statistical Mechanics; Harper and Row: New York, 1976. (24) Gordon, R. G. In AdVances in Magnetic Resonance; Waugh, J. S., Ed.; Academic Press: New York, 1968; p 1. (25) Sceats, M. G.; Rice, S. A. J. Chem. Phys. 1980, 72, 3236. (26) Rovere, M.; Ricci, M. A.; Vellati, D.; Bruni, F. J. Chem. Phys. 1998, 108, 9859. (27) Smirnov, K. S.; Bougeard, D. J. Phys. Chem. 1993, 97, 9434.