Atomistic simulations of liquid water using Lekner electrostatics

Equilibrium molecular dynamics simulations have been performed for liquid water using three different potential models in the NVT and NPT ensembles. The flexible SPC model, the rigid TIP4P model and the rigid/polarizable TIP4P-FQ potential were studied. The Lekner method was used to handle long range electrostatic interactions, and an efficient trivariate cubic spline interpolation method was devised for this purpose. A partitioning of the electrostatic interactions into medium and long range parts was performed, and the concomitant use of multiple timestep techniques led to substantially enhanced computation speeds. The simulations were carried out using 256 molecules in the NVT ensemble at 25°C and 997 kg m^?3 and in the NPT ensemble at 25°C and 1 bar. Various dynamic, structural, dielectric, rotational and thermodynamic properties were calculated, and it was found that the simulation methodologies performed satisfactorily vis-à-vis previous simulation results and experimental observations.     

Effect of electrostatics techniques on the estimation of thermal conductivity via equilibrium molecular dynamics simulation: application to methane hydrate

Equilibrium molecular dynamics (MD) simulations for three system sizes of fully occupied methane hydrate have been performed at around 265 K to estimate the thermal conductivity using the Ewald, Lekner, reaction field, shifted-force and undamped Fennell–Gezelter methods. The TIP4P water model was used in conjunction with a fully atomistic methane potential with which it had been parameterized from quantum simulation. The thermal conductivity was evaluated by integration of the heat flux autocorrelation function (ACF) derived from the Green–Kubo formalism; this approach vas validated by estimation of the average phonon mean free path. The thermal conductivities predicted by non-periodic techniques were in reasonable agreement with the experimental results of 0.62 and 0.68 W/m K, although it was found that the estimates by the non-periodic techniques were up to 25% larger than those of Lekner and Ewald estimates, particularly for larger systems. The results for the Lekner method exhibited the least variation with respect to system size. A decomposition of the heat flux vector into its respective contributions revealed the importance of electrostatic interactions, and how different electrostatic treatments affect the contribution to the thermal conductivity.     

Lekner summations and Ewald summations for quasi-two-dimensional systems

The purpose of this work is to carry out a comparison between the Ewald quasi-2D and Lekner summation methods. These methods were derived to treat the long-range electrostatic interactions in systems periodic in two directions, but bound in the third. The comparison is performed by Monte Carlo simulations on a very simple system, a bilayer of point ions; samplings of the phase space, average energies and structure functions are compared. When correctly implemented, the Lekner summation method is found to be in close agreement with the Ewald quasi-2D method; otherwise, a very complicated bias may plague computations.     

Transport properties of bulk water at 243–550 K: a Comparative molecular dynamics simulation study using SPC/E,TIP4P,and TIP4P/2005 water models

ABSTRACT

Molecular dynamics simulations of various water models – SPC/E (extended simple point charge), TIP4P (transferable intermolecular potential 4 points), and TIP4P/2005 – have been carried out in the canonical (NVT fixed) ensemble over the range of temperatures 243–550?K with Ewald summation. The transport properties (self-diffusion coefficients D, viscosities η, and thermal conductivities λ) of SPC/E, TIP4P, and TIP4P/2005 water were evaluated at 243–550?K and compared with experimental data. The temperature dependence of transport properties of SPC/E, TIP4P and TIP4P/2005 water was discussed to determine how reliable the models are over this temperature range.     

Yukawa potentials in systems with partial periodic boundary conditions. II. Lekner sums for quasi-two-dimensional systems

Yukawa potentials may be long-ranged when the Debye screening length is large. In computer simulations, such long-ranged potentials have to be taken into account with convenient algorithms to avoid systematic bias in the sampling of the phase space. Recently, we provided Ewald sums for quasi-two-dimensional systems with Yukawa interaction potentials [J. Chem. Phys. 126, 056101 (2007); Molec. Phys. paper I of this series]. Sometimes, Lekner sums are used as an alternative to Ewald sums for Coulomb systems. In the present work, we derive the Lekner sums for quasi-two-dimensional systems with Yukawa interaction potentials and we give numerical tests for practical implementations. The main result of this paper is to emphasize that Lekner sums cannot be considered as an alternative to Ewald sums for Yukawa potentials. As a conclusion to this work: Lekner sums should not be used for quasi-two-dimensional systems with Yukawa interaction potentials.     

