Quantum effects in ice Ih

Quantum and classical simulations are carried out on ice Ih over a range of temperatures utilizing the TIP4P water model. The rigid-body centroid molecular dynamics method employed allows for the investigation of equilibrium and dynamical properties of the quantum system. The impact of quantization on the local structure, as measured by the radial and spatial distribution functions, as well as the energy is presented. The effects of quantization on the lattice vibrations, associated with the molecular translations and librations, are also reported. Comparison of quantum and classical simulation results indicates that shifts in the average potential energy are equivalent to rising the temperature about 80 K and are therefore non-negligible. The energy shifts due to quantization and the quantum mechanical uncertainties observed in ice are smaller than the values previously reported for liquid water. Additionally, we carry out a comparative study of melting in our classical and quantum simulations and show that there are significant differences between classical and quantum ice.     

Isotope effects in ice Ih: a path-integral simulation

Ice Ih has been studied by path-integral molecular dynamics simulations, using the effective q-TIP4P/F potential model for flexible water. This has allowed us to analyze finite-temperature quantum effects in this solid phase from 25 to 300 K at ambient pressure. Among these effects we find a negative thermal expansion of ice at low temperatures, which does not appear in classical molecular dynamics simulations. The compressibility derived from volume fluctuations gives results in line with experimental data. We have analyzed isotope effects in ice Ih by considering normal, heavy, and tritiated water. In particular, we studied the effect of changing the isotopic mass of hydrogen on the kinetic energy and atomic delocalization in the crystal as well as on structural properties such as interatomic distances and molar volume. For D(2)O ice Ih at 100 K we obtained a decrease in molar volume and intramolecular O-H distance of 0.6% and 0.4%, respectively, as compared to H(2)O ice.     

Theoretical study on isotope and temperature effect in hydronium ion using ab initio path integral simulation

We have applied the ab initio path integral molecular dynamics simulation to study hydronium ion and its isotopes, which are the simplest systems for hydrated proton and deuteron. In this simulation, all the rotational and vibrational degrees of freedom are treated fully quantum mechanically, while the potential energies of the respective atomic configurations are calculated "on the fly" using ab initio quantum chemical approach. With the careful treatment of the ab initio electronic structure calculation by relevant choices in electron correlation level and basis set, this scheme is theoretically quite rigorous except for Born-Oppenheimer approximation. This accurate calculation allows a close insight into the structural shifts for the isotopes of hydronium ion by taking account of both quantum mechanical and thermal effects. In fact, the calculation is shown to be successful to quantitatively extract the geometrical isotope effect with respect to the Walden inversion. It is also shown that this leads to the isotope effect on the electronic structure as well as the thermochemical properties.     

Quantum simulation of a hydrated noradrenaline analog with the torsional path integral method

An extended version of the torsional path integral Monte Carlo (TPIMC) method is presented and shown to be useful for studying the conformation of flexible molecules in solvated clusters. The new technique is applied to the hydrated clusters of the 2-amino-1-phenyl-ethanol (APE) molecule. APE + nH2O clusters with n = 0-4 are studied at 100 and 300 K using both classical and quantum simulations. Only at the lower temperature is the hydration number n found to impact the conformational distribution of the APE molecule. This is shown to be a result of the temperature-dependent balance between the internal energy and entropy contributions to the relative conformer free energies. Furthermore, at 100 K, large quantum effects are observed in the calculated conformer populations. A particularly large quantum shift of 30% of the total population is calculated for the APE + 2H2O cluster, which is explained in terms of the relative zero point energy of the lowest-energy hydrated structures for this cluster. Finally, qualitative agreement is found between the reported calculations and recent spectroscopy experiments on the hydrated clusters of APE, including an entropically driven preference for the formation of AG-type hydrated structures and the formation of a water "droplet" in the APE + 4H2O cluster.     

Melting temperature of ice Ih calculated from coexisting solid-liquid phases

