Molecular Dynamics Studies of the Interactions of Water and Amino

Nov 24, 2008 - Rebecca Notman* and Tiffany R. Walsh. Department of Chemistry and Centre for Scientific Computing, UniVersity of Warwick,. CoVentry CV4...
0 downloads 0 Views 1MB Size
1638

Langmuir 2009, 25, 1638-1644

Molecular Dynamics Studies of the Interactions of Water and Amino Acid Analogues with Quartz Surfaces Rebecca Notman* and Tiffany R. Walsh Department of Chemistry and Centre for Scientific Computing, UniVersity of Warwick, CoVentry CV4 7AL, United Kingdom ReceiVed October 8, 2008. ReVised Manuscript ReceiVed NoVember 24, 2008 The interactions of silica surfaces with water and biomolecules are of considerable significance in bio- and nanotechnology and in geochemistry. An important goal in the fields of biomineralization and biomimetics is to fine-tune these interactions for the control, e.g., of assembly of materials at the nanoscale. Here we report molecular dynamics simulations of fully hydroxylated R-quartz (101j0), (0001), and (011j1) surfaces in explicit water. We also present free energy estimates of adsorbing water and analogues of amino acid side chains onto the quartz (101j0) surface. We find that at least two layers of structured water form on the surface, which is driven by the formation of a strong hydrogen bond network at the interface. Interestingly, we find that the free energy change to move methane (analogue of the side chain of alanine) from bulk water to the (101j0) interface is favorable. We ascribe this to the presence of microscopic voids on the surface, which can accommodate small hydrophobic moieties and shield them from the solvent. These observations draw some useful insights into the possible mechanisms by which biomolecules, in particular peptides and proteins, bind to quartz and other silica surfaces.

1. Introduction The structure of water at aqueous interfaces plays a key role in many aspects of technology, geochemistry, and biological systems.1,2 At the solid-liquid interface, the properties of both substances deviate from the characteristics of the bulk phases, giving rise to a multitude of interesting behaviors.3 In particular, the silica (SiO2)-water interface is important in molecular biology4 and microtechnologies.5 Silica surfaces are also of special interest in the fields of biomineralization and biomimetics. Silica biomineralization occurs in diatoms, silica sponges, and some plants, and much effort has been directed toward understanding the assembly of the silica nanostructures that are formed in these organisms.6-8 In addition, extensive experimental work has been carried out to artificially synthesize silica platelets,9 nanospheres,10 and nanotubes11-13 using, e.g., synthetic peptides,9,10,14 polyamines,15 and dendrimers16 as templates. In recent years, combinatorial biology techniques have been used to select peptide sequences that preferentially bind to one material, such as quartz, over a range of other materials.17-23 In this emerging field of biomi* To whom correspondence should be addressed. E-mail: r.notman@ warwick.ac.uk. (1) Henderson, M. A. Surf. Sci. Rep. 2002, 46, 1. (2) Li, Z.; Lazaridis, T. Phys. Chem. Chem. Phys. 2007, 9, 573. (3) Barnes, G. T.; Gentle, I. R. Interfacial Science: An Introduction; Oxford University Press: Oxford, U.K., 2005. (4) Sasaki, Y. C. Biochem. Soc. Trans. 2004, 32, 761. (5) Papakonstantinon, P.; Vainos, N. A.; Fotakis, C. Appl. Surf. Sci. 1999, 151, 159. (6) Patwardhan, S. V.; Clarson, S. J.; Perry, C. C. Chem. Commun. 2005, 1113. (7) Sumper, M.; Brunner, E. ChemBioChem 2008, 9, 1187. (8) Gro¨ger, C.; Lutz, K.; Brunner, E. Cell Biochem. Biophys. 2008, 50, 23. (9) Glawe, D. D.; Rodriguez, F.; Stone, M. O.; Naik, R. R. Langmuir 2005, 21, 717. (10) Brott, L. L.; Naik, R. R.; Pikas, D. J.; Kirkpatrick, S. M.; Tomlin, D. W.; Whitlock, P. W.; Clarson, S. J.; Stone, M. O. Nature 2001, 413, 291. (11) Baral, S.; Schoen, P. Chem. Mater. 1993, 5, 145. (12) Shenton, W.; Douglas, T.; Young, M.; Stubbs, G.; Mann, S. AdV. Mater. 1999, 11, 253. (13) Harada, M.; Adachi, M. AdV. Mater. 2000, 12, 839. (14) Yuwono, V. M.; Hartgerink, J. D. Langmuir 2007, 23, 5033. (15) Patwardhan, S. V.; Mukherjee, N.; Clarson, S. J. Silicon Chem. 2002, 1, 47. (16) Knecht, M. R.; Wright, D. W. Langmuir 2004, 20, 4728. (17) Brown, S. Nat. Biotechnol. 1997, 15, 269.

