Role of Presolvation and Anharmonicity in Aqueous Phase Hydrated

Article pubs.acs.org/JPCB

Role of Presolvation and Anharmonicity in Aqueous Phase Hydrated Proton Solvation and Transport Rajib Biswas, Ying-Lung Steve Tse,§ Andrei Tokmakoff, and Gregory A. Voth* Department of Chemistry, James Franck Institute, and Institute for Biophysical Dynamics, The University of Chicago, Chicago, Illinois 60637, United States S Supporting Information *

ABSTRACT: Results from condensed phase ab initio molecular dynamics (AIMD) simulations suggest a proton transfer reaction is facilitated by “presolvation” in which the hydronium is transiently solvated by four water molecules, similar to the typical solvation structure of water, by accepting a weak hydrogen bond from the fourth water molecule. A new version 3.2 multistate empirical valence bond (MS-EVB 3.2) model for the hydrated excess proton incorporating this presolvation behavior is therefore developed. The classical MSEVB simulations show similar structural properties as those of the previous model but with significantly improved diffusive behavior. The inclusion of nuclear quantum effects in the MSEVB also provides an even better description of the proton diffusion rate. To quantify the influence of anharmonicity, a second model (aMS-EVB 3.2) is developed using the anharmonic aSPC/Fw water model, which provides similar structural properties but improved spectroscopic responses at high frequencies.

I. INTRODUCTION Proton solvation and transport is ubiquitous1−9 and is involved in some of the most fundamental and important processes in biological10−13 and energy systems.14−17 A significant effort has been undertaken to understand the proton solvation and transport in different systems by experiment 18−24 and simulation.3−9,13,25−41 Thanks to the development of highresolution spectroscopic measurements18−24 and computer simulation methods,3−9,13,15−17,25−41 the understanding of the proton transport mechanisms, especially for the case in water, has advanced greatly. One reason why the proton transport is different from most other charge transport processes in aqueous systems is due to the Grotthuss hopping mechanism, in which an excess proton is transferred from one water to another by breaking and forming chemical bonds.1−9,11−21,25−33,35−41 Together with vehicular transport (standard diffusion), the hydrated excess proton is the fastest moving ion in water. Because of the Grotthuss mechanism, a reasonable model to study proton transfer between water molecules must be able to correctly describe the dynamical water hydrogen bond (H-bond) network that constantly changes along with the bonding topology. The most natural computational method to study proton solvation and transport is ab initio molecular dynamics (AIMD) simulations in which the Hamiltonian of the system explicitly accounts for the electronic (and nuclear) degrees of freedom.3,6−8,25,27,31,33,35,38−40,42 However, this approach also comes with very high computational cost and, therefore, limits the time and length scales that can be studied. The inherent © 2015 American Chemical Society

approximations associated with the density functional theory (DFT) utilized in AIMD can lead to unphysical results, especially in the case of liquid water.27,40,43,44 However, despite these concerns and because AIMD simulations can naturally account for electronic polarization and charge transfer, they can provide important information for parametrization of more efficient models.34,36 The multistate empirical valence bond (MS-EVB) model has been shown to successfully capture the essential physics and chemistry in different protonated systems.5,13,26,29,37 In this paper, we aim to improve the original version 3 MS-EVB model29 (MS-EVB 3.0) by incorporating the effect of water presolvation,31 which is characterized by the presence of a transient fourth water neighbor of hydronium that is found to influence proton dynamics.40 Since a water molecule is typically 4-fold coordinated, this presolvation allows the hydronium to also be transiently 4-fold coordinated, which creates a similar solvation environment for both the hydronium and its nearby water molecules.31,35,40 Thus, if a proton transfer reaction for this hydronium happens, this newly formed water molecule remains in a typical 4-fold coordinated configuration and can be stabilized. This presolvation idea in which the newly formed water is stabilized can be analyzed by the Marcus picture of proton transfer,45,46 in which the free energy of the system Special Issue: Bruce C. Garrett Festschrift Received: September 28, 2015 Revised: November 16, 2015 Published: November 17, 2015 1793

DOI: 10.1021/acs.jpcb.5b09466 J. Phys. Chem. B 2016, 120, 1793−1804

Article

The Journal of Physical Chemistry B

Figure 1. Potential energy surface scans using different methods for excess proton in (a) water dimer, (b) linear tetramer, and (c) Eigen-like conformations. The numbers 2.2, 2.4, 2.6, and 2.8 represent the different oxygen−oxygen distances dOO in angstroms.

combination of different distinct bonding topologies, which are the “basis states” |i⟩

along the (slow) solvent coordinate is similar before and after a proton transfer reaction so that the transition can happen in the (fast) proton transfer coordinate. Some earlier gas phase studies have also shown that the formation of such a presolvated hydronium ion is the driving force for the dissociative proton transfer reactions in hydrated organic ions.47,48 However, as noted in our previous work, the distance between the hydronium ion and the hydrogen-donating water neighbor can be quite different in the gas and liquid phases.40,49 In the present work, we also investigate the effect of anharmonicity in the water O−H stretching modes on the structural and dynamical properties of the hydrated excess proton, especially on the vibrational spectrum, by using an anharmonic version of the SPC/Fw water model developed by Park et al.29 The outline of the remaining sections of this paper is as follows: In section II, we briefly discuss the method, parametrization of the two new MS-EVB models, and the simulation details. Section III then contains the main results of the present work, and conclusions are given in section IV.

