Quantum effects in liquid water and ice: model dependence

This paper explores the influence of choice of potential model on the quantum effects observed in liquid water and ice. This study utilizes standard rigid models and a more formal context for the rigid-body centroid molecular dynamics methodology used to perform the quantum simulations is provided. Quantum and classical molecular dynamics simulations are carried out for liquid water and ice Ih at 298 and 220 K, respectively, with the simple point charge/extended and TIP4P-Ew water models. The results obtained for equilibrium and dynamical properties are compared with those recently reported on TIP4P [L. Hernandez de la Pena and P. G. Kusalik, J. Chem. Phys. 121, 5992 (2004); L. Hernandez de la Pena et al., J. Chem. Phys 123, 144506 (2005)]. For the liquid, an energy shift of about 8% and an average molecular uncertainty of about 11 degrees were found independently of the water model. The self-diffusion coefficient consistently increases by more than 50% when going from the classical to the quantum system and quantum dynamics are found to reproduce the experimental isotopic shifts with the models examined. The ice results compare remarkably well with those previously reported for the TIP4P water model; they confirm that quantum effects are considerable and that the quantum mechanical uncertainty and the energy shifts due to quantization are smaller in ice than in liquid water. The relevance of these findings in the context of the construction of water models is briefly discussed.     

石英玻璃分子动力学模拟中的原子电荷转移与系综选择

介绍了SiO2体系分子动力学模拟中的Si、O原子电荷转移问题;采用Morse势函数研究了原子电荷转移对石英玻璃模拟的影响,发现原子电荷转移在影响模型密度的同时,还直接影响着原子的最近邻距离.NPT和NVT系综下的模拟结果对比显示,系综对模型中原子最近邻情况影响不大,但在NVT系综下模拟结果表明实际玻璃中存在的较大的空隙结构,找到了以往模拟中密度结果偏高的原因,提出了一种较好的石英玻璃分子动力学建模的方法.该方法不但解决了在调整电荷时维持原子最近邻距离与保证模型密度之间的矛盾,而且可以很好地描述石英玻璃在远程结构上密度不均、存在较大空隙的无序结构.此外,原子自扩散系数的计算结果展示了空隙结构在石英玻璃扩散性质研究中的作用.     

Structure and dynamics of liquid water with different long-range interaction truncation and temperature control methods in molecular dynamics simulations

We have used molecular dynamics simulations to study the physical properties of modified TIP3P water model included in the CHARMM program, using four different methods-the Ewald summation technique, and three different spherical truncation methods-for the treatment of the long-range interactions. Both the structure and dynamics of the liquid water model were affected by the methods used to truncate the long-range interactions. For some of the methods artificial structuring of the model liquid was observed around the cutoff radius. The model liquid properties were also affected by the commonly applied temperature control methods. Four different methods for controlling the temperature of the system were studied, and the effects of these methods on the bulk properties for liquid water were analyzed. The system size was also found to change the dynamics of the model liquid water. Two control simulations with the SPC/E water model were carried out. The self-diffusion coefficient (D), the radial distribution function (g(OO)), the distance dependent Kirkwood G-factor [G(k)(r)] and the intermolecular potential energy (E(pot)) were determined from the different trajectories and compared with the experimental data.     

The ice-vapor interface and the melting point of ice I(h) for the polarizable POL3 water model

We use molecular dynamics simulations to determine the melting point of ice I(h) for the polarizable POL3 water force field (Dang, L. X. J. Chem. Phys.1992, 97, 2659). Simulations are performed on a slab of ice I(h) with two free surfaces at several different temperatures. The analysis of the time evolution of the total energy in the course of the simulations at the set of temperatures yields the melting point of the POL3 model to be T(m) = 180 ± 10 K. Moreover, the results of the simulations show that the degree of hydrogen-bond disorder occurring in the bulk of POL3 ice is larger (at the corresponding degree of undercooling) than in ice modeled by nonpolarizable water models. These results demonstrate that the POL3 water force field is rather a poor model for studying ice and ice-liquid or ice-vapor interfaces. While a number of polarizable water models have been developed over the past years, little is known about their performance in simulations of supercooled water and ice. This study thus highlights the need for testing of the existing polarizable water models over a broad range of temperatures, pressures, and phases, and developing a new polarizable water force field, reliable over larger areas of the phase diagram.     

Development of an improved four-site water model for biomolecular simulations: TIP4P-Ew