metics, the molecular recognition and self-assembly properties of such peptides are used to control nanostructures and assemblies of materials in two and three dimensions.20,24 For these processes to be exploited for potentially wide-reaching applications in bioand nanotechnology, a molecular-level understanding of silica interfaces with water and biomolecules is imperative. Silica surfaces have been extensively studied by experiments. Of significant interest to the present study are experiments on hydrated silica surfaces, which clearly demonstrate that the silica surface perturbs the structure of water at the interface. In many studies R-quartz is the surface of choice as this is the most stable polymorph of silica at ambient conditions. Early nuclear magnetic resonance spectroscopy experiments showed that water is more ordered at the silica surface than in bulk water.25 Sum-frequency generation has indicated that there are several monolayers of water at the R-quartz-water interface if the surface is highly ionized.26 In addition, it was found, using phase-sensitive sumfrequency vibrational spectroscopic techniques, that water in contact with hydrophilic R-quartz may be ordered or disordered, depending on the pH.27 Stålgren et al.28 probed the R-quartz-H2O/ D2O interface using neutron reflectivity and quartz crystal microbalance and provided evidence that the solvent within 5-10 nm of the surface differs from the bulk and that there is preferential enrichment of D2O at the interface. Ordered monolayers of water (18) Whaley, S. R.; English, D. S.; Hu, E. L.; Barbara, P. R.; Belcher, A. M. Nature 2000, 405, 665. (19) Wang, S.; Humphreys, E. S.; Chung, S. Y.; Delduco, D. F.; Lustig, S. R.; Wang, H.; Parker, K. N.; Rizzo, N. W.; Subramonay, S.; Chiang, Y. M.; Jagota, A. Nat. Mater. 2003, 2, 196. (20) Sarikaya, M.; Tamerler, C.; Jen, A. K. Y.; Schulten, K.; Baneyx, F. Nat. Mater. 2003, 2, 577. (21) Sano, K. I.; Shiba, K. J. Am. Chem. Soc. 2003, 125, 14234. (22) Oren, E. E.; Tamerler, C.; Sahin, D.; Hnilova, M.; Seker, U. O. S.; Sarikaya, M.; Samudrala, R. Bioinformatics 2007, 23, 2816. (23) Tamerler, C.; Kacar, T.; Sahin, D.; Fong, H.; Sarikaya, M. Mater. Sci. Eng., C 2007, 27, 558. (24) Evans, J. S.; Samudrala, R.; Walsh, T. R.; Oren, E. E.; Tamerler, C. MRS Bull. 2008, 33, 514. (25) Sermon, P. A. J. Chem. Soc., Faraday Trans. 1980, 76, 885. (26) Du, Q.; Freysz, E.; Shen, Y. R. Phys. ReV. Lett. 1994, 72, 238. (27) Ostroverkhov, V.; Waychunas, G. A.; Shen, Y. R. Phys. ReV. Lett. 2005, 94, 046102. (28) Stålgren, J. J. R.; Boschkova, K.; Ericsson, J.-C.; Frank, C. W.; Knoll, W.; Satija, S.; Toney, M. F. Langmuir 2007, 23, 1943.

10.1021/la803324x CCC: $40.75  2009 American Chemical Society Published on Web 01/06/2009

H2O and Amino Acid Analogue Interactions with Quartz

have also been detected on other silica surfaces, including amorphous silica.29,30 Unfortunately, the experimental data are limited in terms of the molecular-level understanding of silica surfaces. Hence, various computational studies have also been performed to investigate the properties of these systems. Much of this work has focused on the structure of dry surfaces (with the aim of understanding the molecular reconstruction that occurs when the bulk material is cleaved) or on the mechanisms of water dissociation on the surface. For example, the structures and relative energies of different quartz surfaces were modeled atomistically using the ionic core-shell model.31 The authors observed that the dehydrated (0001) surface was the most stable compared with the (101j0), (101j1), and (101j1j) surfaces. Goumans studied the structure and stability of the cleaved, reconstructed, and hydroxylated R-quartz (0001) surface at various thicknesses using periodic density functional theory.32 The hydroxylated R-quartz (0001) surface was found to be much more stable than the bare surfaces, confirming the known strong hydrophilicity of quartz surfaces. In a related study of the same surface by Rignanese et al., the fully hydroxylated surface was found to be energetically favored over the bare or partially hydroxylated surface.33 To fully explore the properties of water at an interface, molecular dynamics (MD) simulations using periodic boundary conditions are necessary. Du and de Leeuw investigated the interactions of bare and hydroxylated R-quartz (0001) surfaces with water by means of MD simulations.34 They observed several ordered monolayers of water on the bare surface and found that, on the hydroxylated surface, interfacial water molecules loosely associate with the surface via hydrogen bonding. Recently, Adeagbo et al. studied the behavior of water confined between two R-quartz (0001) surfaces using Car-Parrinello MD simulations.35 They found that water molecules rapidly reacted with the Si-terminated quartz surface, leading to hydroxylation of both surfaces. In addition, the diffusion of water was found to decrease, the closer its proximity to the surface. MD simulations have also been performed of water on the (111) surface of cristobalite.36 In that study, it was found that water molecules in contact with hydroxylated silica surfaces are oriented in such a way as to maximize hydrogen bonding with the surface. This orientational anisotropy of water on silica was also reported by Lee and Rossky.37 Further conformational ordering of water on R-quartz was reported in an MD study by Yang and Wang, who showed that a monolayer of water on the hydroxylated (0001) surface adopts a flat two-dimensional (2D) structure, where the water molecules are oriented in a H-down configuration to satisfy the requirement of the hydrogen bonding between the water and the surface hydroxyls.38 A 2D icelike structure has also been observed on the hydroxylated β-cristobalite (101j0) surface on the basis of ab initio calculations39 and by classical MD (29) Aarts, I. M. P.; Pipino, A. C. R.; Hoefnagels, J. P. M.; Kessels, W. M. M.; van de Sanden, M. C. M. Phys. ReV. Lett. 2005, 95, 166104. (30) Asay, D. B.; Kim, S. H. J. Phys. Chem. B 2005, 109, 16760. (31) de Leeuw, N. H.; Higgins, F. M.; Parker, S. C. J. Phys. Chem. B 1999, 103, 1270. (32) Goumans, T. P. M.; Wander, A.; Brown, W. A.; Catlow, C. R. A. Phys. Chem. Chem. Phys. 2007, 9, 2146. (33) Rignanese, G.-M.; Charlier, J.-C.; Gonze, X. Phys. Chem. Chem. Phys. 2004, 6, 1920. (34) Du, Z.; de Leeuw, N. H. Dalton Trans. 2006, 2623. (35) Adeagbo, W. A.; Doltsinis, N. L.; Klevakina, K.; Renner, J. Comput. Phys. Commun. 2008, 9, 994. (36) Giovambattista, N.; Debenedetti, P. G.; Rossky, P. J. J. Phys. Chem. B 2007, 111, 9581. (37) Lee, S. H.; Rossky, P. J. J. Chem. Phys. 1994, 100, 3334. (38) Yang, J.; Wang, E. G. Phys. ReV. B 2006, 73, 035406. (39) Yang, J.; Meng, S.; Xu, L. F.; Wang, E. G. Phys. ReV. Lett. 2004, 92, 146102.