N

|Ψ⟩ =

∑ ci|i⟩

(1)

i=1

where the ci’s are the EVB coefficients and N is the total number of basis states. The ci’s are obtained by solving the eigenvalue−eigenfunction problem using the Hamiltonian matrix H consisting of matrix elements hij = ⟨i|H|j⟩, such that Hc = E0c

(2)

Here, E0 is the lowest eigenvalue and the ground state energy, which is a function of the given configuration of the nuclei, and c is the associated eigenvector whose components are the EVB coefficients in eq 1. The diagonal terms hii, which contain terms that are, e.g., (though not necessarily36) described by classical force fields and are expressed as + + hii = VHintra 3O

NH2O

∑ k=1

II. METHODS Theoretical Framework of the MS-EVB Model. It has been well-established that the delocalized nature of the net positive charge defect associated with the hydrated excess proton can be described by the MS-EVB approach in an efficient and accurate way.5,9,26,29 The full description of the MS-EVB methodology can be found in earlier papers.5,26,29 In this paper, we will summarize a few key features and then the present improvements of the MS-EVB model. In the MS-EVB framework, for a given nuclear configuration, the ground state of the system |Ψ⟩ is represented by a linear

k VHintra, + 2O

NH2O

∑ k=1

NH2O k VHinter, + O+,H O 3

2

kk ′ ∑ VHinter, O 2

k 0.5 Å. (b, right panel) The twodimensional δ-dependent RDFs for the region of r = 1.5−2.05 Å. The aMS-EVB 3.2 model also shows similar results.

the ions as in MS-EVB 3.0. Thus, the new MS-EVB models should be more easily applicable to a wider variety of counterions. Hydrated Excess Proton Self-Diffusion Constant. The self-diffusion constant of the hydrated excess proton is calculated from the mean squared displacement (MSD) of the hydrated proton center of excess charge (CEC) defined as26,29

insightful information about the presence of density near 2.0 Å in the O*−Hw RDF. To gauge the influence of the fourth water on the proton transfer events, we studied the conditional O*−Hw RDFs for different values of the proton transfer coordinate δ, defined earlier. Figure 6a clearly shows the presence of small peak around 2.0 Å in the case of the smaller δ, when a proton transfer reaction (δ = 0) can more readily happen. As we mentioned earlier, the presolvation of the hydronium is correlated with the proton transfer events, in which the two water molecules that are solvating the excess proton tend to have a similar solvation environment. When a proton transfer is about to happen (small δ), the presence of the fourth water is observed around 2.0 Å (Figure 6b). Ion Pairing with a Counterion. The interaction of a hydrated proton with a negative counterion is also an important aspect of having a good model. The MS-EVB 3.0 model has a tendency to form strong ion pairing with its counterion. Figure 7 shows one such example with the chloride counterion. To rectify this, one specific LJ interaction between the hydronium hydrogen and the counterion is usually incorporated to resolve this problem.17,73 On the other hand, by using only the standard mixing rules with the chloride ion, the two present MS-EVB 3.2 models do not have such overbinding problems (Figure 7) and do not need a specific LJ interaction between

N

rCEC =

i ∑ ci2 r COC

(4)

i=1

riCOC

where denotes the position of the center-of-charge of hydronium in the ith EVB state. The self-diffusion constant is then obtained from the MSD of the CEC using the standard Einstein relation DCEC = lim

t →∞

⟨|rCEC(t ) − rCEC(0)|2 ⟩ 6t

(5)

The MSDs for the MS-EVB 3.2 and aMS-EVB 3.2 models are shown in Figure 8. The observed diffusion constants for the

Figure 8. Mean squared displacement of the hydrated excess proton center of excess charge as a function of time.

new models are 0.37 ± 0.01 Å2/ps for MS-EVB 3.2 and 0.36 ± 0.01 Å2/ps for aMS-EVB 3.2, whereas that of MS-EVB 3.0 is 0.29 ± 0.03 Å2/ps. The (∼30%) increase in the self-diffusion constant value can be qualitatively justified from the proton transfer barrier obtained from Figure 3. It is also important that the weak interaction of the fourth water molecules with the

Figure 7. Comparison of RDFs between hydronium oxygen O* and the chloride ion (Cl−) for the different MS-EVB models. The overbinding between ions in the MS-EVB 3.0 model is reflected by the large first peak intensity around 3.0 Å. Note the greatly reduced peak intensity in the new models. 1799


Article


