Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results
We have investigated temperature trends of the microscopic structure of the SPC/E and TIP4P-Ew water models in terms of the pair distribution functions, coordination numbers, the average number of hydrogen bonds, the distribution of bonding states of a single molecule as well as the angular distribu...
Gespeichert in:
| Veröffentlicht in: | Condensed Matter Physics |
|---|---|
| Datum: | 2015 |
| Hauptverfasser: | , , |
| Format: | Artikel |
| Sprache: | English |
| Veröffentlicht: |
Інститут фізики конденсованих систем НАН України
2015
|
| Online Zugang: | https://nasplib.isofts.kiev.ua/handle/123456789/153513 |
| Tags: |
Tag hinzufügen
Keine Tags, Fügen Sie den ersten Tag hinzu!
|
| Назва журналу: | Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| Zitieren: | Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results / E. Galicia-Andrés, H. Dominguez, O. Pizio // Condensed Matter Physics. — 2015. — Т. 18, № 1. — С. 13603:1-11. — Бібліогр.: 42 назв. — англ. |
Institution
Digital Library of Periodicals of National Academy of Sciences of Ukraine| id |
nasplib_isofts_kiev_ua-123456789-153513 |
|---|---|
| record_format |
dspace |
| spelling |
Galicia-Andrés, E. Dominguez, H. Pizio, O. 2019-06-14T10:42:47Z 2019-06-14T10:42:47Z 2015 Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results / E. Galicia-Andrés, H. Dominguez, O. Pizio // Condensed Matter Physics. — 2015. — Т. 18, № 1. — С. 13603:1-11. — Бібліогр.: 42 назв. — англ. 1607-324X DOI:10.5488/CMP.18.13603 arXiv:1504.01217 PACS: 61.20.-p, 61.20.Gy, 61.20.Ja, 65.20.Jk https://nasplib.isofts.kiev.ua/handle/123456789/153513 We have investigated temperature trends of the microscopic structure of the SPC/E and TIP4P-Ew water models in terms of the pair distribution functions, coordination numbers, the average number of hydrogen bonds, the distribution of bonding states of a single molecule as well as the angular distribution of molecules by using the constant pressure molecular dynamics simulations. The evolution of the structure is put in correspondence with the dependence of water density on high temperatures down to the region of temperatures where the system becomes supercooled. It is shown that the fraction of molecules with three and four bonds determine the maximum density for both models. Moreover, the temperature dependence of the dielectric constant is obtained and analyzed. Ми дослiдили температурнi залежностi мiкроскопiчної структури моделей води SPC/E i TIP4P-Ew в термiнах парних функцiй розподiлу, координацiйних чисел, середнього числа водневих зв’язкiв, розподiлу зв’язаних станiв однiєї молекули, використовуючи метод молекулярної динамiки при постiйногму тиску. Еволюцiя структури поставлена у вiдповiднiсть до залежностi густини води вiд температури в областi, що знаходиться мiж високими температурами i температурами, де система стає переохолодженою. Показано, що частка молекул з трьома i чотирма зв’язками визначає максимальну густину для обох моделей. Бiльш того, отримано i проаналiзовано температурну залежнiсть дiелектричної сталої. E.G. was supported by CONACyT of Mexico under Ph.D. scholarship. E.G. and O.P. are grateful to Dr. T. Patsahan for very helpful discussions and valuable comments. O.P. is grateful to David Vazquez for technical support of this work. en Інститут фізики конденсованих систем НАН України Condensed Matter Physics Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results Температурна залежнiсть мiкроскопiчної структури i аномалiя густини в моделях води SPC/E та TIP4P-Ew. Результати моделювання методом молекулярної динамiки Article published earlier |
| institution |
Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| collection |
DSpace DC |
| title |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results |
| spellingShingle |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results Galicia-Andrés, E. Dominguez, H. Pizio, O. |
| title_short |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results |
| title_full |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results |
| title_fullStr |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results |
| title_full_unstemmed |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results |
| title_sort |
temperature dependence of the microscopic structure and density anomaly of the spc/e and tip4p-ew water models. molecular dynamics simulation results |
| author |
Galicia-Andrés, E. Dominguez, H. Pizio, O. |
| author_facet |
Galicia-Andrés, E. Dominguez, H. Pizio, O. |
| publishDate |
2015 |
| language |
English |
| container_title |
Condensed Matter Physics |
| publisher |
Інститут фізики конденсованих систем НАН України |
| format |
Article |
| title_alt |
Температурна залежнiсть мiкроскопiчної структури i аномалiя густини в моделях води SPC/E та TIP4P-Ew. Результати моделювання методом молекулярної динамiки |
| description |
We have investigated temperature trends of the microscopic structure of the SPC/E and TIP4P-Ew water models in terms of the pair distribution functions, coordination numbers, the average number of hydrogen bonds, the distribution of bonding states of a single molecule as well as the angular distribution of molecules by using the constant pressure molecular dynamics simulations. The evolution of the structure is put in correspondence with the dependence of water density on high temperatures down to the region of temperatures where the system becomes supercooled. It is shown that the fraction of molecules with three and four bonds determine the maximum density for both models. Moreover, the temperature dependence of the dielectric constant is obtained and analyzed.
Ми дослiдили температурнi залежностi мiкроскопiчної структури моделей води SPC/E i TIP4P-Ew в термiнах парних функцiй розподiлу, координацiйних чисел, середнього числа водневих зв’язкiв, розподiлу
зв’язаних станiв однiєї молекули, використовуючи метод молекулярної динамiки при постiйногму тиску.
Еволюцiя структури поставлена у вiдповiднiсть до залежностi густини води вiд температури в областi, що
знаходиться мiж високими температурами i температурами, де система стає переохолодженою. Показано, що частка молекул з трьома i чотирма зв’язками визначає максимальну густину для обох моделей.
Бiльш того, отримано i проаналiзовано температурну залежнiсть дiелектричної сталої.
|
| issn |
1607-324X |
| url |
https://nasplib.isofts.kiev.ua/handle/123456789/153513 |
| citation_txt |
Temperature dependence of the microscopic structure and density anomaly of the SPC/E and TIP4P-Ew water models. Molecular dynamics simulation results / E. Galicia-Andrés, H. Dominguez, O. Pizio // Condensed Matter Physics. — 2015. — Т. 18, № 1. — С. 13603:1-11. — Бібліогр.: 42 назв. — англ. |
| work_keys_str_mv |
AT galiciaandrese temperaturedependenceofthemicroscopicstructureanddensityanomalyofthespceandtip4pewwatermodelsmoleculardynamicssimulationresults AT dominguezh temperaturedependenceofthemicroscopicstructureanddensityanomalyofthespceandtip4pewwatermodelsmoleculardynamicssimulationresults AT pizioo temperaturedependenceofthemicroscopicstructureanddensityanomalyofthespceandtip4pewwatermodelsmoleculardynamicssimulationresults AT galiciaandrese temperaturnazaležnistʹmikroskopičnoístrukturiianomaliâgustinivmodelâhvodispcetatip4pewrezulʹtatimodelûvannâmetodommolekulârnoídinamiki AT dominguezh temperaturnazaležnistʹmikroskopičnoístrukturiianomaliâgustinivmodelâhvodispcetatip4pewrezulʹtatimodelûvannâmetodommolekulârnoídinamiki AT pizioo temperaturnazaležnistʹmikroskopičnoístrukturiianomaliâgustinivmodelâhvodispcetatip4pewrezulʹtatimodelûvannâmetodommolekulârnoídinamiki |
| first_indexed |
2025-11-27T02:28:44Z |
| last_indexed |
2025-11-27T02:28:44Z |
| _version_ |
1850791477811085312 |
| fulltext |
Condensed Matter Physics, 2015, Vol. 18, No 1, 13603: 1–11
DOI: 10.5488/CMP.18.13603
http://www.icmp.lviv.ua/journal
Temperature dependence of the microscopic
structure and density anomaly of the SPC/E and
TIP4P-Ew water models. Molecular dynamics
simulation results
E. Galicia-Andrés1, H. Dominguez2, O. Pizio1∗
1 Instituto de Química, Universidad Nacional Autónoma de México, Circuito Exterior, 04510, México, D.F., México
2 Instituto de Investigaciones en Materiales, Universidad Nacional Autónoma de México, Circuito Exterior,
04510, México, D.F., México
Received November 6, 2014, in final form January 13, 2015
We have investigated temperature trends of the microscopic structure of the SPC/E and TIP4P-Ew water models
in terms of the pair distribution functions, coordination numbers, the average number of hydrogen bonds, the
distribution of bonding states of a single molecule as well as the angular distribution of molecules by using
the constant pressure molecular dynamics simulations. The evolution of the structure is put in correspondence
with the dependence of water density on high temperatures down to the region of temperatures where the
system becomes supercooled. It is shown that the fraction of molecules with three and four bonds determine
the maximum density for both models. Moreover, the temperature dependence of the dielectric constant is
obtained and analyzed.
Key words: water models, density anomaly, molecular dynamics
PACS: 61.20.-p, 61.20.Gy, 61.20.Ja, 65.20.Jk
1. Introduction
It is our honor and pleasure to contribute to this special issue and to dedicate this paper to Prof.
Dr. Douglas Henderson on the occasion of his 80th birthday. Almost four decades ago, Barker and Hen-
derson published their seminal paper “What is “liquid”? Understanding the states of matter”, that com-
prehensively describes basic theoretical methods to explain the equilibrium properties of simple liquids
[1]. These methods have made substantial progress since that time. In particular, computer simulations
techniques have become essentially powerful and common tools in the theory of liquids, mixtures and so-
lutions. Computer simulations permitted to reach enormous progress in understanding the liquids with
hydrogen bonds, more specifically, water and aqueous solutions, for needs of physical chemistry and for
several inter-disciplinary areas involving complex biological molecules and membranes. In particular,
Henderson with co-workers was actively interested in the properties and behavior of inhomogeneous
aqueous solution, see e.g. [2–4].
The studies of water have long history that has been documented in an enormous number of publica-
tions. Of particular interest are the peculiar thermodynamic and dynamic properties of water manifested
in a set of anomalies. In spite of formal equivalence of the relation between themicroscopic structure and
thermodynamic properties discussed in great detail by Barker and Henderson for simple homogeneous
liquids [1], in the case of water it is inevitable to involve the notion of hydrogen bonding. Many works
have dealt with the exploration of hydrogen bonds network structure, its dynamics and life-time of bonds
at an early stage of computer simulations of aqueous systems, see e.g. a few original papers [5, 6] as well
∗E-mail: pizio@unam.mx
© E. Galicia-Andrés, H. Dominguez, O. PizioE-mail: pizio@unam.mx, 2015 13603-1
http://dx.doi.org/10.5488/CMP.18.13603
http://www.icmp.lviv.ua/journal
E. Galicia-Andrés, H. Dominguez, O. Pizio
as later reviews [7–9], that contain an extensive list of references on the subject. Moreover, it is known
that unique properties of water specifically manifest themselves in the region of states below the freezing
point corresponding to the supercooled liquid, see e.g. a recent review [10] and references therein.
It is seducing to relate at least some of water anomalies with the transformations occuring in the
hydrogen bonds cooperative structure. Therefore, in this short report our focus is in the description of
the evolution of the average number of H-bonds per water molecule and in the distribution of the num-
bers of H-bonds per molecule upon temperature changes. Moreover, we would like to relate the trends of
the behavior of these properties with the density anomaly. While the properties in question describe the
formation of “agglomerates” of bonded molecules at relatively short-range distances between them, we
would like to get additional insights into the simultaneously occuring changes of the long-range correla-
tions between water molecules in terms of the dielectric constant.
One of the most studied anomalies of water is its nonmonotonous temperature dependence of density,
ρ(T ), see e.g. a review providing a nice historical insight [9]. Several non-polarizable water models with
different number of interaction sites are capable of reproducing the presence of maximum density of
water with different precision.
Namely, the temperature of maximum density (TMD) at ambient pressure for most frequently used
models is observed at 253 K (TIP4P), at 278 K (TIP4P/2005), at 274.15 K (TIP4P-EW), at 241 K (SPC/E) and at
182 K (TIP3P). These values for temperature are taken from the table 2 of the review by Vega and Abascal
[11]. The experimental results for these and other properties used by the authors are taken from [12].
The experiment yields 277 K for the TMD of real water. On the other hand, the water density at TMD is
0.988 (TIP4P), 0.993 (TIP4P/2005), 0.994 (SPC/E), 0.98 (TIP3P), while the experimental value is 0.997. Some
deviation from these values reported in different works result from technical details of each simulation.
The origin of density anomaly is commonly attributed to the evolution of hydrogen bond network
with the formation of an open structure as the temperature is lowered or, on the other hand, manifests
itself as the reflection of a possible liquid-liquid phase transition observed by simulations in some water
models [9]. The latter and other scenarios for peculiar properties of supercooled water were documented
in [10].
To summarize this brief introduction, our intention in this work is to to present our own data of
molecular dynamics simulations of the temperature dependence of water density, ρ(T ), in an ample in-
terval of temperature (from 370 K down to 170 K) at ambient pressure by using the SPC/Emodel. Next, we
would like to compare the obtained data with the results previously reported by other authors. However,
the principal focus is to follow changes of water structure in terms of different descriptors involving both
the hydrogen bonds and the ρ(T ) dependence in order to establish a clear relationship between them.
Nevertheless, for the sake of comparison of our observations for the SPC/E, and to make the conclusions
more sound, similar results for the TIP4P-Ew water model [13] are presented as well. Our interest in
the SPC/E model is motivated by its wide use for different purposes, in particular for the description of
thermodynamic properties of mixtures. Extending the previous research in mixtures, having in mind a
more ample project presently being developed in this laboratory, our intention is to use the SPC/E model
as a starting point in exploring water–alcohol and water–DMSO solutions by using molecular dynamics
simulations.
2. Simulation details
All our simulations of water models were performed in the isothermal-isobaric ensemble at ambient
pressure. The initial configuration was constructed with 2000 water molecules placed in a cubic array in
a simulation box. Then, the system was simulated in an ample range of temperatures, from a high tem-
perature 370 K down to 180 K. The DL-POLY Classic package was employed [14]. We used the Berendsen
thermostat and barostat, the running timestep was set to 0.002 ps. Commonly, the periodic boundary con-
ditions were used. The short range interactions were cut-off at 11 Å, whereas the long-range interactions
were handled by the Ewald method with the precision 10−4. After equilibration, several sets of simulation
runs, each for 6–8 ns, with restart from the previous configuration, were performed to obtain averages
for data analysis. The overall length of simulation was dictated by the stability of internal energy, density
and the dielectric constant, all being the functions of time.
13603-2
Density Anomaly of the SPC/E water model. Molecular Dynamics.
3. Results and discussion
The temperature dependence of the density of water and the stability of several ice phases in the
framework of the SPC/E model were studied in several works, see e.g. [15–21]. In figure 1 (a) we show the
results obtained by different authors for the SPC/E water and the results of the present work. The liquid
and glass branches, denominated in figure 1 (a) as Baez et al., were constructed using the points reported
in table II and in table III of reference [18]. The data follow from the simulations of 360 molecules. It
was assumed that the cut-off distance is at 9 Å and the effect of interactions at larger distances was
neglected. The density maximum (ρ = 1.0265 g/cm3) was obtained at 235 K. In the later work from that
laboratory, the density values for disordered hexagonal ice (Ih) in the interval between 170 K and 220 K,
were reported, see reference [20]. These data are also shown in figure 1 (a). However, it is important to
mention that the results refer to the modified SPC/E model, i.e., by multiplying the interaction potential
by a switching function between 6 Å and 9 Å to make the forces and the potential continuous at a cut-off
distance. The modification of the potential and the application of Ewald summation technique to include
the long-range interactions also resulted in the changes of the estimates for TMD, in [19, 20] the improved
values are ρ = 1.003 g/cm3 and T ≈ 260 K. The algorithmic peculiarities of these series of works include a
profound exploration of the effects of Ewald summation on the temperature dependence of density and
other properties of the system as well as calculations of the Gibbs free energy of distinct phases which
requires thermodynamic integration. The line describing the density of stable Ih ice obtained by Gay et
al. is reproduced from [15]. It was mentioned therein that the stability of this specific crystalline phase
crucially depends on the number of molecules in the simulation cell. It is worth mentioning that the
authors performed simulations in isostress, constant temperature ensemble and estimated the melting
point by analyzing the time evolution of the properties of the system upon heating of the solid phase.
Quite recently, Vega et al. [22] evaluated the melting temperature of the most popular non-polarizable
water models. For the SPC/Emodel, the melting temperature value of ice Ih and the densities of liquid wa-
ter and ice are shown in figure 1 (a). Our data for the liquid density of the model starting from a rather
high temperature down to the region of supercooled liquid are compared with the previously published
results by Bryk and Haymet [17] and with very recent results by Fennell et al. [21]. A small difference
between the curve obtained in [17] and the present results in the entire interval of temperatures is pre-
sumably caused by technical details, namely by the distinct cut-off value. On the other hand, it is worth
160 200 240 280 320 360
T (K)
940
960
980
1000
1020
ρ
(k
g/
m
3 )
Present work
Báez et al. (liquid)
Gay et al. (I
h
)
Báez et al. (glass) Fennell et al.
Experiment
Experiment (supercooled)
Bryk et al.
Arbuckle et al. (I
h
)
a
150 200 250 300 350 400
T (K)
920
940
960
980
1000
1020
ρ
(k
g/
m
3 )
SPC/E
TIP4P/2005
Experiment
ice I
h
b
Figure 1. (Color online) Panel a: Temperature dependence of the SPC/E model water density in the liquid
and glassy state as well as of Ih ice at a constant ambient pressure. The results of different authors are
marked in the figure. Panel b: A comparison of the density dependence of the SPC/E and TIP4P type
models. Circles show our results for the TIP4P-Ew model.
13603-3
E. Galicia-Andrés, H. Dominguez, O. Pizio
mentioning that Fennell et al. [21] slightly modified the original SPC/E model by fitting the dielectric con-
stant value at room temperature to the experimental result. Our data locate the TMD at 245 K, i.e., in
between the commonly used value 241 K, see e.g. [11, 17], and the result for the optimized model by Fen-
nell et al. 247 K. As for the value of density at TMD, our result is quite close to the data reported in [17].
In general terms, the water density values in relation to temperature grow in the interval that starts
at high temperatures and is delimited below by the TMD. Next, if the temperature decreases further, the
density decreases until it reaches a minimum value, Tmin, our data locate the minimum at T = 190 K. The
minimum value of density was observed for a set of water models, e.g., for ST2 [9] and TIP4P-2005 [23]. In
the case of the SPC/E model, the melting temperature Tm is located in the interval between the TMD and
Tmin, and the region between the melting point and Tmin actually corresponds to the supercooled states.
Finally, if the temperature decreases below the Tmin, the density grows again as in a normal fluid, but this
growth is not as rapid as in the interval of above the TMD.
In our opinion, a plausible explanation of the presence of the density minimum that involves struc-
tural changes in water was given in [9]. It refers to the behavior of a simple liquid that shrinks upon
cooling due to the increasing effects of attraction between particles with a decreasing temperature. How-
ever, in a certain temperature interval, this normal behavior alters due to augmenting effects of hydro-
gen bonding between molecules and the resulting formation of “complexes” or “agglomerates”. Coop-
erative behavior of these entities result in the growth of an expanded, open structure. The effects of
directional hydrogen bonding overcome the effects of an attractive interaction. Consequently, the fluid
density decreases in a certain interval of temperatures, evidently specific for each model. However, the
directional interactions are characterized by saturation. Therefore, upon further cooling, when changes
of the structural units become weaker, compared with the previous regime (as documented for example
by the changes of the average number of hydrogen bonds per molecule), then the fluid density grows
again similar to the temperature interval above the TMD.
However, the above discussion is incomplete and some thermodynamic arguments must be taken
into account. The liquid is the stable equilibrium phase for temperatures above the melting temperature
and it is unstable below a stability temperature, whereas in between, the liquid is metastable with respect
to crystal ice [24]. In general, metastable liquids and specifically water, can survive even at temperatures
well below the melting point until homogeneous nucleation takes place and water converts into ice [25].
Therefore, a certain part of the curve ρ(T ) presented in figure 1 (a) for the SPC/E model (below the melt-
ing point) characterized by a descending density can correspond to metastable, “long-life” states (on the
simulation time scale). Similar comments concern the behavior of the TIP4P-2005 model, see figure 4 of
reference [23] and the ST2 model [9, 24].
Due to finite-size effects and consequent details in implementation of Ewald sums as well as due to
finite sampling time, a possible spontaneous crystallization, however, is not observed in many works
for a set of water models, as well as in our present simulations. Debates concerning this problem just
started very recently, see e.g. [24, 25]. These authors provide a comprehensive discussion of various issues
relevant to the problem using their own data coming from the Monte Carlo and molecular dynamics
simulations and analyzing the previously published results for different watermodels from the literature.
Some other papers from the special issue of the Journal of Chemical Physics devoted to confined and
interfacial water published in November 2014 can be useful for the reader as well.
Our final comments concerning a set of curves shown in figure 1 (a) refer to the glassy water in
the framework of the SPC/E model. Actually, water, if considered below the homogeneous nucleation
temperature, should be in a glassy form, see e.g., a rather comprehensive discussion of the issue in [26].
In the case of the SPC/E model, the glass transition temperature was located around 177 K [18]. At this
temperature, the dependence of internal energy on temperature shows a kink (we reproduced the results
of that work by performing long simulations with 2000molecules and located the kink around T ≈ 185 K).
Thus, the augmenting density region at low temperatures refers in part to the amorphous states that
cannot be distinguished in terms of the structural characteristics such as the pair distribution functions
coming from the present simulations.
In figure 1 (b), we present the temperature dependence of density of the SPC/E model and compare
it with our results obtained for the TIP4P-Ew model performed in the framework of the above described
methodology. In addition, the results by H. Pi et al. [23] for the liquid and ice Ih phases are given as
well. The TIP4P-Ew model reproduces the experimental result for the TMD and an overall dependence of
13603-4
Density Anomaly of the SPC/E water model. Molecular Dynamics.
density at ambient pressure [27]much better compared to the SPC/E. All the models involved in this figure
(SPC/E, TIP4P-Ew, TIP4P-2005) predict the existence of the density minimum and augmenting density
at lower temperatures. However, it is worth mentioning that the density of ice Ih substantially differs
from the density of supercooled water (and of its glassy form) for the SPC/E model [see figure 1 (a)],
whereas in the case of TIP4P-2005 model, the liquid phase density and that of ice are close to each other
as it follows from the simulations by H. Pi et al. [23], see figure 1 (b). Similar behavior concerning the
density dependences has been documented recently by J. Alejandre et al. [28] for the model from the
TIP4P family developed by fitting the dielectric constant, the maximum density and the equation of state
at low pressure. Unfortunately, we are unaware of the simulation results for the Ih ice phase of TIP4P-Ew
model at room pressure to make our discussion more general at this point. However, it seems that the
differences between liquid and ice phase densities in the region of low temperatures are very sensitive
to the minor changes of parameters describing different models.
The microscopic structure of a set of water models coming out from MD simulations with respect to
the experimental diffraction data at room temperature and ambient pressure was discussed in great de-
tail by Pusztai et al. [29] by using the protocol [30] combining the experimental total scattering structure
factor and partial radial distribution functions from simulations. It was shown that the TIP4P-2005 model
exhibits the best consistency with the experimental structure compared to other commonly used mod-
els. Nevertheless, the SPC/E, ST2 and TIP4P results are not seriously inconsistent at ambient conditions.
Unfortunately, such an analysis of the structure of water has not been performed for other temperatures
and pressures due to experimental difficulties.
A few examples of the radial distribution functions, gi j , where i , j stand for O and H, obtained in
our NPT simulations of the SPC/E model are presented in two panels of figure 2. It can be seen that
the height of the first and the subsequent maxima of both functions, gOO(r ) and gOH, increases with a
decreasing temperature in the entire interval of temperatures studied. These trends are well pronounced
upon cooling in the range starting from a high temperature down to around 190 K, further. At even lower
temperatures, the overall shape of gi j (r ) is similar to that observed at 190 K. However, changes of the
structure are much less pronounced if onemoves to 180 K and 170 K.We omit these results for the sake of
brevity. The functions gOO(r ) and gOH(r ) are given mostly as an illustration, because their characteristics
are necessary to introduce the concept of hydrogen bonds between molecules.
In terms of the coordination number, nOO , it is appreciable that the cooling causes more pronounced
shells around the chosen molecule, cf. figure 3 (a). While at high temperatures even the first coordina-
2 4 6 8 10
r(Å)
0
1
2
3
4
5
g O
O
(r
)
SPC/E
T = 248 K
T = 190 K
T = 330 K
a
2 4 6 8
r(Å)
0
0.5
1
1.5
2
2.5
g O
H
(r
)
T = 298.15 K
SPC/E
T = 240 K
T = 200 K
T = 330 K
b
Figure 2. (Color online) Oxygen–oxygen (panel a) and oxygen–hydrogen (panel b) pair distribution func-
tions of the SPC/E model at different temperatures.
13603-5
E. Galicia-Andrés, H. Dominguez, O. Pizio
2.5 3 3.5 4
r(Å)
0
2
4
6
8
n O
O
(r
)
T = 200 K
T = 330 K
T = 298.15 K
T = 240 K
SPC/E
a
First coordination shell
2 4 6 8 10
r(Å)
0
10
20
30
40
dn
O
O
(r
)/
dr
T = 190 K
T = 330 K
T = 248 K
SPC/E
b
Figure 3. (Color online) Oxygen–oxygen coordination number of the SPC/E model at different temper-
atures (panel a). The derivative of the oxygen–oxygen coordination number at different temperatures
(panel b).
tion shell is not very well pronounced, it is seen perfectly well at low temperatures, and the coordination
number is around four as expected. The derivative of the coordination number is more illustrative [fig-
ure 3 (b)]. It shows the location of the first, the second shell and even the third one, the latter being
developed at very low temperatures and serving as a manifestation of a spatial order that develops in the
system and involves a larger number of particles than the first neighbors.
To proceed, it is necessary to introduce the concept of H-bonds. Energetic and geometric criteria, un-
avoidably involving certain degree of arbitrariness, are used in the literature to interpret the structure
of water and aqueous solutions. A comprehensive analysis of several geometric definitions of hydrogen
bonding, leading to different values for the average number of H-bonds per molecule, has been per-
formed in reference [31].
In this work, we use a single geometric criterion, see e.g. Zhang et al. [32] for two different water mod-
els. Three conditions determine the existence of the H-bond, namely the distance ROO between oxygen
atoms belonging to two water molecules should not exceed the threshold value RC
OO
; the distance ROH
between the acceptor oxygen atom and the hydrogen atom connected to the donor oxygen should not
exceed the threshold value RC
OH
. Finally, the angle O–H· · ·O should be smaller than the threshold value,
see e.g. references [33–37], concerning the analysis of the structure of liquid methanol as well as of su-
percritical water. If the coordinates of two molecules fulfill the above conditions, they are considered as
H-bonded. Zhang et al. proposed constant threshold values for the distances, 3.5 Å and 2.45 Å, for RC
OO
and RCOH, respectively. However, the values used in this work are the ones coming from the correspond-
ing radial distribution functions minima for O–O and O–H at each calculated temperature. In contrast to
[32], we used the threshold constant value for the angle θ, complementary to the angle α between the
vectors OH and rOH (see figure 1 of reference [31]), that is equal to 150◦. With this definition, we obtained
the average number of hydrogen bonds per water molecule similar to the value by Bakó et al. [38] who
used the energetic criterion for bonding.
The temperature dependence of the number of H-bonds averaged upon simulation time and normal-
ized by the number of water molecules in the system, 〈nHB〉Nw , following from our calculations for two
water models, the SPC/E and TIP4P-Ew, is shown in figure 4 (a). This figure explains that a decreasing
temperature results in a larger average number of H-bonds. It can be seen that there is a change of the
slope at T ≈ 180 K for the SPC/E model where the transition into a glassy state is reported. On the other
hand, the change of the slope is observed at a higher temperature, around T = 200 K, for the TIP4P-Ew
13603-6
Density Anomaly of the SPC/E water model. Molecular Dynamics.
200 240 280 320
T (K)
2.6
2.8
3
3.2
3.4
3.6
3.8
4
<
n H
B
>
N
w
SPC/E
a Bonds (O
w
1
H
w
2
+ H
w
2
O
w
1
)
TIP4P-Ew
TMD = 248 K
T
m
= 215 K TMD = 275 K
T
m
= 245.5 K
40 60 80 100 120 140 160
θ
0
0.01
0.02
0.03
0.04
0.05
P(
co
s
θ)
T = 200 K
T = 220.6 K
T = 240
T = 258 K
T = 298.15 K
T = 360 KSPC/E
107
o
53
o
b
0 1 2 3 4 5
bonding state of a water molecule
0
0.2
0.4
0.6
0.8
fr
ac
tio
n
of
d
if
fe
re
nt
ly
b
on
de
d
w
at
er
m
ol
ec
ul
es
T = 200 K
T = 240 K
T = 298.15 K
T = 330 K
SPC/E
c
160 200 240 280 320 360
T(K)
0
0.2
0.4
0.6
0.8
fr
ac
tio
ns
o
f
m
ol
ec
ul
es
in
d
if
fe
re
nt
b
on
di
ng
s
ta
te
s
4 bonds
3 bonds
2 bonds
1 bond
circles - TIP4P-Ew
triangles - SPC/E
d
TMD
TIP4P-Ew
TMD
SPC/E
Figure 4. (Color online) Temperature dependence of the average number of H-bonds per water molecule
(panel a). Distribution function of the angles formed by three oxygens on different molecules at different
temperatures (panel b). Fractions of water molecules in different bonding states on temperature (panel
c). A comparison of fractions of differently bonded molecules on temperature for the SPC/E and TIP4P-Ew
models (panel d).
model. However, this latter model is characterized by a higher number of H-bonds permolecule in the en-
tire temperature interval studied. In both cases we see that the average number of bonds does not reach
four, or, in other words, even at low temperatures, perfect tetrahedral configuration is not attained.
The same conclusion can be reached by considering the “bond-angle distribution function” already
studied by Bryk and Haymet [17]. It describes the probability to find certain configurations involving
three oxygens belonging to different water molecules. The threshold distance is assumed to be at the first
minimum of the pair distribution funcion, gOO(r ) at each temperature studied. The data were normalized
by the number of molecules that one finds in the shell limited by the threshold distance. All the curves
given in figure 4 (b) are characterized by twomaxima. One of them at θ ≈ 53◦ describes the configurations
characteristic of a simple liquid [16]. On the other hand, the second maximum at θ ≈ 107◦ is close to the
angle θ ≈ 109.5◦ describing a perfect tetrahedral configuration. At a high temperature, T = 360 K, a simple
liquid type structure prevails. While the system is cooled, the first maximum decreases and becomes very
13603-7
E. Galicia-Andrés, H. Dominguez, O. Pizio
small at T = 200 K, whereas themaximum describing the configurations close to tetrahedral substantially
grows in value. It is worth mentioning that both maxima slightly shift to higher angles with a decreasing
temperature. Still, even at low temperatures, the shape of the curves is not delta-like, witnessing a certain
degree of disorder in the system. Moreover, we do not observe any peculiarities of this property either at
TMD or at the temperature where the density is at minimum.
In order to obtain a better insight into the changes of the structure in terms of bonds with temper-
ature, we calculate the time-averages of fractions of differently bonded molecules. This representation
was already used by Bakó et al. [38] in their analysis of the simulated structures of water-methanol mix-
tures at room temperature. Our results concerning the SPC/E model at different temperatures are shown
in figure 4 (c). We observe that at a high temperature, T = 330 K, the fractions of molecules participating
in two and three bonds are the largest. Comparing this and room temperature, it can be seen that the
fraction of molecules with four bonds grows at the expense of a decreasing fraction of molecules with
a single bond and two bonds. The fraction of particles with three bonds remains practically unchanged.
A further decrease of temperature down to T = 240 K (note that it is very close to the TMD) leads to a
substantial increase of four-bonded particles again at the expense of singly-bonded and doubly-bonded
species. The fraction of triple-bonded species remains intact. If the temperature falls down to 200 K, not
only the fraction of 4-bonded particles grows substantially, but the fraction of triple-bonded species ex-
hibits a strong decrease in value. To summarize, the growth of density is accompanied by the growth of
the fraction of 4-bonded particles with almost constant fraction of triple-bonded species. On the other
hand, the descending density region, below the TMD, is characterized by the growth of 4-bonded species
with a simultaneous decrease of fractions of triple-bonded, double-bonded and single bonded-molecules.
Interestingly, the changes of the structure at even lower temperature are very small (we do not show
these curves to avoid overload of the figure).
The dependence of fractions of differently bonded molecules on temperature, in the entire tempera-
ture range studied, is shown in figure 4 (d). Here, we present our results for two different models, namely
for the SPC/E and TIP4P-Ew. The pattern is very similar in both cases. The TMD for each model obtained
from the calculations of density on temperature is located very close to the crossing points between
the curves describing the dependence of triple-bonded and 4-bonded molecules. The TMD values are
marked in the figure as well. It is worth mentioning that the curves describing the fraction of triple-
bonded molecules is not sensitive to the temperature variable in the region between high temperatures
and T ≈ 240 K. Moreover, the values for this fraction coming from the calculations of two models in ques-
tion practically coincide in an ample interval of temperatures. However, the difference between these
two models is well observed due to their capability of forming structures describing 4-bonded species.
In order to complete our discussion of the results given in all panels of figure 4, we would like to make
the following comments. As mentioned above, the geometric definition of H-bonds permits flexibility in
the sense that one can attempt to change the values and the set of parameters of H-bond definition. In this
respect, the energetic definition is more restrictive. Besides the choice of a set (ROO, rOH, θ), we also tested
the sets (ROO, rOH) and (ROO,β, cf. reference [31] figure 1). These less restrictive geometric definitions lead
to the results qualitatively similar to the ones discussed above. Namely, the fraction of quadruply-bonded
species grows with a decreasing temperature whereas all other fractions decrease. Unfortunately, we
have not observed any characteristic crossing attributable to the TMD with these definitions. Therefore,
we hope to extend our study of the models in question in the framework of energetic definition of H-
bonding as well.
The properties described above refer to the structure of the systems in question mostly at short inter-
particle separation. On the other hand, the long-range, asympthotic behavior of correlations between
molecules possessing dipole moment is described by the dielectric constant, ε, see e.g. [39]. In general
terms, the calculation of the dielectric constant from simulations is a difficult task [28, 40, 41], especially
at low temperatures. Several factors can influence the result, such as e.g., the number of molecules, the
type of thermostat and barostat, precision of the long-range interactions summation. Actually, very long
runs are necessary to obtain reasonable estimates for this property, as it is usually calculated from the
time-average of the fluctuations of the total dipole moment of the system [42] as follows:
ε= 1+
4π
3kBT V
(
〈M
2
〉−〈M〉
2
)
, (1)
13603-8
Density Anomaly of the SPC/E water model. Molecular Dynamics.
0 5 10 15 20 25
t (ns)
40
60
80
100
ε
330 K
298.15 K 268 K
240 K
SPC/E
a
0 5 10 15 20 25 30 35
t (ns)
40
60
80
100
ε
T = 260 K
T = 250 K
T = 275 K
TIP4P-Ew
T = 298.15 K
T = 270 K
b
0 4 8 12 16
t (ns)
-20
-10
0
10
20
Mα
M
x
M
y M
z
T = 298.15 K
SPC/E
c
240 280 320 360
T (K)
60
80
ε
SPC/E present work
Experiment
Fennell et al.
TIP4P-Ew
d
Figure 5. (Color online) Time dependence of the dielectric constant of the SPC/E and TIP4P-Ew water
model at different temperatures from our NPT molecular dynamics simulations, panels a and b, respec-
tively. Time dependence of the Cartesian projections of the total dipole moment of the system at room
temperature for the SPC/E model (panel c). Temperature dependence of the dielectric constant for the
SPC/E and TIP4P-Ew models from our calculations, from the fitting procedure by Fennell et al. [21] and
the experimental data [27].
where kB is the Boltzmann constant and V is the simulation cell volume. Some of our results at different
temperatures for the SPC/E model and for TIP4P-Ew are shown in figure 5 (a) and 5 (b), respectively. We
would like to recall here that all the data given in figure 5 were obtained for the systems having two
thousand particles. Practically all the runs lasted more than 20 ns.
Still, as it can be seen, there is room for experimenting in order to get better results below the room
temperature by extending the length of the runs or considering larger systems. In addition, we plan to ex-
plore the effects of changing the cutoff distance, of including the long-range corrections and the precision
of Ewald summation as well as to test other thermostats. All the results for each model, however, refer to
the temperatures above the respectivemelting point. Still, there is one delimiting factor. It manifests itself
in the polarization of the samples at low temperatures. At room temperature, figure 5 (c), the projections
of the total dipole moment of the system tend to zero at large times, whereas at lower temperatures, in
13603-9
E. Galicia-Andrés, H. Dominguez, O. Pizio
the supercooled regime and close to the formation of glassy states, the samples inevitably fall into con-
figuration with non-zero one or more Cartesian projections of M. Then, the calculations of the dielectric
constant require different algorithms. A summarizing insight into our data for ε(T ) for two models in
question are given in figure 5 (d). In general terms, the inclination of the temperature dependence of the
dielectric constant is correctly described by twomodels, compared to the experimental results. The fitting
of the dielectric constant at room temperature does not garantee the correct temperature dependence, as
we can see from the results of the modified model [21].
To summarize this brief report, we have used molecular dynamic simulations in the NPT ensemble
to study the density anomaly of water for the SPC/E and TIP4P-Ew models and to establish the relation
between this anomaly with the microscopic structure of the system. The structure is described in terms of
radial distribution functions, oxygen-oxygen coordination number and its derivative, the average num-
ber of H-bonds per molecule, angular distributions and fractions of molecules characterized by a differ-
ent number of bonds. We observed that at temperature of maximum density, the fractions of molecules
with three and with four bonds are almost the same, or in other words it is a reasonably good structural
estimator for the TMD. In addition, we presented our preliminary results concerning the temperature
dependence of the dielectric constant and compared them with experimental data. While the slope of
this dependence is well described, the values are underestimated compared to the experimental results.
It would be desirable in a future work to extend our structural analysis to find estimates for the cooper-
ative bonding of molecules and relate them with the density and other anomalies of water. On the other
hand, it would be interesting to establish the relationship between the structural properties describing
the electric properties of the system and the dielectric constant.
Acknowledgements
E.G. was supported by CONACyT of Mexico under Ph.D. scholarship. E.G. and O.P. are grateful to
Dr. T. Patsahan for very helpful discussions and valuable comments. O.P. is grateful to David Vazquez for
technical support of this work.
References
1. Barker J.A., Henderson D., Rev. Mod. Phys., 1976, 48, 587; doi:10.1103/RevModPhys.48.587.
2. Spohr E., Trokhymchuk A., Henderson D., J. Electroanal. Chem., 1998, 450, 281;
doi:10.1016/S0022-0728(97)00645-1.
3. Crozier P.S., Rowley R.L., Holladay N.B., Henderson D., Busath D.D., Phys. Rev. Lett., 2001, 86, 2467;
doi:10.1103/PhysRevLett.86.2467.
4. Crozier P.S., Henderson D., Rowley R.L., Busath D.D., Biophys. J., 2001, 81, 3077;
doi:10.1016/S0006-3495(01)75946-2.
5. Rahman A., Stillinger F.H., J. Am. Chem. Soc., 1973, 95, 7943; doi:10.1021/ja00805a003.
6. Jorgensen W.L., Chem. Phys. Lett., 1980, 70, 326; doi:10.1016/0009-2614(80)85344-9.
7. Lang E.W., Ludemann H.D., Angew. Chem. Int. Ed. Engl., 1982, 21, 315; doi:10.1002/anie.198203153.
8. Guillot B., J. Molec. Liq., 2002, 101, 219; doi:10.1016/S0167-7322(02)00094-6.
9. Brovchenko I., Oleinikova A., ChemPhysChem, 2008, 9, 2660; doi:10.1002/cphc.200800639.
10. Holten V., Bertrand C.E., Anisimov M.A., Sengers J.V., J. Chem. Phys., 2012, 136, 094507; doi:10.1063/1.3690497;
Preprint arXiv:1111.5587v2, 2011.
11. Vega C., Abascal J.L.F., Phys. Chem. Chem. Phys., 2011, 13, 19663; doi:10.1039/c1cp22168j.
12. Saul A., Wagner W., J. Phys. Chem. Ref. Data, 1989, 18, 1537; doi:10.1063/1.555836.
13. Horn H.W., Swope W.C., Pitera J.W., Madura J.D., Dick T.J., Hura G.L., Head-Gordon T., J. Chem. Phys., 2004, 120,
9665; doi:10.1063/1.1683075.
14. Forester T.R., Smith W., DL-POLY Package of Molecular Simulation, CCLRC, Daresbury Laboratory, Daresbury,
Warrington, England, 1996.
15. Gay S.C., Smith E.J., Haymet A.D.J., J. Chem. Phys., 2002, 116, 8876; doi:10.1063/1.1471556.
16. Bryk T., Haymet A.D.J., J. Chem. Phys., 2002, 117, 10258; doi:10.1063/1.1519538.
17. Bryk T., Haymet A.D.J., Mol. Simul., 2004, 30, 131; doi:10.1080/0892702031000152172.
18. Báez L.A., Clancy P., J. Chem. Phys., 1994, 101, 9837; doi:10.1063/1.467949.
13603-10
http://dx.doi.org/10.1103/RevModPhys.48.587
http://dx.doi.org/10.1016/S0022-0728(97)00645-1
http://dx.doi.org/10.1103/PhysRevLett.86.2467
http://dx.doi.org/10.1016/S0006-3495(01)75946-2
http://dx.doi.org/10.1021/ja00805a003
http://dx.doi.org/10.1016/0009-2614(80)85344-9
http://dx.doi.org/10.1002/anie.198203153
http://dx.doi.org/10.1016/S0167-7322(02)00094-6
http://dx.doi.org/10.1002/cphc.200800639
http://dx.doi.org/10.1063/1.3690497
http://arxiv.org/abs/1111.5587v2
http://dx.doi.org/10.1039/c1cp22168j
http://dx.doi.org/10.1063/1.555836
http://dx.doi.org/10.1063/1.1683075
http://dx.doi.org/10.1063/1.1471556
http://dx.doi.org/10.1063/1.1519538
http://dx.doi.org/10.1080/0892702031000152172
http://dx.doi.org/10.1063/1.467949
Density Anomaly of the SPC/E water model. Molecular Dynamics.
19. Báez L.A., Clancy P., J. Chem. Phys., 1995, 103, 9744; doi:10.1063/1.469938.
20. Arbuckle B.W., Clancy P., J. Chem. Phys., 2002, 116, 5090; doi:10.1063/1.1451245.
21. Fennell C.J., Li L., Dill K.A., J. Phys. Chem. B, 2012, 116, 6936; doi:10.1021/jp3002383.
22. Vega C., Sanz E., Abascal J.L.F., J. Chem. Phys., 2005, 122, 114507; doi:10.1063/1.1862245.
23. Pi H.L., Aragones J.L., Vega C., Noya E.G., Abascal J.L.F., Gonzalez M.A., McBride C., Mol. Phys., 2009, 107, 365;
doi:10.1080/00268970902784926.
24. Limmer D.T., Chandler D., J. Chem. Phys., 2013, 138, 214504; doi:10.1063/1.4807479.
25. Espinoza J.R., Sanz E., Valeriani C., Vega C., J. Chem. Phys., 2014, 141, 18C529; doi:10.1063/1.4897524.
26. Mallamace F., Corsaro C., Mallamace D., Vasi S., Vasi C., Stanley H.E., J. Chem. Phys., 2014, 141, 18C504;
doi:10.1063/1.4895548.
27. CRC Handbook of Chemistry and Physics, 95th Edition, Haynes W.R. (Ed.) CRC Press, Boca Raton, 2014.
28. Alejandre J., Chapela G.A., Saint-Martin H., Mendoza N., Phys. Chem. Chem. Phys., 2011, 13, 19728;
doi:10.1039/c1cp20858f.
29. Pusztai L., Pizio O., Sokolowski S., J. Chem. Phys., 2008, 129, 184103; doi:10.1063/1.2976578.
30. Pusztai L., Harsanyi I., Dominguez H., Pizio O., Chem. Phys. Lett., 2008, 457, 96; doi:10.1016/j.cplett.2008.03.091.
31. Kumar R., Schmidt J.R., Skinner J.L., J. Chem. Phys., 2007, 126, 204107; doi:10.1063/1.2742385.
32. Zhang N., Li W., Chen C., Zuo J., Weng L., Mol. Phys., 2013, 111, 939; doi:10.1080/00268976.2012.760050.
33. Padró J.A., Saiz L., Guàrdia E., J. Mol. Struct., 1997, 416, 243; doi:10.1016/S0022-2860(97)00038-0.
34. Guàrdia E., Martí J., Prado J.A., Saiz L., Vanderkooi A.V., J. Mol. Liq., 2002, 96–97, 3;
doi:10.1016/S0167-7322(01)00342-7.
35. Guàrdia E., Martí J., García-Tarrés L., Laria D., J. Mol. Liq., 2005, 117, 63; doi:10.1016/j.molliq.2004.08.004.
36. Skarmoutsos I., Dellis D., Samios J., J. Phys. Chem. B, 2009, 113, 2783; doi:10.1021/jp809271n.
37. Skarmoutsos I., Guàrdia E., J. Chem. Phys., 2010, 132, 074502; doi:10.1063/1.3305326.
38. Bakó I., Megyes T., Bálint S., Grósz T., Chihaia V., Phys. Chem. Chem. Phys., 2008, 10, 5004; doi:10.1039/b808326f.
39. Yukhnovskyj I.R., Holovko M.F., Statistical Theory of Classic Equilibium Systems, Naukova Dumka, Kiev, 1980 (in
Russian).
40. Rami Reddy M., Berkowitz M., Chem. Phys. Lett., 1989, 155, 173; doi:10.1016/0009-2614(89)85344-8.
41. Gereben O., Pusztai L., Chem. Phys. Lett., 2011, 507, 80; doi:10.1016/j.cplett.2011.02.064.
42. Neumann M., Mol. Phys., 1983, 50, 841; doi:10.1080/00268978300102721.
Температурна залежнiсть мiкроскопiчної структури i
аномалiя густини в моделях води SPC/E та TIP4P-Ew.
Результати моделювання методом молекулярної динамiки
E. Галiсiя-Андрес1, Г. Домiнгез2, О. Пiзiо1
1 Iнститут хiмiї, Нацiональний автономний унiверситет Мексики, Мехiко, Мексика
2 Iнститут дослiдження матерiалiв, Нацiональний автономний унiверситет Мексики, Мехiко, Мексика
Ми дослiдили температурнi залежностi мiкроскопiчної структури моделей води SPC/E i TIP4P-Ew в тер-
мiнах парних функцiй розподiлу, координацiйних чисел, середнього числа водневих зв’язкiв, розподiлу
зв’язаних станiв однiєї молекули, використовуючи метод молекулярної динамiки при постiйногму тиску.
Еволюцiя структури поставлена у вiдповiднiсть до залежностi густини води вiд температури в областi,що
знаходиться мiж високими температурами i температурами, де система стає переохолодженою. Показа-
но, що частка молекул з трьома i чотирма зв’язками визначає максимальну густину для обох моделей.
Бiльш того, отримано i проаналiзовано температурну залежнiсть дiелектричної сталої.
Ключовi слова: моделi води, аномалiя густини, молекулярна динамiка
13603-11
http://dx.doi.org/10.1063/1.469938
http://dx.doi.org/10.1063/1.1451245
http://dx.doi.org/10.1021/jp3002383
http://dx.doi.org/10.1063/1.1862245
http://dx.doi.org/10.1080/00268970902784926
http://dx.doi.org/10.1063/1.4807479
http://dx.doi.org/10.1063/1.4897524
http://dx.doi.org/10.1063/1.4895548
http://dx.doi.org/10.1039/c1cp20858f
http://dx.doi.org/10.1063/1.2976578
http://dx.doi.org/10.1016/j.cplett.2008.03.091
http://dx.doi.org/10.1063/1.2742385
http://dx.doi.org/10.1080/00268976.2012.760050
http://dx.doi.org/10.1016/S0022-2860(97)00038-0
http://dx.doi.org/10.1016/S0167-7322(01)00342-7
http://dx.doi.org/10.1016/j.molliq.2004.08.004
http://dx.doi.org/10.1021/jp809271n
http://dx.doi.org/10.1063/1.3305326
http://dx.doi.org/10.1039/b808326f
http://dx.doi.org/10.1016/0009-2614(89)85344-8
http://dx.doi.org/10.1016/j.cplett.2011.02.064
http://dx.doi.org/10.1080/00268978300102721
Introduction
Simulation details
Results and discussion
|