Langmuir, Vol. 25, No. 3, 2009 1639

simulations at low temperatures.40 Recently, Lopes et al. developed an MD force field for silica that is compatible with the CHARMM force field.41 The force field was verified by simulations of hydroxylated and H-terminated quartz (101j0) and (011j1) surfaces in water. The significant feature of this work is the move toward the development of force fields that accurately model the inorganic material, the biomolecule, and the interactions between them. The results of the studies described above clearly show that the structure of quartz affects the structure and chemistry of water near the interface. Nevertheless, a clear molecular-level description of the quartz-water and quartz-biomolecule interface is still lacking. Here, we aim to characterize the interaction of water with the fully hydroxylated R-quartz (101j0), (0001), and (011j1) surfaces, hereafter referred to as the (101j0), (0001), and (011j1) surfaces. We extend the analysis of the (101j0) and (011j1) surfaces that was presented by Lopes et al.41 with an emphasis on the effect that each surface has on water structuring at the interface and the possible implications for biomolecule interactions. Furthermore, we model the fully hydroxylated (0001) surface, which has not previously been explored by MD simulations using a classical force field. The current work is motivated by phage display experiments, which have selected peptides that preferentially bind to the quartz (101j0) surface,22 and our long-term goal is to understand the interactions of biomolecules with inorganic surfaces. Hence, to gain preliminary insights into the mechanisms of adsorption of biomolecules on quartz surfaces, we have calculated the potential of mean force (PMF) free energy profiles of water and model hydrophobic (alanine) and hydrophilic (serine) amino acids as a function of the distance from the (101j0) surface. In this study we only considered the side chains of the residues; hence, alanine was modeled as methane and serine as methanol. These studies will be a precursor for future investigations of the mechanisms of peptides binding to quartz. A clear understanding of the role that water plays in these binding processes is vital for understanding and exploiting the forces that drive such interactions.

2. Methodology 2.1. System Details. We have carried out MD simulations of the (101j0), (0001), and (011j1) surfaces using Gromacs.42 The initial structures were built using the InorganicBuilder plugin in VMD43 from an R-quartz unit cell with lattice dimensions a ) 0.491, b ) 0.491, and c ) 0.540 nm, which was replicated along the a, b, and c directions to yield an 8 × 8 × 8 supercell. The supercell was then cleaved along the required plane and rotated so that the surface of interest lay perpendicular to the z axis of the simulation cell. After cleavage, the dangling Si and O atoms were terminated with OH groups or H atoms, which yielded two fully hydroxylated surfaces per slab. The edges of the resultant slabs were also trimmed to generate orthorhombic simulation cells that were fully periodic in the x and y directions. This resulted in (101j0), (0001), and (011j1) surfaces with areas of 5.91, 4.96, and 7.65 nm2 and containing 720, 672, and 1920 SiO2 atoms, respectively. Each slab was approximately 2.5 nm thick (21-27 atomic layers). The surfaces were modeled using the silica force field developed by Lopes et al.41 This force field has been parametrized to be compatible with the CHARMM44,45 force field. One of the significant challenges in modeling bioinorganic interactions is the availability (40) Lu, Z.-Y.; Sun, Z.-Y.; Li, Z.-S.; An, L.-J. J. Phys. Chem. B 2005, 109, 5678. (41) Lopes, P. E. M.; Murashov, V.; Tazi, M.; Demchuk, E.; MacKerell, A. D. J. Phys. Chem. B 2006, 110, 2782. (42) Lindahl, E.; Hess, B.; van der Spool, D. J. Mol. Model. 2001, 7, 306. (43) Humphrey, W.; Dalke, A.; Schulten, K. J. Mol. Graphics 1996, 14, 33. (44) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; Sates, D. J.; Swarninathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187.

1640 Langmuir, Vol. 25, No. 3, 2009

Notman and Walsh