this study, D0/DPBC. The correction factors are 1.18, 1.19, 1.22, and 1.24 for MS-EVB 3.2, aMS-EVB 3.2, SPC/Fw, and aSPC/ Fw, respectively. As shown in Figure 10, even though these correction factors are greater than 1, the ratio of the CEC selfdiffusion constant to that of water remains roughly constant as 1/L goes to 0. For the MS-EVB 3.2 and aMS-EVB 3.2 models, combining both the finite size effect correction and the estimated CMD error yields estimated excess proton self-diffusion constants of ∼0.7 Å2/ps for both new MS-EVB models, which is quite close to the experimental result (within a factor of ∼1.3). This level of agreement is obtained without the neglect of the essential underlying physics of the hydrated proton, something that is often done in simpler and more ad hoc reactive MD models for the hydrated excess proton that claim (or are fit to achieve) better agreement with the experimental self-diffusion constant. The remaining disagreement of the present MS-EVB models with experimental results may be due to a small underestimation of the contribution of concerted proton hops,38 although this is speculation at this point and the exact role of those concerted effects remains unclear and dependent on both the specific DFT used in the AIMD simulation as well as the particular AIMD trajectory that is analyzed.40 Linear Infrared Spectrum of Excess Proton in Water. Building on earlier work with the earlier MS-EVB models,77,78 the infrared (IR) spectra of the protonated systems with the MS-EVB 3.2 and aMS-EVB 3.2 models were studied by simulating 2 excess protons and 2 chloride ions in a box of 256 water molecules (∼0.43 M acid concentration). The linear IR response can be obtained from the Fourier transform of the dipole autocorrelation function:

hydronium oxygen introduces a presolvated environment, which facilitates the proton transport.40 This effect is further quantified by studying the self-diffusion constant in absence of the extra LJ interaction between hydronium oxygen and water hydrogen. The calculated CEC self-diffusion constants in absence of the fourth water are 0.30 Å2/ps for MS-EVB 3.2 and 0.29 Å2/ps for aMS-EVB 3.2 (shown in the Supporting Information), which are similar to that of MS-EVB 3.0. The observed self-diffusion constant for aMS-EVB 3.2 is slightly lower than that of aMS-EVB 3.0 developed by Park et al.37 This might be due to the slightly higher proton transfer barrier observed in the present anharmonic model. For both of the present models, the observed self-diffusion constant is also smaller than the experimentally measured proton self-diffusion constant value of 0.94 ± 0.01 Å2/ps.74 This is partially due to the fact that the present models are developed specifically for classical trajectories, and it has been shown earlier that the inclusion of nuclear quantum effects leads to a ∼70% increment of the calculated proton self-diffusion constant. To estimate the influence of the nuclear quantum effects in these new models, the self-diffusion constant has also been calculated using the CMD method. The CMD proton self-diffusion coefficient for MS-EVB 3.2 is 0.55 ± 0.06 Å2/ps and that of aMS-EVB 3.2 is 0.51 ± 0.06 Å2/ps. The remaining difference may come from an underestimated contribution to the quantum tunneling by CMD, which can cause underestimation of the true quantum diffusion constant by up to a factor of 1.3.66,67 To estimate the finite size effect on the self-diffusion constants of both the excess proton CEC and the liquid water, a series of MD simulations was carried out with different box sizes.75 The self-diffusion constants as a function of inverse box length are shown in Figure 9. The finite size corrected

ℏβ 2 ω 2π

I(ω) =

∫ dt e−iωt ⟨μ(0)·μ(t )⟩

(6)

Here ℏ is Planck’s constant, β = 1/kBT, and μ(t) is the system dipole moment at time t. For faster convergence behavior, we instead use the autocorrelation of the time derivative of the dipole moment.77,79 ℏβ 2π

I(ω) =

∫ dt e−iωt ⟨μ̇ (0)·μ̇ (t )⟩

(7)

Here the time derivative of the dipole moment within the point charge approximation is computed as N atom

μ̇(t ) =

∑ i=1

[qi̇ (t )ri(t ) + qi(t )ri̇ (t )]

(8)

where ri(t) and ṙ(t) are the nuclear coordinates and velocity of ith atom at time t, respectively. The charge defect delocalization of the hydrated ion is taken into account by computing the state average charge using the following expression

Figure 9. Self-diffusion coefficients of CEC and water as a function of inverse simulation box length (1/L).

NEVB

diffusion constants are obtained by interpolating the data in the limit of infinite box size, 1/L = 0. The corrected diffusion constants of CEC for MS-EVB 3.2 and aMS-EVB 3.2 models are 0.435 and 0.429 Å2/ps, respectively. The size corrected diffusion constants for SPC/Fw and aSPC/Fw water models are 0.284 and 0.289 Å2/ps, respectively. Similar to results from the previous studies75,76 the system size effect on the diffusion constants shows a systematic increase of about 20% from 1/L 0.051 to 0 Å−1. To obtain an estimate of the correction factor, we calculated the ratio of diffusion constant at infinite box length (D0) to the diffusion constant from the box size used in

qi =

∑ ck2qki k=1

(9)