Molecular dynamics computation of the work of ion solvation: comparison of two models of water

Clusters of the TIP4P water model with 125 and 256 molecules and with the K⁺ or K^? ions were simulated by the molecular dynamics method. The radial profiles of local density, energy, normal pressure and molecular orientation were obtained. The work of cluster formation was calculated from the normal pressure profile. The results for the TIP4P and SPC/E model of water are discussed in comparison with each other. A divergence in behaviour of water molecules in the two models was discovered.     

Fragile to strong crossover and Widom line in supercooled water: A comparative study

The aim of this paper is to discuss the relationship between the dynamics and thermodynamics of water in the supercooled region. Reviewed case studies comprehend bulk water simulated with the SPC/E, TIP4P and TIP4P/2005 potentials, water at protein interfaces, and water in solution with electrolytes. Upon supercooling, the fragile to strong crossover in the α-relaxation of water is found to occur when the Widom line emanating from the liquid-liquid critical point is crossed. This appears to be a general characteristic of supercooled water, not depending on the applied interaction potential and/or different local environments.     

硝酸钾电解质溶液水化结构的分子动力学模拟

应用分子动力学方法对硝酸钾溶液中离子团簇的结构和离子的水化性质进行模拟研究．水分子采用简单点电荷模型，钾离子被看作带电硬球，硝酸根离子采用刚性四节点模型，同时考虑了节点间的库仑长程作用和L-J相互作用，库仑长程作用采用EWALD求和方法处理，得到了溶液中各种离子对的微观构型和径向分布函数，考察了溶液浓度对离子水化性质的影响．研究表明, 在KNO3溶液中存在一定程度的离子缔合，在高离子浓度情形，可以观察到二聚体、三聚体、溶剂分离阴阳离子对以及其它更复杂的离子团簇构型；离子水化数随离子浓度的升高而降低，对不同浓度的溶液得到的K+的水化数为5?7，NO3-的水化数为3.5?4.7，与蒙特卡罗模拟结果和飞行时间中子衍射实验的测量结果一致．     

Effects of an external electromagnetic field on rutile Tio2: A molecular dynamics study

Non-equilibrium molecular dynamics simulations of the bulk rutile phase of TiO₂ have been carried out in an intense external electromagnetic field with frequency in the microwave to far-infrared range. Simulations were performed in constant-volume ensembles with and without a thermostat coupled to the atomic degrees of freedom, at 298 K and from 298 K to well above the crystal melting temperature, respectively. Fields were applied along the a and c crystallographic directions. Both the Lekner and Ewald techniques were used to handle long-range electrostatics and the impact of different choices of Ewald parameters on the results has been evaluated. It was found that the ions respond rapidly to the field. The peaks and troughs in the Ti-Ti radial distribution functions were sharpened relative to the zero-field case at the instants of maximum electric field intensity, but little evidence of this was found for the Ti-O or O-O distributions. For pure Newtonian dynamics, the 500 GHz field excited the vibrational modes to the greatest extent, raising the system temperature at the fastest rate. For temperatures above 2000-2500 K, the crystal structure was found to melt. In the canonical ensemble, the 50 GHz field led to the greatest enhancement of ionic translational mobility.     

Force-field dependence on the interfacial structure of oil–water interfaces

We investigate the performance of different force-fields for alkanes, united (TraPPE) and all atom (OPLS-AA) models, and water (SPC/E and TIP4P-2005), in the prediction of the interfacial structure of alkane (n-octane, and n-dodecane)–water interfaces. We report an extensive comparison of the interfacial thermodynamic properties as well as the interfacial structure (translational and orientational). We use the recently introduced intrinsic sampling method, which removes the averaging effect of the interfacial capillary waves and provides a clear view of the interface structure. The alkane interfacial structure is sensitive to the environment, i.e. alkane–vapour or alkane–water interfaces, showing a stronger structure when it is in contact with the water phase. We find that this structure is fairly independent of the level of detail, full or united atom, employed to describe the alkane phase. The water surface properties show a small dependence on the water model. The dipole moment of the SPC/E model shows asymmetric fluctuations, with a tendency to point both towards the alkane and water phases. On the other hand the dipole moment of the TIP4P-2005 model shows a tendency to point towards the water phase only. Analysis of the intrinsic electrostatic field indicates that the surface water potential is confined to an interfacial region of about 8 Å. Overall we find that the intrinsic structure of alkane–water interfaces is a robust interfacial property, which is independent of the details of the force-field employed. Hence, it should provide a good reference to interpret experimental data.     