of force fields that adequately describe the inorganic surface, the biomolecule, and the cross-terms. As this work is a precursor to further studies of peptides interacting with quartz, this CHARMMcompatible force field was a natural choice. Each surface was solvated with 580 (101j0), 500 (0001), or 1002 (011j1) explicit water molecules, and the thickness of the resulting water layer was 2.5-3.5 nm. Water molecules were added using the genbox code in Gromacs,42 where the system is immersed in an equilibrated box of water and then water molecules are removed from the box where the distance between any water atom and any silica atom is less than the sum of the van der Waals radii of both atoms. We chose the modified TIP3P potential46,47 as our water model as it is compatible with both the silica force field and the CHARMM force field. For the free energy calculations, amino acid analogues were modeled using the CHARMM force field.44,45 Analogues of polar serine (methanol) and nonpolar alanine (methane) were constructed by replacing the backbone part of each amino acid with a hydrogen and adjusting the adjacent atom type to achieve neutrality. 2.2. Simulation Parameters. MD simulations were carried out using a time step of 1 fs for a total of 10 ns per surface. Coordinates were saved every 1 ps (1000 steps). The first 3 ns of each simulation was used as an equilibration period and removed before subsequent analysis. To calculate diffusion coefficients and the hydrogen bond correlation function, we extended each simulation by 1 ns and saved the coordinates more frequently (every 10 fs). Nonbonded Lennard-Jones interactions were cut off at 1.0 nm, and the electrostatics were treated using particle mesh Ewald, with a real-space cutoff of 1.0 nm. Full periodic boundary conditions were employed in all directions so that the system was modeled as an infinite slab in the xy plane. The simulations were carried out in the NPT ensemble; the temperature of the system was coupled to 300 K by means of the Nose´-Hoover thermostat,48,49 and all dimensions of the simulation cell were coupled independently (anisotropic scaling) to reference pressures of 1 bar using the Parrinello-Rahman scheme.50,51 For the (0001) surface, we applied bond constraints to the silanol O-H bonds using the LINCS algorithm52 as implemented in Gromacs. This prevented unfavorable contacts from occurring between adjacent hydroxyl groups on this more dense surface. 2.3. Free Energy Calculations. The adsorption of water and of model hydrophilic and hydrophobic amino acids on the (101j0) surface was investigated by calculating the free energy profiles of the molecules as a function of their distance to the surface. We used the potential of mean constraint force method, which enables one to calculate the free energy change of any given process as a function of some reaction coordinate. In this case, our reaction coordinate, z, was simply the distance between the center of mass of the R-quartz slab and the center of mass of the water or amino acid analogue. We performed 30-40 MD runs per molecule, where the distance z was constrained to a different value in each run. Each run was 4 ns, which included 1 ns of equilibration. This corresponds to moving the molecule from bulk water onto the R-quartz surface. The Gromacs pull code42 was used to keep z constant. By integrating the resulting constraint force, λ, over all values of z, the position-dependent PMF of free energy change, F(z), is obtained:

F(z) )

∫ 〈λ〉z dz

(1)

2.4. Analysis. The majority of the analyses were carried out using the tools available within Gromacs.42 To obtain the density profiles (45) MacKerell, A. D.; Brooks, B.; Brooks, C. L.; Milsson, L.; Roux, B.; Won, Y.; Karplus, M. Encyclopedia of Compuational Chemistry; John Wiley & Sons: Chichester, U.K., 1998; Vol. 1, p 271. (46) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (47) Neria, E.; Fischer, S.; Karplus, M. J. Chem. Phys. 1996, 105, 1902. (48) Nose´, S. Mol. Phys. 1984, 52, 255. (49) Hoover, W. G. Phys. ReV. A 1985, 31, 1695. (50) Parrinello, M.; Rahman, A. J. Appl. Phys. 1981, 52, 7182. (51) Nose´, S.; Klein, M. L. Mol. Phys. 1983, 50, 1055. (52) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. J. Comput. Chem. 1997, 18, 1463.

of the components of the system as a function of the distance to the surface, the simulation cell was divided into slices in the z direction and the number of silica or water atoms in each slice was averaged over time and normalized with respect to the volume of the slice. The lateral density of the first layer of water was analyzed by further dividing the first layer of water into a grid and averaging the number of water molecules present in each area of the grid over time. The orientational distribution of the water within a given layer, moving upward from the surface, was determined by dividing the simulation box into slices and then calculating 〈cos(θ)〉 for all water molecules in each slice over time, where θ is the angle between the dipole vector of water and the z axis of the box. In this analysis the orientation of water relative to the surface can vary from cos(θ) ) -1, where the dipole of water lies perpendicular to the surface and both H atoms of the water are pointing toward the surface, through cos(θ) ) 0, where the dipole vector lies parallel to the surface plane and one H atom points toward the surface and one points away, to cos(θ) ) +1, where both the dipole vector and the H atoms point away from the surface. Hydrogen bond probability density profiles were determined by averaging the number of quartz-quartz, quartz-water, and water-water hydrogen bonds in slices above the surface in the same way as above. To define hydrogen bonds in our system, we used geometric criteria53 whereby two water molecules are considered to be hydrogen bonded only if the inter-oxygen distance is less than 0.35 nm and simultaneously the hydrogen-oxygen (hydrogen bonded) distance is less than 0.245 nm and the H-O · · · O angle is less than 30°. We used the same criteria to define water-water, water-quartz, and quartz-quartz hydrogen bonds. The hydrogen bond probability density profiles were normalized with respect to the bulk water hydrogen bond density to account for the different areas of the surfaces. The dynamics of water-water (bulk), water-water (interface), and water-quartz hydrogen bonds were characterized by the hydrogen bond time correlation function, C(τ), which is defined as

C(τ) ) 〈si(t) si(t + τ)〉

(2)