A re-parameterization of the standard TIP4P water model for use with Ewald techniques is introduced, providing an overall global improvement in water properties relative to several popular nonpolarizable and polarizable water potentials. Using high precision simulations, and careful application of standard analytical corrections, we show that the new TIP4P-Ew potential has a density maximum at approximately 1 degrees C, and reproduces experimental bulk-densities and the enthalpy of vaporization, DeltaH(vap), from -37.5 to 127 degrees C at 1 atm with an absolute average error of less than 1%. Structural properties are in very good agreement with x-ray scattering intensities at temperatures between 0 and 77 degrees C and dynamical properties such as self-diffusion coefficient are in excellent agreement with experiment. The parameterization approach used can be easily generalized to rehabilitate any water force field using available experimental data over a range of thermodynamic points.     

An optimized molecular potential for carbon dioxide

An optimized molecular potential model for carbon dioxide is presented in this paper. Utilizing the established techniques of molecular-dynamics and histogram reweighting grand canonical Monte Carlo simulations, this model is demonstrated to show excellent predictability for thermodynamic, transport, and liquid structural properties in a wide temperature-pressure range with remarkable accuracies. The average deviations of this new model from experimental data for the saturated liquid densities, vapor densities, vapor pressures, and heats of vaporization are around 0.1%, 2.3%, 0.7%, and 1.9%, respectively. The calculated critical point is almost pinpointed by the new model. The experimental radial distribution functions ranging from 240.0 to 473.0 K are well reproduced as compared to neutron-diffraction measurements. The predicted self-diffusion coefficients are in good agreement with the nuclear-magnetic-resonance measurements. The previously published potential models for CO2 are also systematically evaluated, and our proposed new model is found to be superior to the previous models in general.     

A general purpose model for the condensed phases of water: TIP4P/2005

A potential model intended to be a general purpose model for the condensed phases of water is presented. TIP4P/2005 is a rigid four site model which consists of three fixed point charges and one Lennard-Jones center. The parametrization has been based on a fit of the temperature of maximum density (indirectly estimated from the melting point of hexagonal ice), the stability of several ice polymorphs and other commonly used target quantities. The calculated properties include a variety of thermodynamic properties of the liquid and solid phases, the phase diagram involving condensed phases, properties at melting and vaporization, dielectric constant, pair distribution function, and self-diffusion coefficient. These properties cover a temperature range from 123 to 573 K and pressures up to 40,000 bar. The model gives an impressive performance for this variety of properties and thermodynamic conditions. For example, it gives excellent predictions for the densities at 1 bar with a maximum density at 278 K and an averaged difference with experiment of 7 x 10(-4) g/cm3.     

Assessing group-based cutoffs and the Ewald method for electrostatic interactions in clusters and in saturated,superheated, and supersaturated vapor phases of dipolar molecules

Monte Carlo simulations in the canonical, isobaric-isothermal, grand canonical, and Gibbs ensembles were used to assess whether the computationally expensive Ewald summation method for the computation of the first-order electrostatic energy can be replaced with a simpler truncation approach for accurate simulations of the saturated, superheated, and supersaturated vapor phases of dipolar and hydrogen-bonding molecules. Rotationally averaged hydrogen fluoride dimer and trimer energies, thermophysical properties and aggregation in the superheated vapor phase of hydrogen fluoride, nucleation free energy barriers for water, and the vapor–liquid coexistence properties of hydrogen fluoride and water were investigated over a wide range of state points. We find that for densities not too close to the critical density, results obtained from simulations using a spherical potential truncation based on neutral groups (molecules or fragments) for the Coulomb interactions are statistically identical to those obtained using the Ewald summation method. Use of the simpler spherical truncation results in a significant reduction of the computational effort for simulations employing molecular mechanics force fields and also allows for straightforward implementation of many-body expansion methods to compute the potential energy from electronic structure calculations of subsystems of the entire vapor-phase system.     

Modified free volume theory of self-diffusion and molecular theory of shear viscosity of liquid carbon dioxide

In previous work on the density fluctuation theory of transport coefficients of liquids, it was necessary to use empirical self-diffusion coefficients to calculate the transport coefficients (e.g., shear viscosity of carbon dioxide). In this work, the necessity of empirical input of the self-diffusion coefficients in the calculation of shear viscosity is removed, and the theory is thus made a self-contained molecular theory of transport coefficients of liquids, albeit it contains an empirical parameter in the subcritical regime. The required self-diffusion coefficients of liquid carbon dioxide are calculated by using the modified free volume theory for which the generic van der Waals equation of state and Monte Carlo simulations are combined to accurately compute the mean free volume by means of statistical mechanics. They have been computed as a function of density along four different isotherms and isobars. A Lennard-Jones site-site interaction potential was used to model the molecular carbon dioxide interaction. The density and temperature dependence of the theoretical self-diffusion coefficients are shown to be in excellent agreement with experimental data when the minimum critical free volume is identified with the molecular volume. The self-diffusion coefficients thus computed are then used to compute the density and temperature dependence of the shear viscosity of liquid carbon dioxide by employing the density fluctuation theory formula for shear viscosity as reported in an earlier paper (J. Chem. Phys. 2000, 112, 7118). The theoretical shear viscosity is shown to be robust and yields excellent density and temperature dependence for carbon dioxide. The pair correlation function appearing in the theory has been computed by Monte Carlo simulations.     