We carried out molecular-dynamics simulations by using the two-phase coexistence method with the constant pressure, particle number, and enthalpy ensemble to compute the melting temperature of proton-disordered hexagonal ice I(h) at 1-bar pressure. Four models of water were considered, including the widely used TIP4P [W. L. Jorgensen, J. Chandrasekha, J. D. Madura, R. W. Impey, and M. L. Klein, J. Chem. Phys.79, 926 (1983)] and TIP5P [M. W. Mahoney and W. L. Jorgensen J. Chem. Phys.112, 8910 (2000)] models, as well as recently improved TIP4P and TIP5P models for use with Ewald techniques-the TIP4P-Ew [W. Horn, W. C. Swope, J. W. Pitera, J. C. Madura, T. J. Dick, G. L. Hura, and T. Head-Gordon, J. Chem. Phys.120, 9665 (2004)] and TIP5P-Ew [S. W. Rick, J. Chem. Phys.120, 6085 (2004)] models. The calculated melting temperature at 1 bar is T(m) = 229 +/- 1 K for the TIP4P and T(m) = 272.0 +/- 0.6 K for the TIP5P ice I(h), both are consistent with previous simulations based on free-energy methods. For the TIP4P-Ew and TIP5P-Ew models, the calculated melting temperature is T(m) = 257.0 +/- 1.1 K and T(m) = 253.9 +/- 1.1 K, respectively.     

Path integral calculation of free energies: quantum effects on the melting temperature of neon

The path integral formulation has been combined with several methods to determine free energies of quantum many-body systems, such as adiabatic switching and reversible scaling. These techniques are alternatives to the standard thermodynamic integration method. A quantum Einstein crystal is used as a model to demonstrate the accuracy and reliability of these free energy methods in quantum simulations. Our main interest focuses on the calculation of the melting temperature of Ne at ambient pressure, taking into account quantum effects in the atomic dynamics. The free energy of the solid was calculated by considering a quantum Einstein crystal as reference state, while for the liquid, the reference state was defined by the classical limit of the fluid. Our findings indicate that, while quantum effects in the melting temperature of this system are small, they still amount to about 6% of the melting temperature, and are therefore not negligible. The particle density as well as the melting enthalpy and entropy of the solid and liquid phases at coexistence is compared to results obtained in the classical limit and also to available experimental data.     

The melting point of ice Ih for common water models calculated from direct coexistence of the solid-liquid interface

In this work we present an implementation for the calculation of the melting point of ice I(h) from direct coexistence of the solid-liquid interface. We use molecular dynamics simulations of boxes containing liquid water and ice in contact. The implementation is based on the analysis of the evolution of the total energy along NpT simulations at different temperatures. We report the calculation of the melting point of ice I(h) at 1 bar for seven water models: SPC/E, TIP4P, TIP4P-Ew, TIP4P/ice, TIP4P/2005, TIP5P, and TIP5P-E. The results for the melting temperature from the direct coexistence simulations of this work are in agreement (within the statistical uncertainty) with those obtained previously by us from free energy calculations. By taking into account the results of this work and those of our free energy calculations, recommended values of the melting point of ice I(h) at 1 bar for the above mentioned water models are provided.     

The melting point of hexagonal ice (Ih) is strongly dependent on the quadrupole of the water models

It is shown that quadrupolar interactions play a determinant role in the melting temperatures of common water models and that there is a simple relationship between the strength of the quadrupolar forces and the position of the negative charge; our conclusion is that acceptable predictions for the melting temperature can only be obtained when the negatively charged site is shifted from the oxygen atom towards the hydrogens.     

Quantum effects and isotope effects in chemical kinetics

The relations of quantum and kinetic isotope effects are investigated using the exact quantum corrections to the collision and activated complex theories. The latter are computed for the collinear three-atomic reaction H₂ + H → H + H₂ and the related isotopic reactions using realistic potential energy surfaces. Taking into account the bent configurations of the collision complex H-H-H gives a very good agreement between the quantum collision theory and the experimental data for the absolute values and the isotopic ratios of rate constants. Classical trajectory calculations yield considerably lower results.     

Molecular-dynamics studies of surface of ice Ih