where si(t) ) 1 when a particular hydrogen bond i, present at time t, also exists at time t + τ according to the definition of a hydrogen bond used above and si(t) ) 0 otherwise. Thus, the correlation function C(τ) describes the probability that a particular hydrogen bond is intact at time t, given that it was intact at time zero. As noted above, the hydrogen bond correlation function was calculated from a trajectory where coordinates were saved every 10 fs. Lateral diffusion coefficients were calculated from the lateral mean squared displacement of water using the Einstein relation. To calculate the lateral diffusion coefficient of a particular layer (interface or bulk), we only considered the diffusion of molecules that were in that layer at time t ) 0 ps and at t ) 10 ps. Although this approach might bias the results by selecting water molecules which have a lower initial velocity, and are therefore more likely to remain in a particular layer, it should be noted that we used the same approach to estimate the lateral diffusion in the bulk. Therefore, any differences we have observed between diffusion at the interface and diffusion in the bulk are real differences and not an artifact of the methodology.

3. Results 3.1. Snapshots of the Quartz Surfaces. Snapshots of each of the surfaces studied here are presented in Figure 1. Each slab presents a different pattern of hydroxyl groups on the surface, in terms of both the lateral arrangement of the hydroxyl groups and the extent to which the hydroxyl groups protrude from the surface. The consequences of this are that each surface presents a different physical environment to the water, which influences the adsorption of water onto each surface (see below). The silanol groups on the (101j0) and (0001) surfaces are geminal silanols, (53) Jedlovszky, P.; Brodholt, J. P.; Bruni, F.; Ricci, M. A.; Vallauri, R. J. Chem. Phys. 1998, 108, 8528.

H2O and Amino Acid Analogue Interactions with Quartz

Figure 1. Side-on and plan views of the (101j0), (0001), and (011j1) surfaces. Silicon atoms are shown in yellow, oxygen atoms in red, and hydrogen atoms in white. The snapshots were generated using VMD.43

Figure 2. Density profiles of water atoms and surface silanol groups on each of the quartz surfaces as a function of the distance perpendicular to the surface (z distance).

while the (011j1) surface is terminated with vicinal silanol groups. It can be seen from the side-on views of the surfaces that the (101j0) and (0001) surfaces have pronounced indentations, which are due to the Si-O-Si triplets in the crystal structure lying along a direction out of the surface plane. The (011j1) surface, on the other hand, lacks these pronounced indentations as the uppermost layer of Si-O-Si atoms lies in approximately the same plane as the surface. 3.2. Structure of Water Perpendicular to the Surface. The underlying crystal structure of each surface yields a particular arrangement of quartz atoms at the interface, which is expected to influence the molecular arrangements of water molecules in the direction perpendicular to the surface. Density profiles of

Langmuir, Vol. 25, No. 3, 2009 1641

each surface, which show the density of quartz atoms and water as a function of the distance to the surface, are presented in Figure 2. The zero point of reference for each surface is the first peak in the Si atom density distribution. On the (011j1) surface, where the uppermost Si atoms are not perfectly in line with each other, the mean Si atom position is taken as the zero point of reference. First, we consider the density profile on the (101j0) surface. Beginning at z ) 0, as we move along the z axis the water density increases to a small (relative to the bulk density) peak at 0.14 nm, which corresponds to water molecules occasionally penetrating into the void spaces between the hydroxyl groups on the surface. At z ) 0.36 nm there is a peak in the water density which corresponds to the first layer of adsorbed water. Note that the water density of this layer is greater than that of the bulk solvent, which indicates that the water undergoes structural changes at the interface. We then observe a depletion region at z ) 0.47 nm where the water density is reduced, followed by an additional layer of structured water at z ) 0.69 nm. After the end of the third peak at z ) 0.80 nm there are some weak features in the water density, before it reaches the bulk value at approximately z ) 1.12 nm. The water density on the (0001) surface shows a similar variation with the z distance as that on the (101j0) surface. There are two noteworthy variations, however. First, there is no peak at z ) 0.14 nm as water molecules do not penetrate into the cavities on the surface, and second, the maxima and minima in the water density are more pronounced than on the (101j0) surface. Both of these differences may be explained by the smaller hydroxyl-hydroxyl distance on the (0001) surface compared with the (101j0) surface: more closely spaced hydroxyl groups may prevent water molecules from penetrating into the interstices and are likely to have a greater ordering effect on the interfacial water. As with the (101j0) surface, the water density appears to reach the bulk density by z ) 1.12 nm. On the (011j1) surface, there is a peak in the water density at z ) 0.36 nm, due to the first layer of adsorbed water. The density of adsorbed water on the (011j1) surface is significantly higher than on the other surfaces. The water molecules in the first layer also extend beyond the hydroxyl hydrogen atoms and approach the first layer of Si atoms more closely. As well as interactions between water and surface hydroxyl groups, the increase in water density on the (011j1) surface compared with the (101j0) and (0001) surfaces may be due to additional interactions between water and the bridging Si-O-Si oxygen atoms, which lie in the plane of the surface. Another possibility is the presence of vicinal silanols on the (011j1) surface rather than geminal silanols. Geminal silanols show some silanol-silanol hydrogen-bonding interactions (see below) which are in competition with the water-silanol hydrogen bonds. The second peak in the water density also appears closer to the surface than for the (101j0) and (0001) systems. Note that the water density profiles for the (101j0)

Figure 3. Lateral density profiles of the first layer of structured water on the (101j0), (0001), and (011j1) surfaces.

1642 Langmuir, Vol. 25, No. 3, 2009

Figure 4. Orientation of water with respect to the surface normal for each of the R-quartz surfaces, where θ is the angle between the dipole vector of water and the z axis of the simulation box.