Asymptotic analysis of nonlinear nonaxisymmetric waves on the surface of a neutral dielectric jet in a longitudinal electrostatic field

The problem of calculation of finite-amplitude waves on the cylindrical surface of an ideal incompressible dielectric liquid jet in a uniform electrostatic field collinear to the unperturbed jet axis is solved using a second-order asymptotic analytic procedure in ratio of the wave amplitude to the jet radius. Nonlinear corrections to the jet profile, velocity field potential, and electrostatic potentials inside and outside the jet are of resonant nature. The degenerate resonant interaction between the wave determining the initial strain and the waves excited due to nonlinearity of the hydrodynamic equations can take place for waves with different symmetries (different azimuth numbers).     

Comparative atomistic and coarse-grained study of water: What do we lose by coarse-graining?

We employ the inverse Boltzmann method to coarse-grain three commonly used three-site water models (TIP3P, SPC and SPC/E) where one molecule is replaced with one coarse-grained particle with isotropic two-body interactions only. The shape of the coarse-grained potentials is dominated by the ratio of two lengths, which can be rationalized by the geometric constraints of the water clusters. It is shown that for simple two-body potentials either the radial distribution function or the geometrical packing can be optimized. In a similar way, as needed for multiscale methods, either the pressure or the compressibility can be fitted to the all atom liquid. In total, a speed-up by a factor of about 50 in computational time can be reached by this coarse-graining procedure.     

Effect of short- and long-range forces on the structure of water: temperature and density dependence

The individual effects of short-range and long-range forces on the structure of water, a prerequisite for developing a perturbation theory, are assessed using a decomposition of realistic water—water potential models into trial potentials. Computer simulations for one typical liquid density and a number of temperatures ranging from the freezing temperature up to supercritical ones, and for several densities on a supercritical isotherm were performed. The trial potentials were constructed from the ST2 and TIP4P potentials and it is shown that for both potentials the results are practically identical. It is shown that (i) regardless of the thermo-dynamic conditions and potential models used, the structure of water and the mutual orientational arrangement of water molecules, given by a set of site—site correlation functions, are determined nearly exclusively only by the short range forces, and (ii) for high density states the effect of the short range electrostatic part of the intermolecular potential on the spatial arrangement of the water molecules rapidly decreases with increasing temperature but does not become negligible.     

Static dielectric properties of Stockmayer fluids

Dielectric constants for the Stockmayer fluid are computed both from the equilibrium fluctuations of the polarization and from the polarization response to an applied field in a molecular dynamics simulation with periodic boundary conditions and with the dipole interactions treated by the Ewald method. The dielectric constant at finite wavelength is obtained from a Kirkwood-like expression, while for uniform polarization fluctuations a different expression applies due to the absence of surface effects in the Ewald calculation. The results are nearly independent of the number (up to 500) of particles investigated, even at large values of the dipole moment. The dielectric constants obtained from the approximate solution of the hypernetted chain integral equation are considerably larger than the molecular dynamics results except at low dielectric constants. The external field simulation permits study of the non-linear dependence of the dielectric constant on the magnitude of the applied electric field.     

Water and other tetrahedral liquids: order, anomalies and solvation