We performed molecular dynamics calculations of surface of ice Ih in order to investigate formation mechanism of melting layer on the surface. The results showed that the vibrational amplitude of the atoms in the surface layer greatly depends on the crystal orientation, whereas that in the ice bulk is isotropic. The anisotropy of the vibration is due to a dangling motion of the free O-H bonds exist at the surface layer. The dangling motion enhances the rotational motion of the water molecules. The vibrational density of state showed a coupling between the rotational vibration and the lattice vibration of the water molecules in the surface layer. The coupling of the vibrations causes a distortion of ice lattice. Through the hydrogen-bonding network, the distortion transmits to the interior of the crystal. We conclude that the dangling motion of the free O-H bonds exist at the surface layer is one of the dominant factors governing the surface melting of ice crystal.     

Two-dimensional infrared spectroscopy of neat ice Ih

The OH stretch line shape of ice Ih exhibits distinct peaks, the assignment of which remains controversial. We address this longstanding question using two dimensional infrared (2D IR) spectroscopy of the OH stretch of H(2)O and the OD stretch of D(2)O of ice Ih at T = 80 K. The isotropic response is dominated by a 2D line shape component which does not depend on the pump pulse frequency. The decay time of the component that does depend on the pump frequency is calculated using singular value decomposition (bi-exponential decay H(2)O: 30 fs, 490 fs; D(2)O: 40 fs, 690 fs). The anisotropic contribution exhibits on-diagonal peaks, which decay on a very fast timescale (H(2)O: 85 fs; D(2)O: 65 fs), with no corresponding anisotropic cross-peaks. Both isotropic and anisotropic results indicate that randomization of excited dipoles occurs with a very rapid rate, just like in neat liquid water. We conclude that the underlying mechanism relates to the complex interplay between exciton migration and exciton-phonon coupling.     

Two-dimensional infrared spectroscopy of isotope-diluted ice Ih

We present experimental 2D IR spectra of isotope diluted ice Ih (i.e., the OH stretch mode of HOD in D(2)O and the OD stretch mode of HOD in H(2)O) at T = 80 K. The main spectral features are the extremely broad 1-2 excited state transition, much broader than the corresponding 0-1 groundstate transition, as well as the presence of quantum beats. We do not observe any inhomogeneous broadening that might be expected due to proton disorder in ice Ih. Complementary, we perform simulations in the framework of the Lippincott-Schroeder model, which qualitatively reproduce the experimental observations. We conclude that the origin of the observed line shape features is the coupling of the OH-vibrational coordinate with crystal phonons and explain the beatings as a coherent oscillation of the O···O hydrogen bond degree of freedom.     

Transient absorption of vibrationally excited ice Ih

The ultrafast dynamics of HDO:D2O ice Ih at 180 K is studied by midinfrared ultrafast pump-probe spectroscopy. The vibrational relaxation of HDO:D2O ice is observed to proceed via an intermediate state, which has a blueshifted absorption spectrum. Polarization resolved measurements reveal that the intermediate state is part of the intramolecular relaxation pathway of the HDO molecule. In addition, slow dynamics on a time scale of the order of 10-100 ps is observed, related to thermally induced collective reorganizations of the ice lattice. The transient absorption line shape is analyzed within a Lippincott-Schroeder model for the OH-stretch potential. This analysis identifies the main mechanism behind the strong spectral broadening of the v(OH)=1-->2 transition.     

Simulations of proton order and disorder in ice Ih

Computer simulations of ice Ih with different proton orientations are presented. Simulations of proton disordered ice are carried out using a Monte Carlo method which samples over proton degree of freedom, allowing for the calculation of the dielectric constant and for the examination of the degree of proton disorder. Simulations are also presented for two proton ordered structures of ice Ih, the ferroelectric Cmc2(1) structure or ice XI and the antiferroelectric Pna2(1) structure. These simulations indicate that a transition to a proton ordered phase occurs at low temperatures (below 80 K). The symmetry of the ordered phase is found to be dependent on the water potential. The stability of the two proton ordered structures is due to a balance of short-ranged interactions which tend to stabilize the Pna2(1) structure and longer-range interactions which stabilize the Cmc2(1) structure.     

HD isotope effect on the dihydrogen bond of NH+4...BeH2 by ab initio path integral molecular dynamics simulation