and (011j1) surfaces agree with those presented by Lopes et al.41 using the same force field. 3.3. Water Ordering in the Plane of the Surface. The lateral density of the first layer of structured water on each of the surfaces is shown in Figure 3. We observe that water is structured differently on each surface depending on the positions of the hydroxyl groups. On the (101j0) surface, water forms striations that run diagonally across the surface, along the direction of the [011] vector. Water on the (0001) surface is arranged in rings that correspond to the underlying Si-O rings of the surface (see Figure 2). It is clear that the positions of the water molecules are dominated by the positions of the silanol hydroxyl groups but that they are also influenced by the underlying layers of the crystal. We find that the most pronounced lateral ordering occurs in the (011j1) surface. This observation correlates with the higher water density for the first layer on the (011j1) surface compared to the (101j0) and (0001) surfaces (Figure 2). As discussed above, a contributing factor may be strong interactions between the bridging Si-O-Si oxygen atoms, which are present in the plane of the surface. The lateral patterning on the (011j1) surface supports this idea, as it is clear from Figure 3 that the regions of high density correspond to the positions of in-plane bridging oxygen atoms on the surface as well as the surface hydroxyl groups. 3.4. Water Orientation with Respect to the Surface Normal. To gain insight into the alignment of water molecules on the surface, we calculated the orientation of water as a function of the distance to the surface, as shown in Figure 4. Here, the z ) 0 nm baseline is in the same position as in the density profiles (see above). For all three surfaces, the water shows anisotropic orientational ordering in line with the water density. That is, on average, the adsorbed water in each layer has a tendency to be aligned in a particular way relative to the surface. On the (101j0) surface there is a small peak close to zero, which corresponds to those water molecules that penetrate the surface hydroxyl groups. These water molecules are largely restricted to an upright position with one H atom pointing toward the surface and one H atom pointing toward the bulk water, due to spatial restraints in this region. The first peak at z ) 0.23 nm corresponds to water molecules located in the area of low water density adjacent to the surface. The dip at z ) 0.4 nm is at the position of the first layer of adsorbed water. Water in this layer tends to be arranged with one hydrogen atom pointing toward the surface where it is able to interact with the surface silanol groups. This weak orientational ordering of water tails toward zero with an increase in distance to the surface. Other than the first peak, the orientational ordering on water on the (0001) surfaces is essentially the same as on the (101j0) surface. On the (011j1) surface, there is a minimum at z ) 0.36 nm which is associated with water molecules oriented with a hydrogen atom pointing down toward the surface and there are some weak features extending out to about z ) 1.2 nm,

Notman and Walsh

Figure 5. Probability density of water-water (ww), quartz-water (qw), and quartz-quartz (qq) hydrogen bonds as a function of the distance to the R-quartz surfaces. The zero point of reference is the position of the topmost layer of Si atoms. The (101j0) surface has a small peak of qq hydrogen bonds, which is not visible on this scale, whereas the 011 surface has an insignificant number of qq hydrogen bonds.

indicating that there is some orientational ordering associated with each of the structured water layers. 3.5. Hydrogen Bond Structure and Dynamics. Silanol hydroxyl groups are known to form strong hydrogen bond complexes with water37 and biological molecules.54 It follows that hydrogen bonding with the surface is likely to be a driving force for adsorption of water and other molecules and will influence the structure of water close to the surface. The probability densities of water-water and water-quartz hydrogen bonds as a function of the distance from the surface are given in Figure 5. For the (0001) surface, quartz-quartz hydrogen bonds are also shown. The number of quartz-quartz hydrogen bonds for the (101j0) and (011j1) surfaces were also calculated but found to be insignificant. In the same way as above, the zero point of reference is the position of the topmost layer of Si atoms. As predicted, all three surfaces form hydrogen bonds with water at the interface. In addition, there are peaks in the hydrogen bond probability density, which correlate with the structuring of water in the direction normal to the surface and the orientation ordering of water in each of the layers. This suggests that the structured water layers are stabilized by hydrogen bonds both directly with the surface and with themselves. An understanding of the dynamics of hydrogen bonding is also important for elucidating the role that hydrogen bonds may play at an interface, as the lifetime of a hydrogen bond is an indicator of its stability and strength.55 To probe the hydrogen bond dynamics on quartz, we calculated the hydrogen bond time correlation function for water-water hydrogen bonds in the bulk water, water-water hydrogen bonds at the interface, and quartz-water hydrogen bonds, as shown in Figure 6. In general, for all three surfaces, the correlation function decays fastest for the hydrogen bonds in bulk water, followed by the water-water hydrogen bonds at the interface and then the quartz-water hydrogen bonds. This is evidence that the hydrogen bonds formed at the interface are stronger than those formed in the bulk and supports the idea that hydrogen bond formation at the interface plays a key role in the adsorption of water and hydrophilic solutes onto quartz surfaces. 3.6. Lateral Diffusion Coefficients of the Water Layers. To ascertain the effects of the quartz surfaces on the mobility of water, the lateral diffusion coefficients of water molecules in the first layer of structured water were calculated and compared with that of bulk water. It is clearly shown in Table 1 that the lateral diffusion of structured water on all three surfaces is reduced to about 50% of the lateral diffusion in the bulk. As discussed (54) Murashov, V. V.; Leszczynski, J. J. Phys. Chem. A 1999, 103, 1228. (55) Chanda, J.; Bandyopadhyay, S. J. Phys. Chem. B 2006, 110, 23443.

H2O and Amino Acid Analogue Interactions with Quartz