A classical polarizable model for simulations of water and ice

We develop a classical rigid polarizable model of water for molecular simulations of water and ice. The model uses the Rowlinson five-site geometry: oxygen bearing the Lennard-Jones interaction and linearly polarizable point dipole, two positively charged hydrogens, and two massless negative charges placed symmetrically off oxygen so that the experimental dipole moment is reproduced. The target properties are the densities of water and ice Ih, diffusivity, enthalpies of fusion and vaporization, and the ice Ih melting point. The surface tension at lower temperatures is by 7% underestimated whereas the dielectric constant by 6% overestimated. Diffusivity and viscosity worsen at higher temperatures, although the Stokes radius is overestimated only by 2-7%. The ice Ih melting temperature is 260 K and the temperature of maximum density is 269 K. Rescaling the charges by a factor of 1.01 and Lennard-Jones energy by 1.0201 improves the melting point and energy-related quantities but shifts the agreement of kinetic properties to higher temperatures. For the model we propose abbreviation POL4D.     

Computer simulation of the 13 crystalline phases of ice

As a reference for follow-up studies toward more accurate model parametrizations, we performed molecular-dynamics and Monte Carlo simulations for all known crystalline phases of ice, as described by the simple point-charge/extended and TIP4P water models. We started from the measured structures, densities, and temperatures, and carried out classical canonical simulations for all these arrangements. All simulated samples were cooled down close to 0 K to facilitate the comparison with theoretical estimates. We determined configurational internal energies as well as pressures, and monitored how accurately the measured configurations were preserved during the simulations. While these two models predicted very similar thermophysical and structural properties for water at ambient conditions, the predicted features for the corresponding ice polymorphs may differ significantly.     

Transport properties of the rough hard sphere fluid

Results are presented of a systematic study of the transport properties of the rough hard sphere fluid. The rough hard sphere fluid is a simple model consisting of spherical particles that exchange linear and angular momenta, and energy upon collision. This allows a study of the sole effect of particle rotation upon fluid properties. Molecular dynamics simulations have been used to conduct extensive benchmark calculations of self-diffusion, shear and bulk viscosity, and thermal conductivity coefficients. As well, the validity of several kinetic theory equations have been examined at various levels of approximation as a function of density and translational-rotational coupling. In particular, expressions from Enskog theory using different numbers of basis sets in the representation of the distribution function were tested. Generally Enskog theory performs well at low density but deviates at larger densities, as expected. The dependence of these expressions upon translational-rotational coupling was also examined. Interestingly, even at low densities, the agreement with simulation results was sometimes not even qualitatively correct. Compared with smooth hard sphere behaviour, the transport coefficients can change significantly due to translational-rotational coupling and this effect becomes stronger the greater the coupling. Overall, the rough hard sphere fluid provides an excellent model for understanding the effects of translational-rotational coupling upon transport coefficients.     

Comparison of rigid and flexible simple point charge water models at supercritical conditions

This study investigates the differences between the predictions of various properties of rigid and flexible simple point charge water models at supercritical conditions. Molecular dynamics simulations were conducted for supercritical water in a temperature range of 773–1073 K and densities in the range 115–659 kg/m³. We present thermodynamic data, pair correlation functions, self-diffusivity, power spectra, dielectric constants, and variaous measures of hydrogen bonding at different state conditions. The flexible water model performs better in predicting the pressures along the supercritical isotherms simulated. Agreement between experimental and calculated dielectric constants is superior for the flexible water model, particularly at high densities. The flexible model exhibits a greater degree of hydrogen bonding and more persistent hydrogen bonds than does the rigid model. The structural features of supercritical water at high densities are identical for the two water models. At low densities, however, the flexible potential exhibits pair correlation functions with enhanced peaks. Inclusion of flexibility in the potential model does not result in a significant shift in the position of the rotational/librational peak in the power spectrum. The self-diffusivities obtained from the simulations are within the accuracy of the experimental values for both the rigid and flexible models. On balance the inclusion of flexibility improves agreement with the properties of real supercritical water while incurring little or no additional computational burden. © 1996 by John Wiley & Sons, Inc.     

Dynamical and structural properties of benzene in supercritical water