The experimental difference spectrum is obtained by subtracting the pure water spectrum from the 0.4 M HCl spectrum obtained using attenuated total reflection (ATR) IR spectroscopy. In Figure S4 (Supporting Information), we show the difference spectra obtained from the classical MS-EVB simulations and experiment. The experimentally observed spectrum for the 0.4 M HCl shows prominent features at 1760 and 3000 cm−1, a broad peak at 1200 cm−1, and the “acid 1800


Article


Figure 10. IR difference spectra obtained from classical and CMD trajectories of (a, left) aMS-EVB 3.2 and (b, right) MS-EVB 3.2 simulations of 1 HCl aqueous system along with the experimental attenuated total reflection (ATR) difference spectrum. Note that in the case of CMD anharmonic spectra the peak position is in better agreement with the experimental result.

continuum” absorption over the 1000−3400 cm−1 range. In spite of using classical trajectories in calculating these MS-EVB spectra, several experimental features are captured by both new excess proton models. The peaks around 1200 and 1760 cm−1 are observed in a slightly red-shifted position. The acid continuum band spanning 1000−3400 cm−1 is also reproduced, although the high-frequency region, especially above 2500 cm −1 for both MS-EVB models, is lacking the precise intensity trends observed from experiment. The anharmonic water used in the aMS-EVB 3.2 model leads to a better description of the spectroscopic response; however, the high-frequency dip in the difference spectrum is too blue-shifted relative to experiment due to the missing quantum red shift of the water O−H vibrations in the classical MS-EVB MD simulations. Although the MS-EVB formalism is able to include electronic polarization close to the excess proton, most of the spectroscopic responses are influenced by the contribution from the (classical) water, which forms the bulk of the material. It is expected that CMD simulations, although very computationally demanding for the calculation of this particular spectrum, will give an accurate representation of the quantum red shift of the high-frequency range for the water O−H stretching mode.80 (Such CMD results for 300 K would be free of the widely mentioned and often overly emphasized “curvature problem” of CMD, which is present mainly at low temperatures.80,81) The IR absorption coefficients from CMD simulations for both proton models were therefore computed using the centroid dipole moment time derivative correlation function (see Figure 10):82,83 α(ω) =

4π 2ℏβ ICMD(ω) 3Vcn(ω)

256 water molecules only. From Figure 10a,b, it is evident that the different proton and water models show a significant quantum red shift of the high-frequency O−H stretch region. Particularly in the case of the anharmonic proton model, the O−H stretch region gives a more accurate peak position with respect to the experimental results (see Figure 10a). The missing intensity in Figure 10 and Figure S4 in the highfrequency range from 2500 to 3000 cm−1 may possibly be attributed to the lack of polarizability of the underlying water models.84 In addition, the current excess proton MS-EVB models might not be able to completely describe the charge transfer behavior of the electronic wave functions for different asymmetric modes of the hydrated proton,32 nor the full nonlinear dipole moment surface (electrical anharmonicity), both of which can be significant for the spectroscopic properties.

IV. CONCLUSIONS Two new MS-EVB hydrated excess proton models, MS-EVB 3.2 and aMS-EVB 3.2, have been developed and presented in this work on the basis of the presolvation concept of the hydronium cation. The hydrated proton solvation structure and dynamical properties of the new models show similar equilibrium properties but with improved dynamical properties over the MS-EVB 3.0 model. The prevalent nature of the Eigen over Zundel population is observed from the hydrated proton solvation free energy. The proton transfer free energy barrier is further reduced, and the presolvation environment around the hydronium ion through a weak interaction with a fourth water molecule is incorporated; its correlation with proton transfer is shown in these models. The hydrated excess proton selfdiffusion constant is improved, and the inclusion of nuclear quantum effects and finite size simulation box corrections is shown to further increase the excess proton diffusion rate to bring it into much closer agreement with experiment. Both new hydrated excess proton MS-EVB models also qualitatively capture several spectroscopic responses observed in experiments. The use of the anharmonic water model in the aMSEVB 3.2 gives an increased intensity in the low-frequency range, although in the high-frequency range the spectra still lack some intensity for both of the new models compared to experiment. The use of CMD gives rise to the quantum red shift for both proton models in the high-frequency range for the water O−H stretching mode, improving the peak position with respect to the experimental result. Future research will therefore focus on

(10)

Here V is the volume of the simulation box, c is the speed of light in vacuum, n(ω) is the refraction coefficient that was set to unity in our calculations, and ICMD(ω) is obtained from the Fourier transform of dipole moment derivative correlation function ICMD(ω) =

1 2π

̇ ̇ (0) ·μCMD (t )⟩ ∫ dt e−iωt ⟨μCMD

(11)

As previously mentioned, the CMD simulations are computationally demanding and become more so with increasing number of reactive centers. Thus, we have calculated all the CMD spectra for either 256 water molecules only or 1 HCl in 1801


Article