In order to understand the common features of tetrahedral liquids with water-like anomalies, the relationship between local order and anomalies has been studied using molecular dynamics simulations for three categories of such liquids: (a)?atomistic rigid-body models for water (TIP4P, TIP4P/2005, mTIP3P, SPC/E), (b)?ionic melts, BeF(2) (TRIM model) and SiO(2) (BKS potential) and (c)?Stillinger-Weber liquids parametrized to model water (mW) and silicon. Rigid-body, atomistic models for water and the Stillinger-Weber liquids show a strong correlation between tetrahedral and pair correlation order and the temperature for the onset of the density anomaly is close to the melting temperature. In contrast, the ionic melts show weaker and more variable degrees of correlation between tetrahedral and pair correlation metrics, and the onset temperature for the density anomaly is more than twice the melting temperature. In the case of water, the relationship between water-like anomalies and solvation is studied by examining the hydration of spherical solutes (Na(+), Cl(-), Ar) in water models with different temperature regimes of anomalies (SPC/E, TIP4P and mTIP3P). For both ionic and nonpolar solutes, the local structure and energy of water molecules is essentially the same as in bulk water beyond the second-neighbour shell. The local order and binding energy of water molecules are not perturbed by the presence of a hydrophobic solute. In the case of ionic solutes, the perturbation is largely localized within the first hydration shell. The binding energies for the ions are strongly dependent on the water models and clearly indicate that the geometry of the partial charge distributions, and the associated multipole moments, play an important role. However the anomalous behaviour of the water network has been found to be unimportant for polar solvation.     

Calculation of the vapour-liquid coexistence curve for a fluctuating point charge water model

Recent studies of the phase behaviour of several polarizable water models have shown that typically the inclusion of polarizability into current non-polarizable, rigid models of water actually decreases the accuracy of their predictions of the liquid-vapour coexistence curve. Most of the polarizable models examined so far have been shown to characteristically over-predict the saturation pressure and correspondingly under-predict the critical temperature. It is the purpose of this Research Note to report the results of a Gibbs—Duhem integration of the liquid—vapour coexistence curve for the TIP4P-FQ water model and to discuss the possible reasons behind the lack of improvement of polarizable water models over their non-polarizable counterparts.     

Hydration of apolar solutes of varying size: a systematic study

We report results of a systematic and extensive study on the orientational structure of water molecules in the vicinity of apolar solutes of varying size in infinitely diluted aqueous solutions at ambient conditions. A wide range of sizes of the solute, modelled as a Lennard–Jones particle, has been considered including also the limiting case of a smooth planar wall (i.e. a solute of an infinite diameter). Both the bivariate and monovariate angle distributions are used to describe the orientation of water molecules. It is demonstrated that the commonly used latter method may lead to misleading or even erroneous conclusions. It is found that all three considered models of water, SPC/E, TIP4P and TIP5P, exhibit, unlike some simplified models, qualitatively the same behaviour, i.e. with increasing size of the solute only the water molecules in closest proximity of the solute adopt orientations which sacrifice a possible hydrogen bond, whereas the preferred orientation of water molecules farther from the surface, but yet belonging to the first solvation shell, is independent of the solute size.     

Thermodynamic and dielectric properties of polar lattices

Monte Carlo results are presented for the free energy and dielectric constant of systems of permanent dipoles disposed on the sites of simple cubic and face-centred cubic lattices as functions of the parameter ? = μ²ρ/kT. Alternative schemes are considered for taking account of the long-range character of the dipolar interaction, based on an Ewald-type sum or a reaction-field approximation; use of such a procedure is essential if the system is to have the correct dielectric properties. Different methods of calculating the dielectric constant are also compared. It is shown that the two methods of treating the long-range contribution to the energy lead to similar results for the dielectric constant, but use of the Ewald method apparently introduces a small but systematic error. The results on both free energy and dielectric constant are discussed in the light of predictions of a number of analytical approaches. In the case of the simple cubic lattice the general trend in the computed dielectric constant as a function of ? is in closer accord with the solution of the mean spherical approximation than with a simple Padé approximant to the Clausius-Mossotti function.     

Ewald method for polytropic potentials in arbitrary dimensionality

The Ewald summation technique is generalized to power-law 1/| r | ^k potentials in three-, two- and one-dimensional geometries with explicit formulae for all the components of the sums. The cases of short-range, long-range and ‘marginal’ interactions are treated separately. The jellium model, as a particular case of a charge-neutral system, is discussed and the explicit forms of the Ewald sums for such a system are presented. A generalized form of the Ewald sums for a non-cubic (non-square) simulation cell for three- (two-) dimensional geometry is obtained and its possible field of application is discussed. A procedure for the optimization of the involved parameters in actual simulations is developed and an example of its application is presented.