We have employed an anisotropic united atom model of benzene (R. O. Contreras, Ph.D. thesis, Universitat Rovira i Virgili 2002) that reproduces the quadrupolar moment of this molecule through the inclusion of seven point charges. We show that this kind of interaction is required to reproduce the solvation of these molecules in supercritical water. We have computed self-diffusion coefficient and Maxwell-Stefan coefficients as well as the shear viscosity for the mixture water-benzene at supercritical conditions. A strong density and composition dependence of these properties is observed. In addition, our simulations are in qualitative agreement with the experimental evidence that, at medium densities (0.6 g/cm(3) and 673 K), almost half of the benzene molecules have one hydrogen bond with water molecules. We also observe that these bonds are longer lived than the corresponding hydrogen bonds between water molecules. Similarly, we obtain an important reduction of the dielectric constant of the mixture with the increment of the amount of benzene molecules at medium and high densities.     

Semiempirical self-consistent polarization description of bulk water, the liquid-vapor interface, and cubic ice

We have applied an efficient electronic structure approach, the semiempirical self-consistent polarization neglect of diatomic differential overlap (SCP-NDDO) method, previously parametrized to reproduce properties of water clusters by Chang, Schenter, and Garrett [ J. Chem. Phys. 2008 , 128 , 164111 ] and now implemented in the CP2K package, to model ambient liquid water at 300 K (both the bulk and the liquid-vapor interface) and cubic ice at 15 and 250 K. The SCP-NDDO potential retains its transferability and good performance across the full range of conditions encountered in the clusters and the bulk phases of water. In particular, we obtain good results for the density, radial distribution functions, enthalpy of vaporization, self-diffusion coefficient, molecular dipole moment distribution, and hydrogen bond populations, in comparison to experimental measurements.     

A transferable classical potential for the water molecule

We developed a new model for the water molecule which contains only three Gaussian charges. Using the gas-phase geometry the dipole moment of the molecule matches, the quadrupole moment closely approximates the experimental values. The negative charge is connected by a harmonic spring to its gas-phase position. The polarized state is identified by the equality of the intermolecular electrostatic force and the spring force acting on the negative charge. In each timestep the instantaneous position of the massless negative charge is determined by iteration. Using the technique of Ewald summation, we derived expressions for the potential energy, the forces, and the pressure for Gaussian charges. The only properties to be fitted are the half-width values of the Gaussian charge distributions and the parameters of the nonelectrostatic repulsion-attraction potential. We determined the properties of gas-phase clusters up to six molecules, the internal energy and density of ambient water and hexagonal ice. We calculated the equilibrium density of ice VII as a function of pressure. As an additional test, we calculated the pair-correlation function, the isotherm compressibility, the heat capacity, and the self-diffusion coefficients for ambient water. As far as we know, this is the first classical model of water which is able to estimate both ends of the phase diagram, the high pressure ice VII, and the gas clusters of water with excellent accuracy.     

A non-polarizable model of water that yields the dielectric constant and the density anomalies of the liquid: TIP4Q

A four-site rigid water model is presented, whose parameters are fitted to reproduce the experimental static dielectric constant at 298 K, the maximum density of liquid water and the equation of state at low pressures. The model has a positive charge on each of the three atomic nuclei and a negative charge located at the bisector of the HOH bending angle. This charge distribution allows increasing the molecular dipole moment relative to four-site models with only three charges and improves the liquid dielectric constant at different temperatures. Several other properties of the liquid and of ice Ih resulting from numerical simulations with the model are in good agreement with experimental values over a wide range of temperatures and pressures. Moreover, the model yields the minimum density of supercooled water at 190 K and the minimum thermal compressibility at 310 K, close to the experimental values. A discussion is presented on the structural changes of liquid water in the supercooled region where the derivative of density with respect to temperature is a maximum.     

An application of flexible constraints in Monte Carlo simulations of the isobaric--isothermal ensemble of liquid water and ice Ih with the polarizable and flexible mobile charge densities in harmonic oscillators model

The method of flexible constraints was implemented in a Monte Carlo code to perform numerical simulations of liquid water and ice Ih in the constant number of molecules, volume, and temperature and constant pressure, instead of volume ensembles, using the polarizable and flexible mobile charge densities in harmonic oscillators (MCDHO) model. The structural and energetic results for the liquid at T=298 K and rho=997 kg m(-3) were in good agreement with those obtained from molecular dynamics. The density obtained at P=1 atm with flexible constraints, rho=1008 kg m(-3), was slightly lower than with the classical sampling of the intramolecular vibrations, rho=1010 kg m(-3). The comparison of the structures and energies found for water hexamers and for ice Ih with six standard empirical models to those obtained with MCDHO, show this latter to perform better in describing water far from ambient conditions: the MCDHO minimum lattice energy, density, and lattice constants were in good agreement with experiment. The average angle HOH of the water molecule in ice was predicted to be slightly larger than in the liquid, yet 1.2% smaller than the experimental value.