(7) Tuckerman, M. E.; Marx, D.; Parrinello, M. The Nature and Transport Mechanism of Hydrated Hydroxide Ions in Aqueous Solution. Nature 2002, 417, 925−9. (8) Marx, D. Proton Transfer 200 Years after Von Grotthuss: Insights from Ab Initio Simulations. ChemPhysChem 2006, 7, 1848−1870. (9) Knight, C.; Voth, G. A. The Curious Case of the Hydrated Proton. Acc. Chem. Res. 2012, 45, 101−109. (10) Wraight, C. A. Chance and Design - Proton Transfer in Water, Channels and Bioenergetic Proteins. Biochim. Biophys. Acta, Bioenerg. 2006, 1757, 886−912. (11) Decoursey, T. E. Voltage-Gated Proton Channels and Other Proton Transfer Pathways. Physiol. Rev. 2003, 83, 475−579. (12) Cukierman, S. Et Tu, Grotthuss! And Other Unfinished Stories. Biochim. Biophys. Acta, Bioenerg. 2006, 1757, 876−885. (13) Swanson, J. M. J.; Maupin, C. M.; Chen, H.; Petersen, M. K.; Xu, J.; Wu, Y.; Voth, G. A. Proton Solvation and Transport in Aqueous and Biomolecular Systems: Insights from Computer Simulations. J. Phys. Chem. B 2007, 111, 4300−4314. (14) Kreuer, K. D.; Paddison, S. J.; Spohr, E.; Schuster, M. Transport in Proton Conductors for Fuel-Cell Applications: Simulations, Elementary Reactions, and Phenomenology. Chem. Rev. 2004, 104, 4637−4678. (15) Jorn, R.; Savage, J.; Voth, G. A. Proton Conduction in Exchange Membranes across Multiple Length Scales. Acc. Chem. Res. 2012, 45, 2002−2010. (16) Tse, Y.-L. S.; Herring, A. M.; Kim, K.; Voth, G. A. Molecular Dynamics Simulations of Proton Transport in 3M and Nafion Perfluorosulfonic Acid Membranes. J. Phys. Chem. C 2013, 117, 8079− 8091. (17) Savage, J.; Tse, Y.-L. S.; Voth, G. A. Proton Transport Mechanism of Perfluorosulfonic Acid Membranes. J. Phys. Chem. C 2014, 118, 17436−17445. (18) Garczarek, F.; Gerwert, K. Functional Waters in Intraprotein Proton Transfer Monitored by FTIR Difference Spectroscopy. Nature 2006, 439, 109−112. (19) Roberts, S. T.; Petersen, P. B.; Ramasesha, K.; Tokmakoff, A.; Ufimtsev, I. S.; Martinez, T. J. Observation of a Zundel-Like Transition State During Proton Transfer in Aqueous Hydroxide Solutions. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 15154−15159. (20) Roberts, S. T.; Ramasesha, K.; Petersen, P. B.; Mandal, A.; Tokmakoff, A. Proton Transfer in Concentrated Aqueous Hydroxide Visualized Using Ultrafast Infrared Spectroscopy. J. Phys. Chem. A 2011, 115, 3957−3972. (21) Mandal, A.; Ramasesha, K.; De Marco, L.; Tokmakoff, A. Collective Vibrations of Water-Solvated Hydroxide Ions Investigated with Broadband 2dir Spectroscopy. J. Chem. Phys. 2014, 140, 204508. (22) Fayer, M. D. Ultrafast Infrared Vibrational Spectroscopy; CRC Press: Boca Raton, FL, 2013. (23) Fournier, J. A.; Johnson, C. J.; Wolke, C. T.; Weddle, G. H.; Wolk, A. B.; Johnson, M. A. Vibrational Spectral Signature of the Proton Defect in the Three-Dimensional H+(H2o)21 Cluster. Science 2014, 344, 1009−1012. (24) Thämer, M.; De Marco, L.; Ramasesha, K.; Mandal, A.; Tokmakoff, A. Ultrafast 2d IR Spectroscopy of the Excess Proton in Liquid Water. Science 2015, 350, 78−82. (25) Tuckerman, M. E.; Marx, D.; Klein, M. L.; Parrinello, M. On the Quantum Nature of the Shared Proton in Hydrogen Bonds. Science 1997, 275, 817−820. (26) Day, T. J. F.; Soudackov, A. V.; Cuma, M.; Schmitt, U. W.; Voth, G. A. A Second Generation Multistate Empirical Valence Bond Model for Proton Transport in Aqueous Systems. J. Chem. Phys. 2002, 117, 5839−5849. (27) Izvekov, S.; Voth, G. A. Ab Initio Molecular-Dynamics Simulation of Aqueous Proton Solvation and Transport Revisited. J. Chem. Phys. 2005, 123, 044505. (28) Lapid, H.; Agmon, N.; Petersen, M. K.; Voth, G. A. A BondOrder Analysis of the Mechanism for Hydrated Proton Mobility in Liquid Water. J. Chem. Phys. 2005, 122, 014506.