Figure 6. Hydrogen bond time correlation function for the hydrogen bonds formed between the water molecules, Cww(t), in the interface and bulk regions and hydrogen bonds between the quartz surfaces and water, Cqw(t). Table 1. Lateral Diffusion Coefficients, Dlateral (1 × 10-5 cm2 s-1), of the First Monolayer of Water on Each of the Surfaces and of Bulk Water system

Dlateral

(101j0) (0001) (011j1) bulk

2.32 ( 0.42 1.65 ( 0.43 2.18 ( 0.44 4.98 ( 1.57

above, a hydrogen bond network between water and quartz is formed at the interface where the hydrogen bonds have a lifetime greater than that of hydrogen bonds in bulk water. This network of hydrogen bonds provides an explanation for the reduced movement of water molecules on the surfaces. In contrast to our results, Adeagbo35 did not observe a slowing of water diffusion due to hydrogen-bonding interactions between the interfacial water molecules and the hydroxyl groups on the surface; however, this might be due to the relatively short time scales sampled by their simulations (4 ps). 3.7. Free Energy Profiles. We have calculated the free energy change of moving a water, methane, and methanol molecule from bulk water onto the (101j0) surface. The free energy profiles are presented in Figure 7. The free energy profile of water as a function of the distance to the (101j0) surface is in close agreement with the water density profile for the same surface (see Figure 2). That is, the minima in the free energy profile correspond to maxima in the water density. This relationship is a good test of convergence in our PMF calculations as the free energy difference at a particular height can also be determined from the water density perpendicular to the surface.56 The free energy profile shows that there is a minimum in the free energy, which corresponds to the first layer of structured water on the surface. Therefore, the results suggest that the loss in conformational entropy that occurs when water molecules move from the bulk phase to the interface is compensated by interactions with the surface. On the basis of the analysis above, we propose that the driving force for the interaction is the formation of additional water-quartz and water-water hydrogen bonds per water molecule at the interface compared with hydrogen bonds in the bulk. As the water molecule approaches the hydroxyl groups, there is a small (4 kJ mol- 1) free energy barrier to cross the structured water, followed by a local minimum between the interstices on the surface created by the surface silanol groups. Visual inspection of the simulation trajectory reveals that the constrained water molecule tends to be oriented with one hydrogen atom pointing down toward the surface. In this position, the water is able to form hydrogen bonds with the surface silanol groups and also with the oxygen (56) Marrink, S. J.; Berendsen, H. J. C. J. Phys. Chem. 1994, 98, 4155.

Langmuir, Vol. 25, No. 3, 2009 1643

Figure 7. PMF free energy curves of water, methanol, and methane adsorption onto the quartz (101j0) surface, where z is the distance between the molecule center of mass and the basal Si atoms.

atoms in the second atomic layer of the quartz slab. As expected, the free energy then shows a steep increase as the water molecule is pushed toward the surface. To gain insight into the adsorption of polar amino acid residues onto quartz, we calculated the free energy profile of methanol, which is an analogue of the side chain of serine. The free energy curve for methanol has two significant features. First, as we move toward the surface from the bulk water, there is a local minimum at 0.39 nm which corresponds to the position of the structured water layer. There is then a small free energy increase of about 1 kJ mol-1 before we reach another, lower, free energy minimum at approximately 0.18 nm above the surface. This corresponds to the methanol center of mass sitting just above the mean position of the surface hydroxyl hydrogens. Within the center of mass constraint, the methanol molecule is free to rotate. To further our understanding of the free energy profile of methanol, we monitored the angle made between the C · · · O vector of methane and the surface normal for the simulations of methanol at 0.18 nm, 0.39, and 0.96 nm (position in bulk water). For the simulations at 0.39 and 0.96 nm, the methanol molecule adopts the full range of orientations, with the most frequent orientation being with the C · · · O vector lying parallel to the surface; we did not observe any significant differences in the orientation of methanol between the two simulations. For the simulation at 0.18 nm, we found that methanol flips between two main orientations. The position observed most frequently is one where the methanol sits in the gap between four surface silanol groups (see Figure 3) with the three methyl hydrogens located just beneath the surface silanols and with the methanol hydroxyl group pointing toward the solvent. In this position the hydrophobic part of the molecule is shielded from the water and the hydroxyl group is available for hydrogen bonding with water, both of which are likely to be key factors for stabilization of this conformation. Less frequently, methanol flips to a conformation where the methanol hydroxyl group is pointing downward, forming a hydrogen bond with a surface silanol group, and the methyl groups are partly exposed to the water and partly interacting with the surface. Multiple modes of binding are likely to make an important contribution to the binding affinity of particular amino acids to specific surfaces. The free energy profile of methane is similar to that of methanol. There is a favorable net free energy change of about 6 kJ mol-1 to move a methane molecule from the bulk water to the quartz surface. As the distance to the surface decreases, there is a minimum at 0.4 nm, which corresponds to the position of the methane molecule in the structured water layer. We propose that this is due to the effect of the surface shielding the methane from the water as the methane approaches the surface. As the water is fairly weakly structured, we assume that the penalty for methane to be in this region is low compared with the advantage of the

1644 Langmuir, Vol. 25, No. 3, 2009

shielding effect. At 0.3 nm there is a local maximum in the free energy, and we associate this with the methane molecule making unfavorable contacts with the surface hydroxyl groups as it approaches them more closely. After this point there is a decrease in the free energy as the methane is moved into the interstices on the surface. The minimum free energy occurs when the center of mass of the methane molecule is at 0.14 nm above the surface Si atoms, slightly lower than the methanol molecule due to the absence of a hydroxyl group on methane. Finally, as the methane approaches the surface even closer, the free energy curve increases steeply.

4. Discussion Our results show that there are at least two layers of structured water on the surface of quartz, with some weak structuring that extends up to approximately 1.2 nm from the surface. Water in contact with the surface is structured laterally and orients in a way which maximizes the formation of hydrogen bonds with the surface hydroxyl groups. The simulations suggest that hydrogen bond formation is a key driving force for water structuring on hydroxylated quartz surfaces. We observe that the interfacial hydrogen bonds have a longer lifetime, and hence are more stable than hydrogen bonds in the bulk water. This is in good agreement with previous MD simulations, which have also provided evidence for the formation of strong hydrogen bonds between water and silica surfaces,37,41 and our observations are supported by experimental studies that suggest that there is a monolayer of ordered water on amorphous silica.25,27-29,57 Structured water layers are a general feature of solid-liquid interfaces.1 In the case of quartz we found that the water is rather weakly ordered; for example, the orientation of water relative to the plane of the bilayer was always close to a random orientation. In addition, the lateral diffusion coefficients of the first layer of water are lower than in bulk water, yet the water molecules still have some degree of mobility. This is in contrast to other well-studied hydrophilic surfaces; for example, the water on the surface of titania forms a highly structured icelike layer with strong orientational ordering.58 This has implications when one considers how biomolecules might bind to surfaces. Weakly structured water is probably easier to displace than strongly structured water for direct interactions with the surface, while the effects of strongly adsorbed water are likely to extend further into the bulk water and have a greater influence on the local environment of bound molecules. The lateral ordering of water on quartz was found to be strongly dependent on the crystallographic orientation of the surface, as the positions of water molecules are closely correlated to the positions of the surface silanol groups. This patterning of water on the surface is likely to affect the diffusion of water and other molecules on the surface. For example, on each surface there are channels of water density along which diffusion would be expected to occur. In addition, the orientation of these water molecules, together with the underlying crystal structure, determines the groups that are presented for interaction with other compounds. It has been shown that anisotropic mobility affects the aggregation and assembly of molecules on surfaces and interfaces59,60 and regions of high or low water density may bring about anisotropic diffusion on quartz. (57) Tarasevich, Y. I. Colloids J. 2007, 69, 212. (58) Skelton, A. A.; Walsh, T. R. Mol. Simul. 2007, 33, 379. (59) D’angelo, M.; Y. Aristov, V.; Soukiassian, P. Phys. ReV. B 2007, 76, 045304. (60) Roos, K. R.; Roos, K. L.; Lohmar, I.; Wall, D.; Krug, J.; von Hoegen, M. H.; zu Heringdorf, F.-J. M. Phys. ReV. Lett. 2008, 100, 016103.

Notman and Walsh

The free energy profiles of water, methanol, and methane as a function of the distance from the surface have provided some useful insights into the mechanisms by which polar and nonpolar amino acids bind to quartz surfaces. Interestingly, we found that methane and the hydrophobic methyl group of methanol interact directly with the surface, and we propose that this is due to the effect of the surface shielding the hydrophobic moieties from the water. This effect may be limited to small hydrophobic moieties that can spatially fit into the voids between the surface silanol groups. A recent study by Ghiringhelli et al.61 showed that methyl groups hinder the adsorption of amino acids onto the Pt(111) surface; however, the surface charge distribution and topography of the Pt(111) surface are remarkably different from those of the quartz (101j0) surface; in particular there are no suitable interstices on the Pt(111) surface to protect the methyl groups. Our surprising result that hydrophobic moieties bind to quartz concurs with previous studies of quartz-binding peptides that were found to contain a high proportion of nonpolar, neutral residues.22 As well as hydrophobic groups, we also observed that water penetrated between the interstices on the (101j0) surface, although the mechanisms driving this are different. This feature of the (101j0) surface is in contrast to the (0001) and (011j1) surfaces. It will be particularly interesting to make a comparison with the (0001) surface which is terminated with geminal silanol groups, yet has a significantly larger water-depleted region between the first layer of structured water and the silanol groups. Another difference between the surfaces, which may have implications for the binding of larger biomolecules, is the density of the hydroxyl groups. One would expect, for example, the surface coverage of peptides that interact with the silanol groups to be lower on the (011j 1) surface than the (101j 0) surface. Indeed, we hope to address these issues in our ongoing work on peptides interacting with these quartz surfaces.

5. Conclusions We have carried out MD simulations to investigate the interactions of water with the (101j0), (0001), and (011j1) surfaces of R-quartz. In addition, we have calculated the free energy change to adsorb water and a model hydrophobic and hydrophilic amino acid onto the (101j0) surface. Our simulations suggest that there are at least two layers of ordered water on the quartz surfaces. We propose that the driving force for the interaction of water with quartz surfaces is the formation of a network of hydrogen bonds between water and silanol groups and also between adjacent water molecules on the surface. The water layers exhibit a lateral patterning on the surface which is expected to have an effect on the binding and aggregation of biomolecules on quartz and other silica surfaces. Our free energy calculations show that the adsorption of small hydrophobic moieties is favorable on the (101j0) surface as these groups are able to penetrate into small voids on the surface where they are shielded from interfacial water. Acknowledgment. We acknowledge the computer facilities of the Centre for Scientific Computing, University of Warwick, and the use of the U.K. National Grid Service in carrying out this work. The funding for this project is provided by the EPSRC (Grant EP/E02095X/1). LA803324X (61) Ghiringhelli, L. M.; Hess, B.; van den Vegt, N. F. A.; Site, L. D. J. Am. Chem. Soc. 2008, 130, 13460.