This is an open access article published under an ACS AuthorChoice License, which permits copying and redistribution of the article or any adaptations for non-commercial purposes.
Feature Article pubs.acs.org/JPCB
Molecular Modeling of Water Interfaces: From Molecular Spectroscopy to Thermodynamics Yuki Nagata,*,† Tatsuhiko Ohto,‡ Ellen H. G. Backus,† and Mischa Bonn*,† †
Max Planck Institute for Polymer Research, Ackermannweg 10, 55128 Mainz, Germany Graduate School of Engineering Science, Osaka University, 1-3 Machikaneyama, Toyonaka, Osaka 560-8531, Japan
‡
ABSTRACT: Understanding aqueous interfaces at the molecular level is not only fundamentally important, but also highly relevant for a variety of disciplines. For instance, electrode−water interfaces are relevant for electrochemistry, as are mineral−water interfaces for geochemistry and air−water interfaces for environmental chemistry; water−lipid interfaces constitute the boundaries of the cell membrane, and are thus relevant for biochemistry. One of the major challenges in these fields is to link macroscopic properties such as interfacial reactivity, solubility, and permeability as well as macroscopic thermodynamic and spectroscopic observables to the structure, structural changes, and dynamics of molecules at these interfaces. Simulations, by themselves, or in conjunction with appropriate experiments, can provide such molecular-level insights into aqueous interfaces. In this contribution, we review the current state-of-the-art of three levels of molecular dynamics (MD) simulation: ab initio, force field, and coarse-grained. We discuss the advantages, the potential, and the limitations of each approach for studying aqueous interfaces, by assessing computations of the sum-frequency generation spectra and surface tension. The comparison of experimental and simulation data provides information on the challenges of future MD simulations, such as improving the force field models and the van der Waals corrections in ab initio MD simulations. Once good agreement between experimental observables and simulation can be established, the simulation can be used to provide insights into the processes at a level of detail that is generally inaccessible to experiments. As an example we discuss the mechanism of the evaporation of water. We finish by presenting an outlook outlining four future challenges for molecular dynamics simulations of aqueous interfacial systems.
I. INTRODUCTION The history of molecular dynamics (MD) simulations of liquid water dates back to 1971.1 Since the first MD simulation by Rahman and Stillinger,1 MD simulations of bulk water have been extensively performed for understanding the hydrogen bond dynamics of water,2 the vibrational signatures of water,3−5 nucleation of ice,6 and excess protons in water,7 to name a few topics. In parallel with the increasing applications of water MD simulations, the technique of MD simulations has been improved in several aspects. Increasingly accurate force field (FF) models of water have been developed, as have simulation techniques to include nuclear quantum effects (highly relevant, given the abundance of light hydrogen atoms in water) for monitoring realtime dynamics.8 Moreover, MD techniques based on the electronic structure theory9−12 have been developed and employed for studying liquid water structure and dynamics. Simulations of water in confined environments and at aqueous interfaces have also been increasingly reported. The application of MD simulations to water at interfaces started in the 1980s,13 although surface tension was the only quantity for which the simulation could be compared with the experimental data. The surface tension values computed with the popular water models agree well with the experimental data.14 Furthermore, it was © 2016 American Chemical Society
demonstrated that the long-range parts of the Lennard-Jones and electrostatic interactions also largely contribute to the surface tension (∼25% contribution for the intermolecular interactions with >6 Å distance).15 Computing the surface tension of water thus provides a compelling means to examine the accuracy of the long-range intermolecular interactions at aqueous interfaces. In turn, short-range intermolecular interactions govern the hydrogen bond network of water. Vibrational spectroscopic techniques probing the O−H stretch mode have been frequently used for examining the hydrogen bond network, because of the high sensitivity of the O−H stretch mode to the hydrogen bond strength. Since the early 1990s the development of the surfacespecific vibrational sum-frequency generation (SFG) technique and its application to the water surface16 have opened up new possibilities to experimentally probe, in a highly selective manner, the vibrational response of the interfacial water molecules and exclude the contribution from the bulk water molecules.17,18 By simulation of the SFG spectra of water and comparison of the simulated spectra with the experimental data, Received: January 29, 2016 Revised: March 19, 2016 Published: March 24, 2016 3785
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
cm−1 to be of different molecular origins, thus providing conflicting interpretations of the spectra.22,32,44,45 We will discuss here recent efforts employing AIMD simulations, which rely less on assumptions than force field models do, to describe the SFG response of interfacial water. Experimentally, before performing SFG measurements, the surface tension is routinely measured to verify the cleanliness of the sample. The surface tension is a useful quantity in this respect, because it is very sensitive to contamination at the surface as well as to the presence of solutes. Likewise, the theoretical calculation of the macroscopic surface tension from the simulation cell with a limited size is also a good test of the simulation protocol, but also a very challenging exercise. Here, we introduce the recent efforts in calculating the surface tensions from AIMD trajectories.15 II.A. Ab Initio, Force Field, and Coarse-Grained MD Approaches. In this subsection, we introduce the three categories of MD simulations. The first category is the AIMD simulation, for which the interatomic forces are calculated on the basis of electronic structure theory.9,46,47 Provided that the forces are calculated from high-level theory, this technique constitutes the most accurate MD simulation. However, for practical purposes, the number of basis sets must be limited; a generalized gradient approximation for the exchange-correlation (XC) functional must be employed, and van der Waals (vdW) interactions must be approximated as corrections.48 These different approximations can affect the calculation results very critically, and AIMD simulations therefore do not necessarily provide an accurate description. Nevertheless, AIMD has been frequently used to simulate the water−air interface,49−53 because AIMD simulation is a parameter-free technique. A typical time step for integrating equations of motion is 0.3−0.4 fs. The second category is that of FFMD simulations. In this type of simulation, the intermolecular interactions are calculated on the basis of the electrostatic interaction between the (partial) charges on the atoms and Lennard-Jones potentials. The FF parameters are often optimized to reproduce the bulk properties. The electrostatic interactions are calculated with the Ewald or related techniques, 54 where the short- and long-range interactions are calculated in real- and k-space, respectively. Once an FF model has been obtained, the computational cost of an FFMD simulation is much lower than that of an AIMD simulation. The maximum time step can be as large as 2−4 fs with the RESPA algorithm.55 Recently, the FF models of water have been developed for describing the properties at the interface more accurately.22,56,57 Moreover, ab initio based FF models have been applied for simulating interfacial water.24,25,58−62 The third category is that of coarse-grained (CG) MD techniques. Although different routes for CGMD simulation have been presented, the concept underlying CG approaches is common: beads representing a collection of atoms are used to reduce the number of particles in the simulation and enlarge the time step. This technique has been increasingly used for understanding the structure and dynamics of larger systems with N > 100 000 such as biological membranes and proteins, where N is the number of atoms. A typical time step used in CGMD simulations is a few tens of femtoseconds. The number of atoms in the simulation cell and the total simulation time which these different MD simulation techniques can deal with are displayed in Figure 1. One can see that, from AIMD to FFMD, and from FFMD to CGMD, the accessible simulation time and number of atoms are dramatically increased. However, these techniques have individual intrinsic challenges.
the structure and dynamics of water in the MD simulation can be connected with the SFG signals. As such, MD simulation can be used for correlating SFG signals to the molecular-level picture.17,19−25 Most of the SFG simulations have been conducted with FF models.17,19−25 The recently developed ab initio MD (AIMD) technique for computing the SFG spectra26−29 allows access to the complex water interface without parameters, because the forces in an AIMD simulation are calculated directly from the electronic structure theory. Despite the substantial progress in both experimental and simulation approaches for the study of aqueous interfaces, simulation and experimental results are not always wellconnected. For instance, the following examples apply: (1) The ∼3100 cm−1 O−H stretch band of water in the imaginary part of the vibrational SFG response for the water−air interface22,25,30−32 and its assignments23,27,33,34 are still under debate. (2) The surface tension of water at the water−air interface has been simulated only very recently with AIMD techniques,15 while the surface tensions of the other aqueous interfaces have not been studied. (3) The molecular mechanisms of the variation of the surface tensions caused by addition of the complex ions/surfactants/lipids and their mixture in water are not clear yet. Unveiling these effects is quite important for engineering the lipid/surfactant and controlling the surface property of aqueous surfaces in the presence of lipids or surfactants. In addition, simulations can be used for predicting the molecular mechanism of the physical phenomena at interfaces. For instance, MD simulations at the water−air interface have been used for predicting the molecular mechanism behind water evaporation at the water−air interface35−37 and the behavior of the water droplet at the surface,38−40 as well as understanding the molecular mechanism of surface phenomena such as the fusion mechanism of the lipid bilayer41,42 and ion transport mechanism in membranes.43 Such simulations motivate new experiments, and, thus, MD simulations can very nicely complement experimental approaches. In this article, we introduce simulation techniques aimed at computing the experimentally accessible observables such as the surface tensions and SFG responses and review their recent applications to the water−air and water−lipid/surfactant interfaces. New insights into water evaporation from the simulation are also presented. Subsequently, we will list some of the challenging topics which simulations should urgently address for a better understanding of experimental data. The organization of this article is as follows. In section II, we introduce the recent progress in SFG and surface tension calculations, as well as the prediction of the molecular mechanism of water evaporation. Relevant simulation techniques and algorithms are also introduced. In section III, four challenges for MD simulations are given by referring to recent experimental findings. Section IV is devoted to the conclusions.
II. MOLECULAR DYNAMICS SIMULATION AT AQUEOUS INTERFACES Since MD simulations can visualize the time-dependent molecular structures in silico, these simulations allow us to connect the molecular structures with the experimental observables. Therefore, the analysis of MD simulations has been developed as a versatile tool to interpret experimental data and gain molecular-level insights into the interfacial molecular organization. On the other hand, previous papers using different FF models have assigned the peak in the SFG response at ∼3100 3786
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
The values for the surface tension of water at the water−air interface obtained with FFMD and CGMD simulations are summarized in the top panel of Figure 2. Except the SPC model,
Figure 1. Number of water molecules and the length of total trajectories used in water−air interface simulations. The filled points are the data used in refs 52, 63, and 64. Here, note that a MARTINI CG water bead contains four water molecules. The boundaries between AIMD and FFMD and between FFMD and CGMD are determined by the fact that the computation time increases with a factor of ∼N3 and N ln N for density functional theory calculation and the particle mesh Ewald method,65 respectively.
In the subsequent sections, we will highlight the strengths and weaknesses of these different MD approaches, by considering the surface tension and vibrational SFG response of water at aqueous interfaces. II.B. Surface Tension of Water. The surface tension of water can be accessed via the stress tensor Figure 2. (Top) Surface tension values simulated with various FF models of water together with the MARTINI CG model. The surface tensions with the SPC/Fw and MARTINI models are from refs 15 and 69, respectively, while the other simulation data are from ref 14. (Bottom) Surface tension values simulated with AIMD trajectories. The error bars represent the 95% confidential intervals. The data are reproduced from ref 15.
∞
γ=
∫−∞ dz(pN (z) − pT (z))
(1)
where the z-direction lies along the surface normal, and pN(z) and pT(z) represent the normal and tangential components of the pressure tensor, respectively. In contrast to the mechanical route for computing the surface tension, the surface tension can be also calculated using thermodynamic considerations. A robust method using the thermodynamic relation is the test area method, for which the surface tension is given by66 γ=−
⎛ ΔU + ⎞ kT ⎛ ⎜⎜ln exp⎜ − ⎟ ⎝ kT ⎠ 2ΔA ⎝
⎛ ΔU − ⎞ ⎟ − ln exp⎜ − ⎝ kT ⎠
all tested FF water models provide less than ∼15% deviation from the experimental result of 72 mN/m.14,15 In contrast, the MARTINI CG water model significantly underestimates the surface tension.64,69 Furthermore, the trend of the surface tension with the simulation cell size is opposite for the MARTINI CG and FF models; in the MARTINI model, a larger cell size results in a lower surface tension value, while the FF models with larger cell size provide larger surface tensions (results not shown). This difference presumably arise from the lack of longrange interactions in the CG models,64 for which only the shortranged interactions (∼10 Å) are considered, unlike the FF model. As such, although the surface tension of the water− dodecane interface with the MARTINI model exhibits good agreement with the experimental data,64 this agreement may be in part fortuitous, as the absence of long-ranged interactions is compensated by overly strong short-ranged interactions. More detailed analyses of the surface tension calculations with the CGMD trajectories would be required to verify this statement. Subsequently, we turn to AIMD simulation. The surface tensions computed from the AIMD trajectories with different levels of theory are scattered from 30 to 170 mN/m (see the bottom panel of Figure 2).15 A small basis set such as the doubleξ basis set plus the polarization function (DZVP) causes substantial Pulay stress due to incompleteness of the basis sets,70 resulting in a very large value for the computed surface tension. Thus, the triple-ξ or more accurate basis set would be required
⎞ ⎟⎟ ⎠ (2)
where ΔU ± represents the conformational energy difference by scaling the surface area from A = 2LxLy to A ± ΔA, and ⟨B⟩ denotes the thermal average of quantity B. In both ways, the surface tension can be calculated from the MD trajectories of the simulation cells with periodic boundary conditions. Care should be taken, however, as the periodic cell limits the size of the interface, causing artificial finite size effects. Finite size effects on the surface tension can be accounted for by using either long-range Lennard-Jones tail corrections14 or the extrapolation method by making use of the ln(L z )/L x 2 dependence.67,68 Our MD simulation at the water−air interface shows that these two approaches provide similar corrections for the finite size effects on the surface tension at the water−air interface.15 However, this does not necessarily mean that these two different concepts lead to similar corrections at other aqueous interfaces such as the water−lipid interface and water− air interface with ions. Further studies should address the importance of finite size effects in calculations of the surface tension. 3787
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
spectra. We have recently presented an efficient calculation algorithm for computing the SFG spectra of the water O−H stretch mode based on the surface-specific velocity−velocity correlation function, by separating degrees of freedom of the nuclei from solvation effects such as the induced dipole and polarizability.26 The SFG response can be expressed as
for AIMD simulation of the surface tension of water. Furthermore, a comparison of the van der Waals (vdW) corrections48 indicates that the choice of the vdW correction critically affects the surface tension values; without the vdW corrections, the calculated surface tension value is a factor of 3 lower than the value with the Grimme’s D3 vdW correction.72 Moreover, the Grimme’s D2 vdW correction71 provides a 40% larger surface tension value than the D3 vdW correction.72 The D3 correction provides the surface tension value closest to the experimental data. In contrast, the value for the calculated surface tension is remarkably insensitive to the choice of the XC functionals among BLYP,73,74 PBE,75 and revPBE.76 Simulating the surface tension using the second-order Møller−Plesset perturbation theory level (MP2) water12,77,78 is also a challenge in the near future. II.C. SFG of Water and Vibrational Modes at the Water− Air Interface. The SFG signal is generated by infrared (IR) and visible (vis) pulses that overlap spatially and temporally at the water interface. Experimentally, the SFG intensity, ISFG, is recorded as a function of IR frequency, and scales as ISFG ∝ |χ(2)(ω)|2, with χ(2)(ω) the second-order nonlinear response function. χ(2)(ω) is a complex function, and its frequencydependent real and imaginary parts can be obtained from phaseresolved measurements.19,79−81 Of particular interest is the IR resonant part of Im[χ(2)(ω)], since that represents the dissipative part of the response function, and is directly comparable to IR absorption and Raman scattering in bulk. Theoretically, the abc component of the resonant part of the SFG signal can be computed from the dipole moment−polarizability (μ−α) time correlation function as (2) (ω) = χabc
1 kT
⎧ Q (ω)μ′(ω)α′(ω) ∞ −iωt dt e ⎪ 0 iω 2 ⎪ OH ⎪ r j⃗ ̇ (t ) · r j⃗OH(t ) (2) χabc (ω) = ⎨ × ∑ r OH ċ , i (0) ⎪ | r j⃗OH(t )| i,j ⎪ ⎪ ⎩0
∫
a=b
a≠b (4)
where μ′(ω) and α′(ω) denote the transition dipole moment and polarizability at frequency ω, respectively, and rOHi(t) represents the intramolecular O−H distance of O−H group i at time t.26 Q(ω) is the quantum correction factor.90 The parametrization of μ′(ω) and α′(ω) is given in ref 91. Figure 3
∞
∫−∞ ddt ⟨Mc(0)Aab(t )⟩e−iωt dt
(3) Figure 3. SFG spectra of the O−H stretch mode at the isotopically diluted water−air interface simulated from AIMD trajectories at different levels of theory. The spectra are normalized to the peak intensity of the dangling O−H peaks centered at ∼3700 cm−1. This figure is reprinted with permission from ref 26. Copyright 2015 American Institute of Physics.
where Mc(0) is the c-component of the dipole moment for the total system at time t = 0 and Aab(t) is the ab component of the polarizability for the total system at time t. So far, in the SFG simulation community, much attention has been paid to the water O−H stretch mode at the water−air interface, partially because of the simple modeling of the water−air interface.19,82−84 Most of the simulated (Im[χ(2)]) SFG spectra of the O−H stretch mode reported prior to 2008 exhibited a positive 3700 cm−1 peak and a negative 3300−3500 cm−1 band, while a positive band centered at ∼3100 cm−1 had not been reported.20,84,85 The presence of a positive 3100 cm−1 SFG band was experimentally reported using the phase-sensitive SFG approaches by several groups.30,81,86,87 Motivated by these experimental observations, several theoretical groups have tried to reproduce this positive 3100 cm−1 band. The first simulation reproducing this positive peak was reported by Morita’s group, where they modified the induced dipole moment calculation used in the time correlation function calculation (eq 3).32,44,88 Subsequently, Skinner’s group also reported the positive band22,45 by including the three-body interaction in the FF model. However, recent SFG measurements question the presence of this 3100 cm−1 SFG feature.31,89 These controversial discussions challenge both the unambiguity of experimental results and the quality of the simulations. To avoid bias in the computational results and interpretation originating from the choices of the details of FF models, applying a parameter-free AIMD simulation technique to the SFG calculation seems to be a promising route. However, the huge computational cost required for AIMD simulation has prohibited the widespread use of AIMD simulations for computing the SFG
displays the computed SFG spectra of the O−H stretch mode at the isotopically diluted (no inter/intra-molecular vibrational mode couplings) water−air interface from the AIMD trajectories at various levels of theory. There is some frequency offset for the PBE based functionals relative to the BLYP functionals, but for all functionals the AIMD data clearly do not show any positive ∼3100 cm−1 peak, consistent with the most recent SFG measurements.31,89 This implies that the AIMD simulation can address the SFG spectra with adequate accuracy. Note that, in this AIMD simulation, the nuclei were treated classically, and thus the O−H frequency would be red-shifted when properly accounting for nuclear quantum effects.92,93 Therefore, good agreement of the SFG spectra between the BLYP simulation and experiment would arise from the cancellation of the error due to the omission of the nuclear quantum effects and the inadequate description of the BLYP function.93 Above, we have focused only on the O−H stretch mode. However, the vibrational modes of water consist not only of the O−H stretch mode (3000−3700 cm−1), but also of the H−O−H bending mode (∼1650 cm−1), the librational mode (400−1000 cm−1), and the hydrogen bond stretch mode (∼200 cm−1). The SFG spectra for these frequency ranges provide another route to examine the interfacial water structure and the consistency 3788
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
water−PC and water−AO interfaces.27 The isotopic dilution serves to remove effects due to intra- and intermolecular coupling. Water at the PC interface is homogeneously oriented up, while at the AO interface it shows both up and down orientations. In both systems the negatively charged part of the headgroup dominates the orientation of water, because the negatively charged part allows a water molecule to interact with it via a hydrogen bond.99 The different orientations of water near AO and PC interfaces can be attributed to the fact that water molecules form a continuous hydrogen bond network from the bulk up to the carbonyl group of the PC molecule.27,100 This ability of the carbonyl group to terminate the hydrogen bonding network of water is essential for understanding the hydration structure of the lipid/surfactant molecules. In fact, the large difference of the surface dipole moment between the 1,2dimyristoyl-sn-glycero-3-phosphatidylcholine (DMPC) and dihexadecylphosphatidylcholine (DHPC)101 can be reasonably attributed to the hydration activity due to the carbonyl group. The precise location of the carbonyl groups in the surfactant headgroups is also expected to affect the dynamic surface tension.102 II.E. Water Evaporation Mechanism. Above, we introduced recent MD simulations of aqueous interfaces, with particular focus on the SFG response and surface tension calculations. MD simulation has proven to be a very powerful tool for interpreting and predicting SFG spectra. Beyond the thermally averaged picture of the water structure and dynamics at the aqueous interfaces, MD simulations allow access to rare events that are very challenging to probe spectroscopically. At aqueous interfaces, such rare events include ion and water transport across a lipid membrane, evaporation of water molecules, and fusion processes of lipid bilayers. MD simulations are able to provide rather unique mechanistic insights into the processes underlying these rare events. Hereafter, we shall briefly introduce the recent study on the water evaporation. It has been known that water molecules which evaporate from the water surface experience collision-like interactions with the interfacial water molecules, just before evaporating.35,103 A key question was how such a desorbing water molecule, which manages to “escape” the water surface, gains its driving force to break the strong hydrogen bond with that surface. To understand this, we traced, using MD trajectories, the time evolutions of the evaporating water molecule (A) and the water molecule (B) that it strongly interacts with, just prior to evaporating. The trajectories of these two molecules are plotted in Figure 5a. From this figure, it is apparent that the time evolution of the trajectory of molecule A at t > 0 is not symmetric with respect to the time evolution of the trajectory of molecule B at t < 0, indicating that this process cannot be described by a two-body collision.37 Note that here time t = 0 was set when the oxygen− oxygen distance of molecules A and B is minimal, right before the evaporation occurs. The MD trajectories of molecules A and B are schematically depicted in Figure 5b. To understand the molecular interactions between molecules A, B, and the other water molecules, we computed the time evolution of the hydrogen bond strengths of the two molecules with each other and their surroundings, through their electronic occupancies.104 This is shown in Figure 5c. The hydrogen bonding of molecule A with the other water molecules, that of molecule B with other water molecules, and the hydrogen bond strength between molecules A and B are represented by the red, blue, and green lines, respectively. This figure indicates that water molecule B is strongly hydrogen bonded to the other water
between SFG experiments and simulation. Recent studies have focused on the bending mode SFG feature at the water−air interfaces.24,94,95 Combining experiment and simulation contributes to a complete understanding of the thermodynamics of the water−air interface as well as establishes a reliable model of water in a heterogeneous environment. II.D. SFG of the O−H Stretch Vibrational Mode at the Water−Surfactant Interface. SFG studies have focused not only on the pure water−air interface but also on electrolyte solution−air interfaces and a variety of aqueous interfaces such as water−surfactant/lipid interfaces. The behavior of ions at the water−air interfaces is important to understand for example the Hofmeister series, while understanding the water−surfactant/ lipid interfaces provides the information on how to control the surface tension of water with those amphiphilic compounds, and the stabilization mechanism of self-assembled lipid bilayers. Experimentally measured phase-sensitive SFG spectra show positive (negative) 3000−3600 cm−1 O−H stretch bands at the water−surfactant interface with negatively (positively) charged headgroups of the surfactants.81,96,97 From this observation, it was concluded that the water dipole moments are pointing uniformly up toward the surfactant layer (down toward the bulk water) when the negatively charged (positively charged) surfactants are present, respectively. In the case of zwitterionic lipid/surfactants, the interpretation of the O−H stretch bands is more complicated. The SFG response of water at the amine oxide (AO) exhibits both negative and positive features, while at the phosphatidylcholine (PC) interface it exhibits only positive features regardless of the presence of the positively charged −N(CH3)3+ group.33,98 The question to answer is thus what the key factor determining the orientation of water near zwitterionic surfactant/lipids is. To address this issue, we computed the SFG spectra of water at these interfaces and disentangled the contributions of water near the positively and negatively charged parts of the SFG spectra.23,27 Figure 4 displays the comparison of the experimentally measured and calculated SFG response at the isotopically diluted
Figure 4. (Left) Comparison of the experimentally measured and simulated O−H stretch SFG spectra at the isotopically diluted water− lipid (PC)/surfactant (AO) interfaces. (Right) The peak assignments are schematically depicted. The figures are reproduced with permission from ref 27. Copyright 2015 American Chemical Society. 3789
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
Figure 5. (a) Time evolutions of the distances between the oxygen atom and M, where M is the midpoint between the oxygen atoms of molecules A and B at t = 0. Mirror images with respect to t = 0 are also plotted as dotted lines, to highlight the asymmetry. (b) Schematics of molecular trajectories of A (blue arrow) and B (green arrow), while the water molecule circled in orange is molecule C. (C) Time evolutions of the electronic occupancies (H bond strength) between molecules A and B (green), between molecule A and other water molecules (red), and between molecule B and other water molecules (blue). This figure is reprinted with permission from ref 37. Copyright 2015 Americal Physical Society.
increases the surface tension.108 Moreover, the addition of tetraphenylarsonium (Ph4As+) and tetraphenylborate (Ph4B−) ions to water, which have similar molecular sizes, structures, and masses, but opposite sign of charge, reduces the surface tension at the water−oil interface in a different manner; the surface tension of the water−oil interface is ∼15% larger for the Ph4B− aqueous solution than for the Ph 4 As + solution, when the ion concentration in the bulk is the same.109 Unveiling the molecular-level physics underlying the surface tensions provides a route to the optimal design rules for controlling surface activities. MD simulation with complex ions/molecules is, however, still challenging, mainly because a proper FF model is sometimes not available and/or has not yet been reported. In particular, the location and orientation of the ions/molecules at the water−air interface are governed by the molecular interactions of these ions/molecules with the interfacial water molecules.110 For an accurate description of the structure and dynamics of these ions/ molecules at the aqueous interface, the quality of the FF model for these ions/molecules should be carefully examined. For example, a FF model for the simple osmolyte TMAO molecule was shown to be insufficiently accurate to reliably describe the hydrogen bond interaction of the TMAO with water molecules; AIMD simulations111 and experiment112 have shown that the FF model113 significantly underestimates the hydrogen bond interaction between the water molecules and oxygen atom of the TMAO molecule, artificially accelerating the rotational motion of water near TMAO molecule. As is clear from this example, combining simulation and experiment can lead to more reliable models for complex molecules/ions in water. III.B. SFG Active Region versus Interfacial Region Contributing to Surface Tension. So far, surface tension has been measured to ensure the surface activity of the solutes in aqueous solutions. Therefore, surface tension measurements are frequently used together with SFG measurement. A recent set of experiments has shown that the surface propensity of ions in solution inferred from SFG does not necessarily correlate with interfacial adsorption concluded from surface tension measurements.109 This raises a fundamental question: Is the SFG active surface area which contributes to the SFG signal equivalent to the area contributing to the surface tension? Although the layer-bylayer contributions of the water molecules to the SFG spectra
molecules (represented by water molecule C) 44 fs before the A−B collision. Through the B−C hydrogen bond interaction, the large potential energy associated with the A−B hydrogen bond at t = 0 is converted into kinetic energy, which is passed to molecule A, resulting in the driving force for molecule A to evaporate (see also Figure 5b). Several hundreds of evaporation occurrences could be traced in time, and all of these follow essentially the same pathway: evaporation follows one specific pathway, involving a delicately timed, concerted motion of several water molecules to “launch” a single molecule from the surface. The notion that one unique pathway exists for water evaporation is incongruous with the previous hypothesis that the evaporation is a purely stochastic process, in which water molecules are thermally activated and those with the high-energy tail of the Boltzmann distribution can evaporate.105,106 Unfortunately, no experiment has yet been performed which can differentiate between these two mechanisms. A possible experiment would be to interrogate the kinetic and rotational energy content of desorbing molecules; the presence of one specific evaporation mechanism indicates the possibility of nonthermal energy distributions. Furthermore, it is also interesting to trace the reverse process of evaporation, i.e., the accommodation process, and investigate whether this process satisfies the requirements of detailed balance.
III. OUTLOOK AND CHALLENGE AT THE AQUEOUS INTERFACE In this section, we discuss four typical challenges related to aqueous interfaces which should be addressed by MD simulations. All of these topics are rather specific, but the required simulation techniques are very general. III.A. Interplay of Complex Ions at Aqueous Interfaces. In contrast to a number of simulation works focusing on surface tension and SFG spectra of pure water and aqueous solution with monatomic ions, the surface tension and SFG spectra of the aqueous solution with complex ions/molecules have not been studied yet in much detail. Complex ions/molecules critically affect the properties of interfacial water, enabling us to optimize the surfactants for specific purposes. For example, the addition of decyldimethylamine N-oxide (DDAO) surfactant dramatically decreases the surface tension of the water,107 whereas adding a trimethylamine N-oxide (TMAO) molecule to water modestly 3790
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
Figure 6. (a) Plots of the depth-dependent contributions to surface tension pN(z) − pT(z) in eq 1 obtained from the real part (top) and reciprocal space part (bottom) potential calculations at 328 K.115 The Gibbs dividing surface is located at |z| ∼20 Å, and vapor regions are located at |z| > 20 Å. (b) Simulated angle profiles formed by the z-axis and O−H (O−D) groups at 300 K.117 The Gibbs dividing surface is located at z = 6.8 Å. (c) Simulated imaginary part of the SFG responses originating from the region of z > zr, for various values of zr at 298 K for pure water (green), NaCl solution (blue), and NaI solution (red).114 The Gibbs dividing surface is located at z = 0 Å, and the region z > 0 is the vapor region. These figures are reprinted with permission from refs 115, 117, and 114. Copyright 1995 American Institute of Physics. Copyright 2012 American Physical Society. Copyright 2007 American Chemical Society.
have been identified at the water−air interface with ions56,114 and at the water−lipid interface,23 the layer-by-layer analyses can be extended to the surface tension,115 allowing us to connect the SFG active region with the surface tension active region. Figure 6 displays the contributions of the layers to the surface tension, the orientation of the water molecule, and the SFG spectra at the water−air interface. Note that the simulations have been performed independently with different simulation techniques and under different conditions. The data indicate that the surface tension originates from the water molecules in the near-surface water layers, ∼6−8 Å from the Gibbs dividing surface (Figure 6a), while the specific orientation of the water molecules persists in the first 4 Å depth from the Gibbs dividing surface (Figure 6b) and equivalently the SFG signal only depends on the ∼4 Å thick water layer (Figure 6c). However, when ions are added into the water, the distributions of the ions can be very inhomogeneous,110 and thus, it is not self-evident that, for electrolyte solutions, the SFG active regions are equally similar to the regions contributing to surface tensions, as first pointed out by Petersen and Saykally.116 This should be surveyed more carefully. III.C. Surface Tension of H2O and D2O. Isotopic dilution is an important tool to avoid the overlaps of different spectral features, to simplify the water vibrational response by reducing inter- and intramolecular couplings, and to test the quantum nature of the water molecules. Many studies have addressed the different behavior of the water molecules with a hydrogen atom and a deuterium atom at the aqueous interfaces: it has been demonstrated, for instance, that the combination of the gasphase H2O and liquid-phase HDO is more stable than the combination of the gas-phase HDO and liquid-phase H2O.118 Moreover, at the water−air interface, HOD molecules straddling the interface preferentially have a free O−H group and a hydrogen bonded OD group, owing to quantum nuclear effects.117,119 This observation indicates that the O−H···O is weaker than the O−D···O hydrogen bond. It seems therefore reasonable to assume that the D2O surface tension is larger than the H2O surface tension. However, this assumption contradicts the experimentally determined surface tensions: the surface
tensions of H2O and D2O at the water−air interface are 71.98 and 71.87 mN/m at 25 °C, respectively. These values apparently arise from a rather complex interplay of the nuclear quantum effects, electrostatic interactions, and van der Waals interactions. Additional light can be shed on these issues using path integral and/or related MD simulation techniques. III.D. Stability and Function of Lipid Bilayer; Atomic vs Coarse-Graining Description. In section II.D, we illustrated the importance of the carbonyl group on the water orientations from FFMD simulation combined with spectroscopy experiments, highlighting the effects of the chemical composition of the lipids on the water. Another example of the theoretical SFG study revealing the effects of the chemical composition is given by ref 99, where it was demonstrated that the size of the positively charged headgroups of the zwitterionic phospholipids critically affects the orientation of the interfacial water. The headgroup of the −NH3+ induces the down-oriented water dipole moment, while near the −N(CH3)3+ headgroup, the orientations of the water molecules are little affected.23,27 These molecular details have been normally excluded in the CG model.120 Moreover, long-ranged interactions are normally neglected in CGMD simulation. It is therefore a question whether or not the hydration of a lipid headgroup can be welldescribed in the CG models. The new findings from the combined SFG experiment and FFMD simulation studies can provide detailed molecular conformations at the aqueous interfaces, which can testify for the reproducibility of the CG models. In this context, it is worth pointing out a technique to bridge the all-atom description in the FFMD simulation and the CG technique, called the adaptive resolution technique.121−124 In this technique, by changing the resolution of the molecular modeling within the geometry-dependent definition, the FFMD and CGMD can be connected smoothly. On the other hand, this adaptive resolution method may generate another artificial “interface” between the regions with the atomic resolution and CG resolution. The effect of these artifacts on the surface tension should also be addressed in more detail. So far, most of the all-atom MD simulations at the water− surfactant/lipid interfaces have been performed using the FF 3791
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B models, while recent ab initio study predicts the strong hydration effects on the electronic structure of the lipid headgroup.125 It is interesting to see how this can be connected with the vibrational signature. This would provide a clear direction for future SFG experiments.
techniques. For the topic discussed in section III.A, probably large-scale AIMD simulation would be required; for the topic discussed in section III.C, the path integral or relevant MD techniques127 are required. For the topic discussed in section III.D, advanced adaptive resolution techniques would be appreciated. As such, addressing these questions also requires the development of new theoretical approaches and algorithms.
IV. CONCLUDING REMARKS In this article, we have attempted to provide an overview of the recent progress on FFMD, AIMD, and CGMD simulations at the aqueous interface reporting on surface tension and the vibrational SFG response, two key quantities characterizing the structure of interfacial water molecules. These three MD simulation techniques have different advantages and intrinsic challenges; AIMD simulations are based on the electrostatic structure theory, and while therefore it is a parameter-free technique, the computational cost is very demanding. Moreover, approximations, such as the generalized gradient approximation, cause serious limitation for reproducing experimental observables. FFMD simulation is much faster than AIMD simulation, yet different force field models may cause different physical pictures. CGMD simulation is even faster than FFMD simulation, but the details of the chemical structure are ignored, making the reproducibility of the experimental data sometimes questionable.126 For the computationally expensive AIMD simulation, we have introduced an efficient algorithm to calculate the SFG responses. AIMD SFG simulations for water−air and water−lipid/ surfactant indicate that the SFG spectra with AIMD simulation agree with the experimental data. On the other hand, the surface tension values show a very large variability from 30 to 170 mN/m at different levels of theory (the experimental value is 72 mN/m), indicating that an improved description of the long-range interactions is required. For FFMD simulation, the surface tension is well-reproduced with popular water models, while the different force field models still have difficulties in reproducing the SFG spectra. For CGMD simulation, more detailed studies to examine the reproducibility of the surface tension would be required. While combining the adaptive resolution technique with the CGMD simulation is promising, calculating the surface tension of the aqueous interface would be also a critical check for the technique. We also illustrated that, with FFMD simulation, we can predict a molecular mechanism of the water evaporation: The welltimed, concerted motion of the interfacial water molecules can bring the kinetic energy to the water molecule, which is the driving force for this molecule to break its hydrogen bond and fly away from the water surface. Here, a challenge is to find the experimental evidence that the water evaporation is controlled by the concerted motions of several water molecules. Despite the substantial progress in simulating aqueous interfaces and the answers that have been obtained for important questions, many challenges remain. Four challenges, also relevant for the experimentalists, which can be answered in the coming years, are in our opinion the following: (1) Can MD simulation estimate and understand the influence of complex ions on the surface tension? (2) How does the region contributing to the surface tension differ from the SFG active region, especially for different systems? (3) Why is the H2O−air surface tension slightly higher than the D2O−air surface tension? (4) To which extent can a CG simulation model the details of the chemical structure of the lipid/surfactant? Although the challenges are not limited to these four topics, these topics include the important ingredients for the state-of-the-art
■
AUTHOR INFORMATION
Corresponding Authors
*E-mail:
[email protected]. *E-mail:
[email protected]. Phone: +49-6131-379161. Notes
The authors declare no competing financial interest. Biographies
Yuki Nagata was born in Kobe, Japan, in 1979. He received his bachelor and master degrees from the University of Tokyo, and his Ph.D. degree from Kyoto University. He worked as a research scientist at BASF and as a postdoctoral researcher at the group of Shaul Mukamel, University of California, Irvine. After that, he joined the Max Planck Institute for Polymer Research, Mainz, Germany, as a group leader of the Surface Spectroscopy Theory group. His research focus is the molecular mechanism of the aqueous interface and the theoretical design of vibrational spectroscopy.
Tatsuhiko Ohto is a research associate at Osaka University, Osaka, Japan. He received his Ph.D. degree from the University of Tokyo in 2013. During his Ph.D. course, he spent 5 months as a visiting student in the Surface Spectroscopy Theory group at the Max Planck Institute for Polymer Research, Mainz, Germany. After postdoctoral research at the Advanced Institute of Science and Technology, Tsukuba, Japan, he joined Osaka University. His research interest is the theoretical modeling, primarily based on first-principles calculations, of the 3792
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B
(4) Møller, K. B.; Rey, R.; Hynes, J. T. Hydrogen Bond Dynamics in Water and Ultrafast Infrared Spectroscopy: A Theoretical Study. J. Phys. Chem. A 2004, 108, 1275−1289. (5) Nagata, Y.; Yoshimune, S.; Hsieh, C.-S.; Hunger, J.; Bonn, M. Ultrafast Vibrational Dynamics of Water Disentangled by Reverse Nonequilibrium Ab Initio Molecular Dynamics Simulations. Phys. Rev. X 2015, 5, 021002. (6) Matsumoto, M.; Saito, S.; Ohmine, I. Molecular Dynamics Simulation of the Ice Nucleation and Growth Process Leading to Water Freezing. Nature 2002, 416, 409−413. (7) Marx, D.; Tuckerman, M. E.; Hutter, J.; Parrinello, M. The Nature of the Hydrated Excess Proton in Water. Nature 1999, 397, 601−604. (8) Habershon, S.; Manolopoulos, D. E.; Markland, T. E.; Miller, T. F. Ring-Polymer Molecular Dynamics: Quantum Effects in Chemical Dynamics from Classical Trajectories in an Extended Phase Space. Annu. Rev. Phys. Chem. 2013, 64, 387−413. (9) Car, R.; Parrinello, M. Unified Approach for Molecular Dynamics and Density-Functional Theory. Phys. Rev. Lett. 1985, 55, 2471−2474. (10) Kühne, T. D.; Krack, M.; Mohamed, F. R.; Parrinello, M. Efficient and Accurate Car-Parrinello-like Approach to Born-Oppenheimer Molecular Dynamics. Phys. Rev. Lett. 2007, 98, 066401. (11) Kühne, T. D. Second Generation Car-Parrinello Molecular Dynamics. WIREs Comput. Mol. Sci. 2014, 4, 391−406. (12) Del Ben, M.; Schönherr, M.; Hutter, J.; Vandevondele, J. Bulk Liquid Water at Ambient Temperature and Pressure from MP2 Theory. J. Phys. Chem. Lett. 2013, 4, 3753−3759. (13) Lee, C.-Y.; McCammon, J. A.; Rossky, P. J. The Structure of Liquid Water at an Extended Hydrophobic Surface. J. Chem. Phys. 1984, 80, 4448−4455. (14) Vega, C.; De Miguel, E. Surface Tension of the Most Popular Models of Water by Using the Test-Area Simulation Method. J. Chem. Phys. 2007, 126, 154707. (15) Nagata, Y.; Ohto, T.; Bonn, M.; Kühne, T. D. Surface Tension of Ab Initio Liquid Water at the Water-Air Interface. arXiv:1603.01658. (16) Du, Q.; Freysz, E.; Shen, Y. R. Surface Vibrational Spectroscopic Studies of Hydrogen Bonding and Hydrophobicity. Science 1994, 264, 826−828. (17) Ishiyama, T.; Imamura, T.; Morita, A. Theoretical Studies of Structures and Vibrational Sum Frequency Generation Spectra at Aqueous Interfaces. Chem. Rev. 2014, 114, 8447−8470. (18) Geissler, P. L. Water Interfaces, Solvation, and Spectroscopy. Annu. Rev. Phys. Chem. 2013, 64, 317−337. (19) Morita, A.; Hynes, J. T. A Theoretical Analysis of the Sum Frequency Generation Spectrum of the Water Surface. Chem. Phys. 2000, 258, 371−390. (20) Morita, A.; Hynes, J. T. A Theoretical Analysis of the Sum Frequency Generation Spectrum of the Water Surface. II. TimeDependent Approach. J. Phys. Chem. B 2002, 106, 673−685. (21) Stiopkin, I. V.; Weeraman, C.; Pieniazek, P. A.; Shalhout, F. Y.; Skinner, J. L.; Benderskii, A. V. Hydrogen Bonding at the Water Surface Revealed by Isotopic Dilution Spectroscopy. Nature 2011, 474, 192− 195. (22) Pieniazek, P. A.; Tainter, C. J.; Skinner, J. L. Surface of Liquid Water: Three-Body Interactions and Vibrational Sum-Frequency Spectroscopy. J. Am. Chem. Soc. 2011, 133, 10360−10363. (23) Nagata, Y.; Mukamel, S. Vibrational Sum-Frequency Generation Spectroscopy at the Water/Lipid Interface: Molecular Dynamics Simulation Study. J. Am. Chem. Soc. 2010, 132, 6434−6442. (24) Nagata, Y.; Hsieh, C.-S.; Hasegawa, T.; Voll, J.; Backus, E. H. G.; Bonn, M. Water Bending Mode at the Water−Vapor Interface Probed by Sum-Frequency Generation Spectroscopy: A Combined Molecular Dynamics Simulation and Experimental Study. J. Phys. Chem. Lett. 2013, 4, 1872−1877. (25) Nagata, Y.; Hasegawa, T.; Backus, E. H. G.; Usui, K.; Yoshimune, S.; Ohto, T.; Bonn, M. The Surface Roughness, but Not the Water Molecular Orientation Varies with Temperature at the Water−air Interface. Phys. Chem. Chem. Phys. 2015, 17, 23559−23564. (26) Ohto, T.; Usui, K.; Hasegawa, T.; Bonn, M.; Nagata, Y. Toward Ab Initio Molecular Dynamics Modeling for Sum-Frequency Gen-
structure and dynamics of molecules at interfaces and electron transport of metal−molecule−metal systems.
Ellen H. G. Backus is the group leader of the Water at Interfaces group in the Molecular Spectroscopy department at the Max Planck Institute for Polymer Research (MPIP), Mainz, supported by a Minerva Grant from the Max Planck Society, and an ERC Starting Grant. She obtained her Ph.D. in 2005 at the Leiden University in the group of Mischa Bonn and Aart Kleyn. After postdoctoral research at the University of Zurich with Peter Hamm, she moved back to The Netherlands for an independent postdoctoral researcher position in the group of Huib Bakker at AMOLF in Amsterdam. In 2012, she joined the MPIP. Her research focuses on the structure and dynamics of water at various interfaces.
Mischa Bonn is a director at the Max Planck Institute for Polymer Research (MPIP), Mainz, Germany. He works on label-free (ultrafast) vibrational spectroscopy and microscopy of biomolecular systems and water in such systems. He received his Ph.D. in 1996 from the University of Eindhoven for research performed at AMOLF in Amsterdam. After postdoctoral research at the Fritz Haber Institute in Berlin and Columbia University in New York, he worked at Leiden University from 1999 as an assistant and associate professor. In 2004 he became a group leader at AMOLF. In 2011 he joined the MPIP. His research interests are the structure and dynamics of molecules at interfaces, and electron transfer across interfaces.
■
REFERENCES
(1) Rahman, A.; Stillinger, F. H. Molecular Dynamics Study of Liquid Water. J. Chem. Phys. 1971, 55, 3336−3359. (2) Luzar, A.; Chandler, D. Hydrogen-Bond Kinetics in Liquid Water. Nature 1996, 379, 55−57. (3) Bakker, H. J.; Skinner, J. L. Vibrational Spectroscopy as a Probe of Structure and Dynamics in Liquid Water. Chem. Rev. 2010, 110, 1498− 1517. 3793
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B eration Spectra; an Efficient Algorithm Based on Surface-Specific Velocity-Velocity Correlation Function. J. Chem. Phys. 2015, 143, 124702. (27) Ohto, T.; Backus, E. H. G.; Hsieh, C.-S.; Sulpizi, M.; Bonn, M.; Nagata, Y. Lipid Carbonyl Groups Terminate the Hydrogen Bond Network of Membrane-Bound Water. J. Phys. Chem. Lett. 2015, 6, 4499−4503. (28) Sulpizi, M.; Salanne, M.; Sprik, M.; Gaigeot, M.-P. Vibrational Sum Frequency Generation Spectroscopy of the Water Liquid−vapor Interface from Density Functional Theory-Based Molecular Dynamics Simulations. J. Phys. Chem. Lett. 2013, 4, 83. (29) Wan, Q.; Galli, G. First-Principles Framework to Compute SumFrequency Generation Vibrational Spectra of Semiconductors and Insulators. Phys. Rev. Lett. 2015, 115, 1−5. (30) Tian, C.; Shen, Y. R. Isotopic Dilution Study of the Water/vapor Interface by Phase-Sensitive Sum-Frequency Vibrational Spectroscopy. J. Am. Chem. Soc. 2009, 131, 2790−2791. (31) Yamaguchi, S. Development of Single-Channel HeterodyneDetected Sum Frequency Generation Spectroscopy and Its Application to the Water/vapor Interface. J. Chem. Phys. 2015, 143, 034202. (32) Ishiyama, T.; Morita, A. Vibrational Spectroscopic Response of Intermolecular Orientational Correlation at the Water Surface. J. Phys. Chem. C 2009, 113, 16299−16302. (33) Mondal, J. A.; Nihonyanagi, S.; Yamaguchi, S.; Tahara, T. Three Distinct Water Structures at a Zwitterionic Lipid/water Interface Revealed by Heterodyne-Detected Vibrational Sum Frequency Generation. J. Am. Chem. Soc. 2012, 134, 7842−7850. (34) Re, S.; Nishima, W.; Tahara, T.; Sugita, Y. Mosaic of Water Orientation Structures at a Neutral Zwitterionic Lipid/Water Interface Revealed by Molecular Dynamics Simulations. J. Phys. Chem. Lett. 2014, 5, 4343−4348. (35) Mason, P. E. Molecular Dynamics Study on the Microscopic Details of the Evaporation of Water. J. Phys. Chem. A 2011, 115, 6054− 6058. (36) Varilly, P.; Chandler, D. Water Evaporation: A Transition Path Sampling Study. J. Phys. Chem. B 2013, 117, 1419−1428. (37) Nagata, Y.; Usui, K.; Bonn, M. Molecular Mechanism of Water Evaporation. Phys. Rev. Lett. 2015, 115, 236102. (38) Giovambattista, N.; Debenedetti, P. G.; Rossky, P. J. Effect of Surface Polarity on Water Contact Angle and Interfacial Hydration Structure. J. Phys. Chem. B 2007, 111, 9581−9587. (39) Taherian, F.; Marcon, V.; Van Der Vegt, N. F. A.; Leroy, F. What Is the Contact Angle of Water on Graphene? Langmuir 2013, 29, 1457− 1465. (40) Rafiee, J.; Mi, X.; Gullapalli, H.; Thomas, A. V.; Yavari, F.; Shi, Y.; Ajayan, P. M.; Koratkar, N. A. Wetting Transparency of Graphene. Nat. Mater. 2012, 11, 217−222. (41) Marrink, S. J.; Mark, A. E. The Mechanism of Vesicle Fusion as Revealed by Molecular Dynamics Simulations. J. Am. Chem. Soc. 2003, 125, 11144−11145. (42) Stevens, M. J.; Hoh, J. H.; Woolf, T. B. Insights into the Molecular Mechanism of Membrane Fusion from Simulation: Evidence for the Association of Splayed Tails. Phys. Rev. Lett. 2003, 91, 188102. (43) Gurtovenko, A. A.; Vattulainen, I. Pore Formation Coupled to Ion Transport through Lipid Membranes as Induced by Transmembrane Ionic Charge Imbalance: Atomistic Molecular Dynamics Study. J. Am. Chem. Soc. 2005, 127, 17570−17571. (44) Ishiyama, T.; Morita, A. Analysis of Anisotropic Local Field in Sum Frequency Generation Spectroscopy with the Charge Response Kernel Water Model. J. Chem. Phys. 2009, 131, 244714. (45) Pieniazek, P. A.; Tainter, C. J.; Skinner, J. L. Interpretation of the Water Surface Vibrational Sum-Frequency Spectrum. J. Chem. Phys. 2011, 135, 044701. (46) Marx, D.; Hutter, J. Ab Initio Molecular Dynamics: Basic Theory and Advanced Methods; Cambridge University Press, 2009. (47) Hutter, J.; Iannuzzi, M.; Schiffmann, F.; VandeVondele, J. CP2K: Atomistic Simulations of Condensed Matter Systems. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2014, 4, 15−25.
(48) Klimeš, J.; Michaelides, A. Perspective: Advances and Challenges in Treating van Der Waals Dispersion Forces in Density Functional Theory. J. Chem. Phys. 2012, 137, 120901. (49) Kuo, I. F.; Mundy, C. J. An Ab Initio Molecular Dynamics Study of the Aqueous Liquid-Vapor Interface. Science 2004, 303, 658−660. (50) Baer, M. D.; Mundy, C. J.; McGrath, M. J.; Kuo, I.-F. W.; Siepmann, J. I.; Tobias, D. J. Re-Examining the Properties of the Aqueous Vapor−liquid Interface Using Dispersion Corrected Density Functional Theory. J. Chem. Phys. 2011, 135, 124712. (51) Baer, M. D.; Kuo, I. F. W.; Tobias, D. J.; Mundy, C. J. Toward a Unified Picture of the Water Self-Ions at the Air-Water Interface: A Density Functional Theory Perspective. J. Phys. Chem. B 2014, 118, 8364−8372. (52) Kühne, T. D.; Pascal, T. A.; Kaxiras, E.; Jung, Y. New Insights into the Structure of the Vapor/Water Interface from Large-Scale FirstPrinciples Simulations. J. Phys. Chem. Lett. 2011, 2, 105−113. (53) Tobias, D. J.; Stern, A. C.; Baer, M. D.; Levin, Y.; Mundy, C. J. Simulation and Theory of Ions at Atmospherically Relevant Aqueous Liquid-Air Interfaces. Annu. Rev. Phys. Chem. 2013, 64, 339−359. (54) Toukmaji, A. Y.; Board, J. A. Ewald Summation Techniques in Perspective: A Survey. Comput. Phys. Commun. 1996, 95, 73−92. (55) Tuckerman, M. E.; Berne, B. J.; Martyna, G. J. Reversible Multiple Time Scale Molecular Dynamics. J. Chem. Phys. 1992, 97, 1990−2001. (56) Ishiyama, T.; Morita, A. Molecular Dynamics Study of Gas-Liquid Aqueous Sodium Halide Interfaces. I. Flexible and Polarizable Molecular Modeling and Interfacial Properties. J. Phys. Chem. C 2007, 111, 721− 737. (57) Tainter, C. J.; Pieniazek, P. A.; Lin, Y.-S.; Skinner, J. L. Robust Three-Body Water Simulation Model. J. Chem. Phys. 2011, 134, 184501. (58) Hasegawa, T.; Tanimura, Y. A Polarizable Water Model for Intramolecular and Intermolecular Vibrational Spectroscopies. J. Phys. Chem. B 2011, 115, 5545−5553. (59) Medders, G. R.; Paesani, F. Infrared and Raman Spectroscopy of Liquid Water through “First-Principles” Many-Body Molecular Dynamics. J. Chem. Theory Comput. 2015, 11, 1145−1154. (60) Bowman, J. M.; Wang, Y.; Liu, H.; Mancini, J. S. Ab Initio Quantum Approaches to the IR Spectroscopy of Water and Hydrates. J. Phys. Chem. Lett. 2015, 6, 366−373. (61) Laury, M. L.; Wang, L.-P.; Pande, V. S.; Head-Gordon, T.; Ponder, J. W. Revised Parameters for the AMOEBA Polarizable Atomic Multipole Water Model. J. Phys. Chem. B 2015, 119, 9423−9437. (62) Fanourgakis, G. S.; Xantheas, S. S. Development of Transferable Interaction Potentials for Water. V. Extension of the Flexible, Polarizable, Thole-Type Model Potential (TTM3-F, v. 3.0) to Describe the Vibrational Spectra of Water Clusters and Liquid Water. J. Chem. Phys. 2008, 128, 074506. (63) Sedlmeier, F.; Horinek, D.; Netz, R. R. Nanoroughness, Intrinsic Density Profile, and Rigidity of the Air-Water Interface. Phys. Rev. Lett. 2009, 103, 136102. (64) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; de Vries, A. H. The MARTINI Force Field: Coarse Grained Model for Biomolecular Simulations. J. Phys. Chem. B 2007, 111, 7812−7824. (65) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N· log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (66) Gloor, G. J.; Jackson, G.; Blas, F. J.; de Miguel, E. Test-Area Simulation Method for the Direct Determination of the Interfacial Tension of Systems with Continuous or Discontinuous Potentials. J. Chem. Phys. 2005, 123, 134703. (67) Schmitz, F.; Virnau, P.; Binder, K. Logarithmic Finite-Size Effects on Interfacial Free Energies: Phenomenological Theory and Monte Carlo Studies. Phys. Rev. E 2014, 90, 012128. (68) Schmitz, F.; Virnau, P.; Binder, K. Determination of the Origin and Magnitude of Logarithmic Finite-Size Effects on Interfacial Tension: Role of Interfacial Fluctuations and Domain Breathing. Phys. Rev. Lett. 2014, 112, 125701. (69) Yesylevskyy, S. O.; Schäfer, L. V.; Sengupta, D.; Marrink, S. J. Polarizable Water Model for the Coarse-Grained MARTINI Force Field. PLoS Comput. Biol. 2010, 6, e1000810. 3794
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B (70) Marx, D.; Hutter, J. Ab Initio Molecular Dynamics: Basic Theory and Advanced Methods; Cambridge University Press, 2009. (71) Grimme, S. Semiempirical GGA-Type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787−1799. (72) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (73) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (74) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (75) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (76) Zhang, Y.; Yang, W. Comment on “Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1998, 80, 890−890. (77) Willow, S. Y.; Zeng, X. C.; Xantheas, S. S.; Kim, K. S.; Hirata, S. Why Is MP2-Water “Cooler” and “Denser” than DFT-Water? J. Phys. Chem. Lett. 2016, 7, 680−684. (78) Willow, S. Y.; Salim, M. A.; Kim, K. S.; Hirata, S. Ab Initio Molecular Dynamics of Liquid Water Using Embedded- Fragment Second-Order Many-Body Perturbation Theory towards Its Accurate Property Prediction. Sci. Rep. 2015, 5, 14358. (79) Ostroverkhov, V.; Waychunas, G.; Shen, Y. R. New Information on Water Interfacial Structure Revealed by Phase-Sensitive Surface Spectroscopy. Phys. Rev. Lett. 2005, 94, 046102. (80) Stiopkin, I. V.; Jayathilake, H. D.; Bordenyuk, A. N.; Benderskii, A. V. Heterodyne-Detected Vibrational Sum Frequency Generation Spectroscopy. J. Am. Chem. Soc. 2008, 130, 2271. (81) Nihonyanagi, S.; Yamaguchi, S.; Tahara, T. Direct Evidence for Orientational Flip-Flop of Water Molecules at Charged Interfaces: A Heterodyne-Detected Vibrational Sum Frequency Generation Study. J. Chem. Phys. 2009, 130, 204704. (82) Perry, A.; Neipert, C.; Kasprzyk, C. R.; Green, T.; Space, B.; Moore, P. B. A Theoretical Description of the Polarization Dependence of the Sum Frequency Generation Spectroscopy of the Water/vapor Interface. J. Chem. Phys. 2005, 123, 144705. (83) Auer, B. M.; Skinner, J. L. IR and Raman Spectra of Liquid Water: Theory and Interpretation. J. Chem. Phys. 2008, 128, 224511. (84) Auer, B. M.; Skinner, J. L. Vibrational Sum-Frequency Spectroscopy of the Water Liquid/Vapor Interface. J. Phys. Chem. B 2009, 113, 4125−4130. (85) Auer, B. M.; Skinner, J. L. Vibrational Sum-Frequency Spectroscopy of the Liquid/vapor Interface for Dilute HOD in D2O. J. Chem. Phys. 2008, 129, 214705. (86) Nihonyanagi, S.; Ishiyama, T.; Lee, T.-K.; Yamaguchi, S.; Bonn, M.; Morita, A.; Tahara, T. Unified Molecular View of the Air/water Interface Based on Experimental and Theoretical χ(2) Spectra of an Isotopically Diluted Water Surface. J. Am. Chem. Soc. 2011, 133, 16875− 16880. (87) Verreault, D.; Hua, W.; Allen, H. C. From Conventional to PhaseSensitive Vibrational Sum Frequency Generation Spectroscopy: Probing Water Organization at Aqueous Interfaces. J. Phys. Chem. Lett. 2012, 3, 3012−3028. (88) Ishiyama, T.; Takahashi, H.; Morita, A. Molecular Dynamics Simulations of Surface-Specific Bonding of the Hydrogen Network of Water: A Solution to the Low Sum-Frequency Spectra. Phys. Rev. B 2012, 86, 035408. (89) Nihonyanagi, S.; Kusaka, R.; Inoue, K.; Adhikari, A.; Yamaguchi, S. Accurate Determination of Complex χ(2) Spectrum of the Air/water Interface. J. Chem. Phys. 2015, 143, 124707. (90) Egorov, S. A.; Skinner, J. L. Semiclassical Approximations to Quantum Correlation Functions. Chem. Phys. Lett. 1998, 293, 469. (91) Corcelli, S. A.; Skinner, J. L. Infrared and Raman Line Shapes of Dilute HOD in Liquid H2O and D2O from 10 to 90 Degrees C. J. Phys. Chem. A 2005, 109, 6154−6165.
(92) Paesani, F.; Iuchi, S.; Voth, G. A. Quantum Effects in Liquid Water from an Ab Initio-Based Polarizable Force Field. J. Chem. Phys. 2007, 127, 074506. (93) Medders, G. R.; Paesani, F. Dissecting the Molecular Structure of the Air/Water Interface from Quantum Simulations of the SumFrequency Generation Spectrum. J. Am. Chem. Soc. 2016, 138, 3912. (94) Vinaykin, M.; Benderskii, A. V. Vibrational Sum-Frequency Spectrum of the Water Bend at the Air/ Water Interface. J. Phys. Chem. Lett. 2012, 3, 3348. (95) Ni, Y.; Skinner, J. L. IR and SFG Vibrational Spectroscopy of the Water Bend in the Bulk Liquid and at the Liquid-Vapor Interface, Respectively. J. Chem. Phys. 2015, 143, 014502. (96) Livingstone, R. A.; Nagata, Y.; Bonn, M.; Backus, E. H. G. Two Types of Water at the Water-Surfactant Interface Revealed by TimeResolved Vibrational Spectroscopy. J. Am. Chem. Soc. 2015, 137, 14912−14919. (97) Wen, Y.; Zha, S.; Liu, X.; Yang, S.; Guo, P.; Shi, G.; Fang, H.; Shen, Y. R.; Tian, C. Unveiling Microscopic Structures of Charged Water Interfaces by Surface-Specific Vibrational Spectroscopy. Phys. Rev. Lett. 2016, 116, 016101. (98) Chen, X.; Hua, W.; Huang, Z.; Allen, H. C. Interfacial Water Structure Associated with Phospholipid Membranes Studied by PhaseSensitive Vibrational Sum Frequency Generation Spectroscopy. J. Am. Chem. Soc. 2010, 132, 11336−11342. (99) Roy, S.; Gruenbaum, S. M.; Skinner, J. L. Theoretical Vibrational Sum-Frequency Generation Spectroscopy of Water near Lipid and Surfactant Monolayer Interfaces. J. Chem. Phys. 2014, 141, 18C502. (100) Ishiyama, T.; Terada, D.; Morita, A. Hydrogen-Bonding Structure at Zwitterionic Lipid/Water Interface. J. Phys. Chem. Lett. 2016, 7, 216−220. (101) Gawrisch, K.; Ruston, D.; Zimmerberg, J.; Parsegian, V. A.; Rand, R. P.; Fuller, N. Membrane Dipole Potentials, Hydration Forces, and the Ordering of Water at Membrane Surfaces. Biophys. J. 1992, 61, 1213−1223. (102) Razafindralambo, H.; Blecker, C.; Mezdour, S.; Deroanne, C.; Crowet, J.-M.; Brasseur, R.; Lins, L.; Paquot, M. Impacts of the Carbonyl Group Location of Ester Bond on Interfacial Properties of Sugar-Based Surfactants: Experimental and Computational Evidences. J. Phys. Chem. B 2009, 113, 8872−8877. (103) Musolino, N.; Trout, B. L. Insight into the Molecular Mechanism of Water Evaporation via the Finite Temperature String Method. J. Chem. Phys. 2013, 138, 134707. (104) Kumar, R.; Schmidt, J.; Skinner, J. Hydrogen Bonding Definitions and Dynamics in Liquid Water. J. Chem. Phys. 2007, 126, 204107. (105) Ward, C. A.; Stanga, D. Interfacial Conditions during Evaporation or Condensation of Water. Phys. Rev. E 2001, 64, 051509. (106) Smith, J. D.; Cappa, C. D.; Drisdell, W. S.; Cohen, R. C.; Saykally, R. J. Raman Thermometry Measurements of Free Evaporation from Liquid Water Droplets. J. Am. Chem. Soc. 2006, 128, 12892− 12898. (107) Tuncay, M.; Cristian, S. D. A Study of the Binding of Dimethldodecylamine Oxide by Beta-Cyclodextrin Using Surface Tension Measurements. J. Colloid Interface Sci. 1994, 167, 181−185. (108) Sagle, L. B.; Cimatu, K.; Litosh, V. A.; Liu, Y.; Flores, S. C.; Chen, X.; Yu, B.; Cremer, P. S. Methyl Groups of Trimethylamine N-Oxide Orient Away from Hydrophobic Interfaces. J. Am. Chem. Soc. 2011, 133, 18707−18712. (109) Carrier, O.; Backus, E. H. G.; Shahidzadeh, N.; Franz, J.; Wagner, M.; Nagata, Y.; Bonn, M.; Bonn, D. Oppositely Charged Ions at WaterAir and Water-Oil Interfaces: Unifying the Molecular Picture with Thermodynamics. J. Phys. Chem. Lett. 2016, 7, 825−830. (110) Jungwirth, P.; Tobias, D. J. Molecular Structure of Salt Solutions: A New View of the Interface with Implications for Heterogeneous Atmospheric Chemistry. J. Phys. Chem. B 2001, 105, 10468−10472. (111) Usui, K.; Hunger, J.; Sulpizi, M.; Ohto, T.; Bonn, M.; Nagata, Y. Ab Initio Liquid Water Dynamics in Aqueous TMAO Solution. J. Phys. Chem. B 2015, 119, 10597−10606. 3795
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796
Feature Article
The Journal of Physical Chemistry B (112) Hunger, J.; Tielrooij, K. J.; Buchner, R.; Bonn, M.; Bakker, H. J. Complex Formation in Aqueous Trimethylamine-N-Oxide (TMAO) Solutions. J. Phys. Chem. B 2012, 116, 4783−4795. (113) Kast, K. M.; Brickmann, J.; Kast, S. M.; Berry, R. S. Binary Phases of Aliphatic N-Oxides and Water: Force Field Development and Molecular Dynamics Simulation. J. Phys. Chem. A 2003, 107, 5342− 5351. (114) Ishiyama, T.; Morita, A. Molecular Dynamics Study of GasLiquid Aqueous Sodium Halide Interfaces. II. Analysis of Vibrational Sum Frequency Generation Spectra. J. Phys. Chem. C 2007, 111, 738− 748. (115) Alejandre, J.; Tildesley, D. J.; Chapela, G. A. Molecular Dynamics Simulation of the Orthobaric Densities and Surface Tension of Water. J. Chem. Phys. 1995, 102, 4574−4583. (116) Petersen, P. B.; Saykally, R. J. On the Nature of Ions at the Liquid Water Surface. Annu. Rev. Phys. Chem. 2006, 57, 333−364. (117) Nagata, Y.; Pool, R. E.; Backus, E. H. G.; Bonn, M. Nuclear Quantum Effects Affect Bond Orientation of Water at the Water-Vapor Interface. Phys. Rev. Lett. 2012, 109, 226101. (118) Markland, T. E.; Berne, B. J. Unraveling Quantum Mechanical Effects in Water Using Isotopic Fractionation. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 7988−7991. (119) Liu, J.; Andino, R. S.; Miller, C. M.; Chen, X.; Wilkins, D. M.; Ceriotti, M.; Manolopoulos, D. E. A Surface-Specific Isotope Effect in Mixtures of Light and Heavy Water. J. Phys. Chem. C 2013, 117, 2944− 2951. (120) Ingólfsson, H. I.; Lopez, C. A.; Uusitalo, J. J.; de Jong, D. H.; Gopal, S. M.; Periole, X.; Marrink, S. J. The Power of Coarse Graining in Biomolecular Simulations. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2014, 4, 225−248. (121) Praprotnik, M.; Site, L. D.; Kremer, K. Multiscale Simulation of Soft Matter: From Scale Bridging to Adaptive Resolution. Annu. Rev. Phys. Chem. 2008, 59, 545−571. (122) Potestio, R.; Fritsch, S.; Español, P.; Delgado-Buscalioni, R.; Kremer, K.; Everaers, R.; Donadio, D. Hamiltonian Adaptive Resolution Simulation for Molecular Liquids. Phys. Rev. Lett. 2013, 110, 108301. (123) Mukherji, D.; Kremer, K. Coil-Globule-Coil Transition of PNIPAm in Aqueous Methanol: Coupling All-Atom Simulations to Semi-Grand Canonical Coarse-Grained Reservoir. Macromolecules 2013, 46, 9158−9163. (124) Fogarty, A. C.; Potestio, R.; Kremer, K. Adaptive Resolution Simulation of a Biomolecule and Its Hydration Shell: Structural and Dynamical Properties. J. Chem. Phys. 2015, 142, 195101. (125) Mashaghi, A.; Partovi-Azar, P.; Jadidi, T.; Nafari, N.; Maass, P.; Tabar, M. R. R.; Bonn, M.; Bakker, H. J. Hydration Strongly Affects the Molecular and Electronic Structure of Membrane Phospholipids. J. Chem. Phys. 2012, 136, 114709. (126) Vaikuntanathan, S.; Rotskoff, G. M.; Geissler, P.; Hudson, A. Necessity of Capillary Modes in a Minimal Model of Nanoscale Hydrophobic Solvation. Proc. Natl. Acad. Sci. U. S. A. 2016. (127) Kessler, J.; Elgabarty, H.; Spura, T.; Karhan, K.; Partovi-Azar, P.; Hassanali, A. A.; Kuhne, T. D. Structure and Dynamics of the Instantaneous Water/Vapor Interface Revisited by Path-Integral and Ab-Initio Molecular Dynamics Simulations. J. Phys. Chem. B 2015, 119, 10079−10086.
3796
DOI: 10.1021/acs.jpcb.6b01012 J. Phys. Chem. B 2016, 120, 3785−3796