achieving spectroscopic accuracy through additional improvements of these new MS-EVB models, while not diminishing their several already optimal or near-optimal properties.

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.5b09466. PESs and comparison of new models RDFs with MSEVB 3.0 RDFs and AIMD RDFs, MSDs in the absence of the fourth water interaction, and classical IR difference spectra in case of 2HCl aqueous system (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*Phone: 773-702-7250. Fax: 773-702-0805. E-mail: gavoth@ uchicago.edu. Present Address

§ Department of Chemistry, Chinese University of Hong Kong, Shatin, New Territories, Hong Kong.

Author Contributions

R.B. and Y.-L.S.T. contributed equally. The project was conceived by and the manuscript was written through the contributions of all authors. R.B. and Y.-L.S.T. carried out the calculations in the research. All authors analyzed the data. All authors have given approval to the final version of the manuscript. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This research was supported in part by the Department of Energy (DOE), Office of Basic Energy Sciences (BES), Division of Chemical Sciences, Geosciences, and Biosciences, through Grants DE-SC0005418 to G.A.V. and DE-SC0014305 to A.T., and by the National Science Foundation (NSF), Grant CHE-1465248 to G.A.V. The Croucher Foundation is acknowledged for a postdoctoral research fellowship to Y.-L.S.T. The computational resources in this work were provided in part by a grant of computer time from the U.S. Department of Defense (DOD) High Performance Computing Modernization Program at the Engineer Research and Development Center (ERDC) and Navy DOD Supercomputing Resource Centers and in part by the University of Chicago Research Computing Center (RCC).

■

REFERENCES

(1) von Grotthuss, C. J. T. Sur la Décomposition de l’Eau et des Corps Qu’Elle Tient en Dissolution à l’Aide de l’É ́ lectricité Galvanique. Ann. Chim. 1806, 58, 54−73. (2) Agmon, N. The Grotthuss Mechanism. Chem. Phys. Lett. 1995, 244, 456−462. (3) Tuckerman, M.; Laasonen, K.; Sprik, M.; Parrinello, M. Ab-Initio Molecular-Dynamics Simulation of the Solvation and Transport of Hydronium and Hydroxyl Ions in Water. J. Chem. Phys. 1995, 103, 150−161. (4) Lobaugh, J.; Voth, G. A. The Quantum Dynamics of an Excess Proton in Water. J. Chem. Phys. 1996, 104, 2056−2069. (5) Schmitt, U. W.; Voth, G. A. The Computer Simulation of Proton Transport in Water. J. Chem. Phys. 1999, 111, 9361−9381. (6) Marx, D.; Tuckerman, M. E.; Hutter, J.; Parrinello, M. The Nature of the Hydrated Excess Proton in Water. Nature 1999, 397, 601−604. 1802


Article

The Journal of Physical Chemistry B (29) Wu, Y. J.; Chen, H. N.; Wang, F.; Paesani, F.; Voth, G. A. An Improved Multistate Empirical Valence Bond Model for Aqueous Proton Solvation and Transport. J. Phys. Chem. B 2008, 112, 467−482. (30) Markovitch, O.; Chen, H.; Izvekov, S.; Paesani, F.; Voth, G. A.; Agmon, N. Special Pair Dance and Partner Selection: Elementary Steps in Proton Transport in Liquid Water. J. Phys. Chem. B 2008, 112, 9456−9466. (31) Berkelbach, T. C.; Lee, H. S.; Tuckerman, M. E. Concerted Hydrogen-Bond Dynamics in the Transport Mechanism of the Hydrated Proton: A First-Principles Molecular Dynamics Study. Phys. Rev. Lett. 2009, 103, 238302. (32) Swanson, J. M. J.; Simons, J. Role of Charge Transfer in the Structure and Dynamics of the Hydrated Proton. J. Phys. Chem. B 2009, 113, 5149−5161. (33) Marx, D.; Chandra, A.; Tuckerman, M. E. Aqueous Basic Solutions: Hydroxide Solvation, Structural Diffusion, and Comparison to the Hydrated Proton. Chem. Rev. 2010, 110, 2174−2216. (34) Knight, C.; Maupin, C. M.; Izvekov, S.; Voth, G. A. Defining Condensed Phase Reactive Force Fields from Ab Initio Molecular Dynamics Simulations: The Case of the Hydrated Excess Proton. J. Chem. Theory Comput. 2010, 6, 3223−3232. (35) Tuckerman, M. E.; Chandra, A.; Marx, D. A Statistical Mechanical Theory of Proton Transport Kinetics in HydrogenBonded Networks Based on Population Correlation Functions with Applications to Acids and Bases. J. Chem. Phys. 2010, 133, 124108. (36) Knight, C.; Lindberg, G. E.; Voth, G. A. Multiscale Reactive Molecular Dynamics. J. Chem. Phys. 2012, 137, 22A525. (37) Park, K.; Lin, W.; Paesani, F. A Refined Ms-Evb Model for Proton Transport in Aqueous Environments. J. Phys. Chem. B 2012, 116, 343−352. (38) Hassanali, A.; Giberti, F.; Cuny, J.; Kuhne, T. D.; Parrinello, M. Proton Transfer through the Water Gossamer. Proc. Natl. Acad. Sci. U. S. A. 2013, 110, 13723−13728. (39) Hassanali, A. A.; Giberti, F.; Sosso, G. C.; Parrinello, M. The Role of the Umbrella Inversion Mode in Proton Diffusion. Chem. Phys. Lett. 2014, 599, 133−138. (40) Tse, Y. L.; Knight, C.; Voth, G. A. An Analysis of Hydrated Proton Diffusion in Ab Initio Molecular Dynamics. J. Chem. Phys. 2015, 142, 014104. (41) Peng, Y.; Swanson, J. M. J.; Kang, S.-g.; Zhou, R.; Voth, G. A. Hydrated Excess Protons Can Create Their Own Water Wires. J. Phys. Chem. B 2015, 119, 9212−9218. (42) Marx, D.; Hutter, J. Ab Initio Molecular Dynamics: Basic Theory and Advanced Methods; Cambridge University Press: New York, 2009. (43) Yoo, S.; Zeng, X. C.; Xantheas, S. S. On the Phase Diagram of Water with Density Functional Theory Potentials: The Melting Temperature of Ice I-H with the Perdew-Burke-Ernzerhof and BeckeLee-Yang-Parr Functionals. J. Chem. Phys. 2009, 130, 221102. (44) Yoo, S.; Xantheas, S. S. Communication: The Effect of Dispersion Corrections on the Melting Temperature of Liquid Water. J. Chem. Phys. 2011, 134, 121105. (45) Ando, K.; Hynes, J. T. Molecular Mechanism of Hcl Acid Ionization in Water: Ab Initio Potential Energy Surfaces and Monte Carlo Simulations. J. Phys. Chem. B 1997, 101, 10464−10478. (46) Ando, K.; Hynes, J. T. Acid-Base Proton Transfer and Ion Pair Formation in Solution. Adv. Chem. Phys. 1999, 110, 381−430. (47) Hamid, A. M.; Sharma, P.; El-Shall, M. S.; Hilal, R.; Elroby, S.; Aziz, S. G.; Alyoubi, A. O. Hydration of the Pyrimidine Radical Cation and Stepwise Solvation of Protonated Pyrimidine with Water, Methanol, and Acetonitrile. J. Chem. Phys. 2013, 139, 084304. (48) Ibrahim, Y. M.; Mautner, M. M. N.; Alshraeh, E. H.; El-Shall, M. S.; Scheiner, S. Stepwise Hydration of Ionized Aromatics. Energies, Structures of the Hydrated Benzene Cation, and the Mechanism of Deprotonation Reactions. J. Am. Chem. Soc. 2005, 127, 7053−7064. (49) Jagoda-Cwiklik, B.; Cwiklik, L.; Jungwirth, P. Behavior of the Eigen Form of Hydronium at the Air/Water Interface. J. Phys. Chem. A 2011, 115, 5881−5886.

(50) Wu, Y.; Tepper, H. L.; Voth, G. A. Flexible Simple Point-Charge Water Model with Improved Liquid-State Properties. J. Chem. Phys. 2006, 124, 024503. (51) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic-Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (52) Lee, C. T.; Yang, W. T.; Parr, R. G. Development of the ColleSalvetti Correlation-Energy Formula into a Functional of the ElectronDensity. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (53) Grimme, S. Accurate Description of Van Der Waals Complexes by Density Functional Theory Including Empirical Corrections. J. Comput. Chem. 2004, 25, 1463−1473. (54) 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. (55) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09; Gaussian, Inc.: Wallingford, CT, 2009. (56) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1−19. (57) Dang, L. X. Mechanism and Thermodynamics of Ion Selectivity in Aqueous-Solutions of 18-Crown-6 Ether - a Molecular-Dynamics Study. J. Am. Chem. Soc. 1995, 117, 6954−6960. (58) Martyna, G. J.; Klein, M. L.; Tuckerman, M. Nose-Hoover Chains - the Canonical Ensemble Via Continuous Dynamics. J. Chem. Phys. 1992, 97, 2635−2643. (59) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, U.K., 1989. (60) Cao, J. S.; Voth, G. A. A New Perspective on Quantum TimeCorrelation Functions. J. Chem. Phys. 1993, 99, 10070−10073. (61) Cao, J. S.; Voth, G. A. The Formulation of Quantum-Statistical Mechanics Based on the Feynman Path Centroid Density. 2. Dynamical Properties. J. Chem. Phys. 1994, 100, 5106−5117. (62) Cao, J. S.; Voth, G. A. The Formulation of Quantum-Statistical Mechanics Based on the Feynman Path Centroid Density. 3. PhaseSpace Formalism and Analysis of Centroid Molecular-Dynamics. J. Chem. Phys. 1994, 101, 6157−6167. (63) Cao, J. S.; Voth, G. A. The Formulation of Quantum-Statistical Mechanics Based on the Feynman Path Centroid Density. 4. Algorithms for Centroid Molecular-Dynamics. J. Chem. Phys. 1994, 101, 6168−6183. (64) Jang, S.; Voth, G. A. Path Integral Centroid Variables and the Formulation of Their Exact Real Time Dynamics. J. Chem. Phys. 1999, 111, 2357−2370. (65) Jang, S.; Voth, G. A. A Derivation of Centroid Molecular Dynamics and Other Approximate Time Evolution Methods for Path Integral Centroid Variables. J. Chem. Phys. 1999, 111, 2371−2384. (66) Hone, T. D.; Voth, G. A. A Centroid Molecular Dynamics Study of Liquid Para-Hydrogen and Ortho-Deuterium. J. Chem. Phys. 2004, 121, 6412−6422. (67) Hone, T. D.; Rossky, P. J.; Voth, G. A. A Comparative Study of Imaginary Time Path Integral Based Methods for Quantum Dynamics. J. Chem. Phys. 2006, 124, 154103. (68) Paesani, F.; Zhang, W.; Case, D. A.; Cheatham, T. E.; Voth, G. A. An Accurate and Simple Quantum Model for Liquid Water. J. Chem. Phys. 2006, 125, 184507. (69) Eigen, M. Proton Transfer Acid-Base Catalysis + Enzymatic Hydrolysis. Part I: Elementary Processes. Angew. Chem., Int. Ed. Engl. 1964, 3, 1−19. (70) Zundel, G.; Metzger, H. Energiebänder Der Tunnelnden Ü berschuß-Protonen in Flüssigen Säuren. Eine Ir-Spektroskopische Untersuchung Der Natur Der Gruppierungen H5o2+. Z. Phys. Chem. 1968, 58, 225−245. (71) Zundel, G. Hydrogen Bonds with Large Proton Polarizability and Proton Transfer Processes in Electrochemistry and Biology. Adv. Chem. Phys. 1999, 111, 1−217. 1803


Article

The Journal of Physical Chemistry B (72) Botti, A.; Bruni, F.; Imberti, S.; Ricci, M. A.; Soper, A. K. Ions in Water: The Microscopic Structure of a Concentrated HCl Solution. J. Chem. Phys. 2004, 121, 7840−7848. (73) Wang, F.; Izvekov, S.; Voth, G. A. Unusual ″Amphiphilic″ Association of Hydrated Protons in Strong Acid Solution. J. Am. Chem. Soc. 2008, 130, 3120−3126. (74) Roberts, N. K.; Northey, H. L. Proton and Deuteron Mobility in Normal and Heavy Water Solutions of Electrolytes. J. Chem. Soc., Faraday Trans. 1 1974, 70, 253−262. (75) Yeh, I. C.; Hummer, G. System-Size Dependence of Diffusion Coefficients and Viscosities from Molecular Dynamics Simulations with Periodic Boundary Conditions. J. Phys. Chem. B 2004, 108, 15873−15879. (76) Figueirido, F.; Delbuono, G. S.; Levy, R. M. On Finite-Size Effects in Computer-Simulations Using the Ewald Potential. J. Chem. Phys. 1995, 103, 6133−6142. (77) Kim, J.; Schmitt, U. W.; Gruetzmacher, J. A.; Voth, G. A.; Scherer, N. E. The Vibrational Spectrum of the Hydrated Proton: Comparison of Experiment, Simulation, and Normal Mode Analysis. J. Chem. Phys. 2002, 116, 737−746. (78) Xu, J.; Zhang, Y.; Voth, G. A. Infrared Spectrum of the Hydrated Proton in Water. J. Phys. Chem. Lett. 2011, 2, 81−6. (79) Thomas, M.; Brehm, M.; Fligg, R.; Vohringer, P.; Kirchner, B. Computing Vibrational Spectra from Ab Initio Molecular Dynamics. Phys. Chem. Chem. Phys. 2013, 15, 6608−6622. (80) Paesani, F.; Voth, G. A. A Quantitative Assessment of the Accuracy of Centroid Molecular Dynamics for the Calculation of the Infrared Spectrum of Liquid Water. J. Chem. Phys. 2010, 132, 014105. (81) Rossi, M.; Liu, H.; Paesani, F.; Bowman, J.; Ceriotti, M. Communication: On the Consistency of Approximate Quantum Dynamics Simulation Methods for Vibrational Spectra in the Condensed Phase. J. Chem. Phys. 2014, 141, 181101. (82) Witt, A.; Ivanov, S. D.; Shiga, M.; Forbert, H.; Marx, D. On the Applicability of Centroid and Ring Polymer Path Integral Molecular Dynamics for Vibrational Spectroscopy. J. Chem. Phys. 2009, 130, 194510. (83) Liu, J.; Miller, W. H.; Paesani, F.; Zhang, W.; Case, D. A. Quantum Dynamical Effects in Liquid Water: A Semiclassical Study on the Diffusion and the Infrared Absorption Spectrum. J. Chem. Phys. 2009, 131, 164509. (84) Brancato, G.; Tuckerman, M. E. A Polarizable Multistate Empirical Valence Bond Model for Proton Transport in Aqueous Solution. J. Chem. Phys. 2005, 122, 224507.

1804


Role of Presolvation and Anharmonicity in Aqueous Phase Hydrated

Recommend Documents