In order to investigate the HD isotope effect on a dihydrogen bonded cation system, we have studied NH+4...BeH2 and its isotopomers by ab initio path integral molecular dynamics. It is found that the dihydrogen bond can be exchanged by NH+(4) rotation. The deuterated isotopomer (ND+(4)...BeD(2); DD) can exchange the dihydrogen bond more easily than other isotopomers such as (NH+4...BeH2; HH). This unusual isotope effect is ascribed to the "quantum localization" which occurs when the effective energy barrier for the rotational mode becomes higher by the zero point energy of other modes. We also found that the binding energy of dihydrogen bonds for DD species is the smallest among the isotopomers.     

Quantum mechanical single molecule partition function from path integral Monte Carlo simulations

An algorithm for calculating the partition function of a molecule with the path integral Monte Carlo method is presented. Staged thermodynamic perturbation with respect to a reference harmonic potential is utilized to evaluate the ratio of partition functions. Parallel tempering and a new Monte Carlo estimator for the ratio of partition functions are implemented here to achieve well converged simulations that give an accuracy of 0.04 kcal/mol in the reported free energies. The method is applied to various test systems, including a catalytic system composed of 18 atoms. Absolute free energies calculated by this method lead to corrections as large as 2.6 kcal/mol at 300 K for some of the examples presented.     

A van der Waals density functional study of ice Ih

Density functional theory with the van der Waals density functional (vdW-DF) is used to calculate equilibrium crystal structure, binding energy, and bulk modulus of ice Ih. It is found that although it overestimates the equilibrium volume, vdW-DF predicts accurate binding energy of ice Ih, as compared with high level quantum chemistry calculations and experiment. Inclusion of the nonlocal correlation, i.e., van der Waals interaction, leads to an overall improvement over the standard generalized gradient approximation in describing water ice.     

Direct dynamics study of ultrafast vibrational energy relaxation in ice Ih

Car-Parrinello molecular dynamics (CPMD) and a previously developed wave packet model are used to study ultrafast relaxation in water clusters. Water clusters of 15 water molecules are used to represent ice Ih. The relaxation is studied by exciting a symmetric or an asymmetric stretch mode of the central water molecule. The CPMD results suggest that relaxation occurs within 100 fs. This is in agreement with experimental work by Woutersen and Bakker and the earlier wave packet calculations. The CPMD results further indicate that the excitation energy is transferred both intramolecularly and intermolecularly on roughly the same time scale. The intramolecular energy transfer occurs predominantly between the symmetric and asymmetric modes while the bend mode is largely left unexcited on the short time scale studied here.     

Investigation of the hydrogen bonding in ice Ih by first-principles density function methods

It is a well recognized difficult task to simulate the vibrational dynamics of ices using the density functional theory (DFT), and there has thus been rather limited success in modelling the inelastic neutron scattering (INS) spectra for even the simplest structure of ice, ice Ih, particularly in the translational region below 400 cm(-1). The reason is partly due to the complex nature of hydrogen bonding (H-bond) among water-water molecules which require considerable improvement of the quantum mechanical simulation methods, and partly owing to the randomness of protons in ice structures which often requires simulation of large super-lattices. In this report, we present the first series of successful simulation results for ice Ih using DFT methods. On the basis of the recent advancement in the DFT programs, we have achieved for the first time theoretical outcomes that not only reproduce the rotational frequencies between 500 to 1200 cm(-1) for ice Ih, but also the two optic peaks at ～240 and 320 cm(-1) in the translational region of the INS spectra [J. C. Li, J. Chem. Phys 105, 6733 (1996)]. Besides, we have also investigated the impact of pairwise configurations of H(2)O molecules on the H-bond and found that different proton arrangements of pairwise H(2)O in the ice Ih crystal lattice could not alter the nature of H-bond as significantly as suggested in an early paper [J. C. Li and D. K. Ross, Nature (London) 365, 327 (1993)], i.e., reproducing the two experimental optic peaks do not need to invoke the two H-bonds as proposed in the previous model which led to considerable debates. The results of this work suggest that the observed optic peaks may be attributed to the coupling between the two bands of H-O stretching modes in H(2)O. The current computational work is expected to shed new light on the nature of the H-bonds in water, and in addition to offer a new approach towards probing the interaction between water and biomaterials for which H-bond is essential.