Isothermal isobaric molecular dynamics simulation of water
Introduction
Because of its importance in chemistry and biology, water and aqueous mixtures remain the subject of intense experimental and theoretical studies over the past few years [6], [7], [8], [9], [10], [11], [12], [13]. The cohesive energy of liquid water is important since water molecules can be trapped with their nearest neighbors through a hydrogen bond network. The general admitted behavior of the molecular liquid is a pseudo short range or local order and a long range disorder. Such scheme projects our vision of liquid structure toward an intermediate profile between the regular ice structure and the chaotic gas phase structure. From the molecular point of view, melting is associated with a partial loss of the crystal lacunars structure [5], water molecules being maintained by hydrogen bonds. An increase of the temperature will amplify the molecular rotational and translational motions.
Since water is a highly polarizable molecule, several alternatives to the commonly used fixed-point charge models [14], [15], [16], [17], [18] were proposed, such as the dipole-polarizable models [19] or the “fluctuating charge” models [20], [21], [22]. Such highly polarizable molecules suggest a systematic analysis of the kinetic aspect (dynamics) with respect to the charge distribution (dipole) and vice versa [23], [14], [24].
The purpose of the following investigations is to put into evidence the effect of the respective permanent dipole moments variations of the trapped molecules (guest or inner molecule type A within the tetrahedral cluster) and the molecular traps (hosts or outer molecules type B) upon the hydrogen bond network construction which characterizes liquid water. Experiments such as Raman spectrometry give evidence of the existence in liquid water of two hydrogen bond types; those remaining unbroken or SHB (Strong Hydrogen Bond) and those broken WHB (Weak Hydrogen Bond) subject to the association and dissociation process, where each broken connection is rebuilt after a short time period of about 200 fs. With the help of the X-ray diffraction technique [1], [2], [3], a team of researchers focused on the liquid water profile at the femtosecond (fs) timescales. The results obtained at room temperature relate to the proportions of molecules having a four-coordination or pseudo crystal (4 SHB) arrangement making 20% of the entire system and those having two connections (2 WHB) for the remainder.
Molecular dynamics simulations performed by Dang and Chang on liquid water [25] have shown that water molecules at the liquid–gas interface have dipole moment values close to that of the gas phase molecule (isolated). They found that the coordination hydrogen number (n) per molecule is (n = 2) in the area close to this interface, whereas (n = 4) in the bulk liquid. In a recent quantum calculation study [26] we observe that the mean charge on oxygen atoms per molecule increases with the cluster size, a result found to be consistent with the Morokuma’s [27], [28], [29] analysis of the total interaction energy. Indeed, it is observed that the main components responsible in the stability of large clusters are the electrostatic energy (ΔEelectrostatic) and charge transfer energy (ΔEchargetransfer) and these components decrease when the cluster size increases. Moreover, an analysis [30] of the hydrogen bond network and strength with the topology and arrangement of water molecules type of large clusters appear to follow a D–A (i.e. 1 Donor molecule “D”–1 Acceptor molecule “A”) and DD–AA (i.e. 2 Donor molecules “DD”–2 Acceptor molecules “AA”) combinations for three-dimensional cluster. Such scheme is in perfect concordance with the tetrahedral water edifice with a central molecule in a DD–AA (4 SHB) state and surrounding molecules in a D–A (2 WHB) state. Fractures and recombination of these hydrogen bond connections simultaneously led us to assume the concept of temporal configuration whose lifespan is of the femtosecond order of magnitude. This perpetual molecular arrangement factor will vary according to the studied system temperature and can be correlated to the diffusion coefficient value D. Experiments give evidence that its value is specific to each temperature and pressure applied [31]. A reduction in the value of D involves an increase of the viscosity η and vice versa.
At the melting point, experimental investigations [5], [32], [33], [34], [35], [36], [37], [38] have shown that the molecular motions in the liquid state are characteristic of such environment and the structural data obtained from X-ray diffraction techniques, Raman spectroscopy and neutron scattering give evidence of the tetrahedral form of liquid water structure. These techniques are able to give a rather realistic image of the distribution in space of the hydrogen and oxygen atoms through the measured structure factor S(q) and the corresponding radial distribution function g(r). Spectroscopic measurements by neutron scattering of ice [39] confirmed the hypothesis of Bernal and Fowler [40] who proposed that molecules are maintained in the condensed phase (ice and bulk liquid water) with their small distortion geometry. The liquid water dielectric constant is very high compared to the associated liquids whose permanent dipole moments are comparable [41], [42], [11], [43], [15], [44], [45], [46], [47]. This fact is a consequence of the correlation between the water molecule dipole moment and the hydrogen bond network characteristic of water. Collective polarization is an additional factor which contributes to increase the dielectric constant value. The experiments showed that the dipole moment value [11], [13], [14] passes from 1.85 D for an isolated water molecule to 2.6 D and up to 3.0 D in the condensed phase.
Attempts made to perform numerical simulations by molecular dynamics or Monte Carlo calculations showed the efficiency and limits of several potential models [48], [49], [50]. Molecular interactions from gas phase models have shown their abilities to give good structural and thermodynamic results for free dimmers. Actually, the practice showed a great difficulty in optimizing a series of parameters relevant of the particular behavior of water molecules, such as the charge distribution, the permanent dipole, the geometry… This optimization allows, however, a better reproduction of the structural and thermodynamic properties of water. Numerous physical models and potential interaction functions are quoted in the scientific literature for Monte Carlo and molecular dynamics simulations, among the most frequently used, one can mention the ST2 model [48] or the simple point charge developed by Stillinger and Rahman for the study of liquid water, the single point charge model SPC and SPCE [49] often used for mixtures and aqueous solution, and the Jorgensen’s model known as the transferable interaction potential to four polarizable sites TIP4P [50].
In our previous work [4], the calculated density plotted as a function of temperature compares favorably well with experiments particularly at low temperatures. The plotted temperature-density profile displays a similar behavior to experimental data in the temperature range of 273–280 K; beyond this range the differences are more significant. The thermal expansion of the box might be responsible of these significant changes and should affect all points since the simulation is undertaken within the same conditions (pressure control, piston mass [35], integration step, temperature control…). As a consequence, our previous model can not be sufficient to take into account every parameter for a better correlation between the molecular behavior and the calculated density. The later will be taken as a validity test for any intermolecular potential model, for more realistic isothermal isobaric molecular dynamics simulation. Therefore, for every given distance of the charge delocalization on the oxygen atom, an isothermal isobaric molecular dynamics simulation is carried out. Each delocalized charge is a point with a given value of the dipole moment μA for the trapped molecule A whose influence in the liquid is reflected on the self-diffusion coefficient. The melting process is followed by a steady variation of D.
Section snippets
Computational aspects and the numerical model
Assuming that water structure reflects a continuous fracture and recombination process of water tetrahedral clusters, the central molecule is in a DD–AA continuum and differences between water molecules are made as either an inner molecule (type A) or an outer molecule (type B) are considered. Oxygen–hydrogen equilibrium distances RO–H are found to be different whether the molecule is occupying a central position or is in the first coordination shell. The complete electronic and structural data
Results and discussion
The permanent dipole moment difference between the two water molecular species A and B reported in Fig. 1 is mainly the consequence of the charge asymmetry on oxygen atoms for A and B type water molecules. The negative charge originally located on the oxygen atom and subsequently delocalized on a dummy site M at a distance dOM gives evidence of the existence of two distinct zones: zone A occupied by a central molecule A (trapped) and zone B formed by the outer molecules B (traps). A range of dOM
Conclusion
In agreement with the experimental results obtained from Raman spectroscopy, X-ray, and neutron scattering, our molecular dynamics simulation shows the necessary existence of two hydrogen bond species (Strong and Weak OH bonds). Their calculated proportions (XSOH and XWOH) give a detailed description of the molecular environment characteristic for each charge delocalization supported by the oxygen atom, when running LPE and HPE molecular dynamics simulations. The pressure increase (2100 bar)
References (57)
- et al.
Chem. Phys.
(1989) - et al.
Food Chemistry
(2008) - et al.
Fluid Phase Equilibr.
(2001) - et al.
J. Chem. Phys.
(1986) - et al.
Science
(2004) - et al.
J. Phys. Chem. B
(2006) - et al.
Science
(2004) - et al.
J. Soc. Algiers Chim.
(2003) - et al.
The Structure and the Properties of Water
(1968) - et al.
J. Chem. Phys.
(1971)
J. Phys. Chem. Ithaca.
Phys. Rev. A
Nature
J. Chem. Phys.
J. Chem. Phys.
Mol. Phys.
J. Chem. Phys.
Mol. Phys.
J. Phys. Chem. B
J. Chem. Phys.
J. Chem. Phys.
J. Chem. Phys.
Mol. Phys.
Chem. Phys. Lett.
J. Chem. Phys.
Cited by (5)
Polarization response of methane encapsulated in water cages
2013, Computational and Theoretical ChemistryCitation Excerpt :The structural models of sI and sII are therefore used in this work. Some classical molecular dynamics simulations have been performed to study the structures and properties [8,19–28] of methane hydrates. A number of ab initio studies [29–36] have been carried out to investigate the electronic structures, stability, guest diffusions, etc, in recent years.
Structure and hydrogen bonding in aqueous sodium chloride solutions using theoretical water model AB<inf>4</inf>: Effects of concentration
2011, Computational and Theoretical ChemistryCitation Excerpt :Several works on dynamic properties of the aqueous solutions were carried out using molecular dynamics simulation [17–20]. In order to supplement recent work on NaCl aqueous solution [21–23], our study will concern the structure of this solution at different concentrations by using the AB4 model [7]. Series of molecular dynamics simulations (DM) were carried out under ambient conditions and (NPT) ensemble following a concentration range from 1.32 to 6.61 mol kg−1.
Structure and dynamics of the peptide strand KRFK from the thrombospondin TSP-1 in water
2018, Journal of Molecular ModelingStructure and stability of hydrogen bonds under conditions of heating in nanopores
2015, High TemperatureSilver nanoparticles behave as hydrophobic solutes towards the liquid water structure in the interaction shell. A Raman study in the O-H stretching region
2009, Physical Chemistry Chemical Physics