Time-dependent density fluctuations in liquid water

Chemical Physics Letters 649 (2016) 119–122

Contents lists available at ScienceDirect

Chemical Physics Letters journal homepage: www.elsevier.com/locate/cplett

Time-dependent density fluctuations in liquid water Conor J. Waldron 1 , Niall J. English ∗ School of Chemical and Bioprocess Engineering, University College Dublin, Belfield, Dublin 4, Ireland

a r t i c l e

i n f o

Article history: Received 24 November 2015 In final form 16 February 2016 Available online 22 February 2016

a b s t r a c t Temporal system-mass-density fluctuation analysis was performed on liquid−water moleculardynamics simulations at ambient pressure and 200 and 300 K, in three increasingly-large systems. A prominent mode in system-density fluctuations was observed at molecular-librational frequencies of ∼600−800 cm−1 (with pronounced temperature dependence). This mode displayed marked system-size dependence, disappearing for larger systems. Persistent system-density fluctuations were clearly evident at 10−11 cm−1 for all systems and temperatures, with lower-amplitude ‘overtones’ evident only in larger systems. It is conjectured that this reflects ∼3 ps timescales observed in earlier studies for dissipation of local-density fluctuations in liquid water in this 200−300 K temperature range. © 2016 Elsevier B.V. All rights reserved.

1. Introduction Water has various anomalous, intriguing characteristics which cannot be explained easily from a simple liquid’s perspective [1]. The mechanisms underpinning mass-density fluctuations constitute one such rather formidable challenge, where we refer both to spatial and temporal variations in density, and density can be both ‘local’ (in the sense of local regions of liquid within a ‘correlation length’) or defined for the system as a whole (i.e. the aggregate mass of molecules per unit of system volume, m/V, where m is the collective mass of the total number of molecules in the system, and V is the system volume). Here, a correlation length refers to the distance over which spatial fluctuations in local density are correlated, in the sense that local density has both spatial and temporal fluctuations as characterised by normal thermodynamic fluctuation-dissipation in a liquid [2]. This correlation length can be evaluated from the spatial autocorrelation function of local density, sampled over a sufficiently long period of time to gain adequate statistics to overcome temporal fluctuations [2]. The properties of water are explained in large part by the intricacies of the behaviour of hydrogen bonds [3–8], and molecular dynamics (MD) has contributed to our understanding of this in the ambient liquid and supercooled state, and in terms of strain in hydrogen-bonded polygons [3,9]. Hydrogen-bond orientational

∗ Corresponding author. E-mail address: [email protected] (N.J. English). 1 Present address: Department of Chemical Engineering, University College London, London, United Kingdom. http://dx.doi.org/10.1016/j.cplett.2016.02.038 0009-2614/© 2016 Elsevier B.V. All rights reserved.

mechanisms have been studied explicitly by MD in terms of a ‘jump’ mechanism explaining water reorientation via hydrogenbond cleavage and molecular orientation occurring concertedly [6,7]. Indeed, these fundamental properties of collective hydrogenbond reorientation and jumps have been linked to variations in local density [8]. In any event, quite apart from the question of local, or indeed system, mass density per se (and its intriguing properties such as its ambient-pressure maximum at 4 ◦ C [10]), the temporal fluctuations thereof in liquid water offer an intriguing glimpse into the rich tapestry of phenomena determining water behaviour. In particular, English and Tse [2] have studied local-density fluctuations in liquid water by million-molecule MD for both ambient and supercooled conditions, establishing spatial correlation length for local density and timescales of local high- and low-density regions; compared with smaller, historic sizes, they found finite-size effects are problematic in describing density fluctuations in small systems, due to phonon wavelengths permitted by artificially small size and periodic boundary condition (PBC) restrictions [2,11,12]. Local-density fluctuations were found to occur in liquid water within timescales of ∼3 ps, largely independent of temperature and system size [2]. English et al. have studied system-size effects via MD on methane hydrates [13], planar liquid methane−water interfaces [14] and ice crystallisation and melting [15] and early-onset nucleation [16] on systems with up to ∼8 million molecules, finding omission of (lower-frequency) vibrations in smaller systems led to artificial suppression of methane-hydrate precursor [14] and ice crystal [15] and pre-cursor [16] formation. In all cases [2,11–16], the role of collective lower-frequency modes (the manifold of which becomes more accessible in larger systems) was highlighted, i.e., essentially

120

C.J. Waldron, N.J. English / Chemical Physics Letters 649 (2016) 119–122

system-density fluctuations were suggested as being influential on these other phenomena, in an important interplay with largely local properties and behaviour. Given the importance of density fluctuations in determining water’s unique behaviour, it is clear that a more rigorous investigation of system-density fluctuations in liquid water would certainly elucidate unambiguously the nature and behaviour of these modes, and the effect of system size thereon, and clarify effects on localdensity fluctuations. In any event, owing to previous SPC-type models in local-density fluctuation analysis in liquid water [2], we use SPC [17] for liquid water at two temperatures, 200 and 300 K, above the melting point of 190 K [18], for a variety of increasing system sizes (188, 4096 and 39788 molecules) so as to study temperature and size effects, but not in the supercooled region. We do so to avoid any potential ‘entangling’ of density fluctuations with the development of ice-like pre-cursors, given that these are separate and distinct physical processes [16]. We also wish to avoid examining potential ‘two-liquid’, supercooled water, given that this is outside the ambit of this study, with our focus on equilibrium time-dependent system-density fluctuations in (single-phase) liquid water. 2. Methodology Equilibrium MD was performed in cubic simulation boxes for 1 ns, for the various temperatures and system sizes, under constantpressure, constant-temperature (NPT) conditions at ambient (1 bar) pressure. However, in an effort to mitigate any possibility for the choice of barostat or inertia parameters affecting the volumedilation dynamics of the system, we employed a formulation of Gauss’s principle of least constraint applied to both temperature and pressure, so that these are inherent constants of motion [19]. Briefly, the Newtonian equations of motion for r (position) and p (momentum) of site i in this formulation is altered as follows [19]: r˙ i =

p˙ i ˙ i + εr mi

(1)

V˙ = 3V ε˙

(2)

˙ i − ςpi p˙ i = fi − εp

(3)

where the force is fi and the thermostat’s inertia parameter, ␨, is given by N 

ς = −ε˙ +

f i · pi

i=1 N 

(4) p2i

i=1

˙ is specwhilst the system-volume dilation (or contraction) rate, ε, ified by 1 m

ε˙ =

N   

rij · pij



 

ij + ij / rij  



i=1 j>i N   

r2ij

  ij + ij / rij  + 9PV 

(5)



i=1 j>i

In Eqs. (4) and (5), it is understood that the summations run over   molecules i and j, rather than individual sites, while ij and ij refer, respectively, to the first and second derivatives of the potential. The mass m refers to that of the molecule, while P is a pre-defined constant (1 bar) of the motion, but V varies dynamically (according to Eqn. (2)). Also, rij = ri − rj and pij = pi − pj , and periodic boundary

conditions and holonomic constraints are applied to all position calculations. The box lengths of initial systems were such to conform to a system density of 1.00 g/cm3 . The Particle-Mesh Ewald method was used to handle long-range electrostatic interactions [20]. Normalised autocorrelation functions (ACFs) of the time-derivatives of the system mass density were computed:

 ˛˙ i (t)˛˙ i (0) c (t) =  ˛

˛˙ i (0)˛˙ i (0)

(6)

where ␣ denotes density (m/V, i.e. total mass molecules per unit of system volume). The use of time-derivatives, ˛, ˙ in the definition has the advantage of being acutely sensitive to temporal fluctuations (e.g. periodic oscillations or ‘vibrations’) in the system density, allowing this to be probed conveniently and straightforwardly via power spectra (Fourier-cosine transformation) of their ACFs. 3. Results and discussion Typical behaviour of the system density is depicted in Fig. 1, displaying prominent oscillations with a period of ∼3−3.5 ps, regardless of the temperature or system size. Resultant ACFs of the system-density time derivative (from Eqn. (6)) displayed strong, undamped and sustained harmonic oscillations with corresponding periods of ∼3.0−3.5 ps, and this was found in all of the ACFs, regardless of the temperature or system size. However, superimposed on these oscillations, particularly for smaller systems, were much shorter-period vibrations, characteristic of individual molecular librations (rotation oscillations) [21,22]. The peaks present in the power spectra of the system-density-derivative ACFs, depicting the dominant modes, are summarised in Table 1. From Table 1, it is clear that a prominent peak in density occurs at molecularlibrational frequencies of ∼600−800 cm−1 , with higher frequencies at 200 K owing to slower hydrogen-bonding kinetics [4] leading to a greater degree of strain in hydrogen-bonding networks [3] that dominate (local) molecular libration [9]. However, this ‘librational’ density-fluctuation mode displayed a marked system size dependence, disappearing for larger system sizes, as individual molecular librations no longer affect fluctuations in system density, given the diminishing effect of local, near-neighbour phenomena on the collective nature (i.e. system density) of increasingly larger systems. The persistent system-density fluctuations evident in Fig. 1 show up clearly at circa 10−11 cm−1 in power spectra for all system sizes and essentially independent of temperature, characteristic of a ∼3.0−3.5 ps system-density oscillation. However, given the consistency of the peak positions (both in this general ∼10−11 cm−1 range, this suggests that artificial coupling of the system density with the external reservoirs is not a serious problem. In particular, lower-amplitude, double, triple and quadruple ‘overtones’ are evident at up to ∼40 cm−1 only in larger systems (cf. Fig. 2) as the manifold of available collective-vibrational mode ‘echoes’ becomes more extensive with increasing simulation-box lengths, as argued in ref. 2, and weak-amplitude standing waves develop. It is conjectured that the ‘fundamental’ mode at ∼10−11 cm−1 reflects the timescale of ∼3−3.5 ps observed in earlier MD studies of for the dissipation of local- (as opposed to system-) density fluctuations in liquid water in this 200−300 K temperature range [2]. This provides evidence of the strong coupling of collective system-density vibrations with more local ones in affecting the frequency and timescales of local-density fluctuations (studied, inter alia, in ref. 2). In addition, more recently, Santra et al. [23] have uncovered ∼4 ps persistence times in local structural quantities in ab initio MD simulation of liquid water, which they noted is similar to relaxation times for local density fluctuations in ref. 2, and we note here that this is effectively identical to those noted for system-density

C.J. Waldron, N.J. English / Chemical Physics Letters 649 (2016) 119–122

121

Fig. 1. Instantaneous densities of 39,788-molecule system at 200 K after ∼1.2 ns, showing the ‘beating’ pattern with circa 10−11 cm−1 frequency (or ∼3 ps period). Density of 39 788-molecule system at 200 K.

Table 1 Typical peak positions of power spectra of the system-density-derivative ACF (Eqn. 6). Temperature K

System size No. molecules

System-density peak cm−1

Libration peak cm−1

200 300 200 300 200 300

188 188 4096 4096 39788 39788

9 9 9 9 9, 17, 27, 35 9, 18, 28, 36

759 585 732 587 none none

Fig. 2. Power spectra (from Fourier-cosine transformation of density-derivative autocorrelation functions) for 39,788-molecule system at 200 K. Modes at ∼10−11, 20 30 and 40 cm−1 frequency are present, with the 40 cm−1 less evident/visible (without magnification). Spectra of density ACF for 39,788-molecule system at 200 K. Modes at ∼10−11, 20 30 and 40 cm−1 frequency are present.

(as opposed to local density) in the present work. Further, Ruocco and Sette [24] have reviewed and discussed with some acuity the origin of collective dynamics in liquid water, in view of MD and inelastic x-ray scattering (IXS) studies; they have described two collective modes in the dynamical structure factor which correspond to the apparent fast-sound velocity of ∼3200 ms−1 and another at constant, lower energy of ∼30−50 cm−1 (essentially independent of temperature), reflective of ∼1500 ms−1 sound velocity. In particular, Rahman and Stillinger [25] determined the lower-, constant- energy mode at ∼30−40 cm−1 from MD simulation with the ST2 potential: this is dissimilar to the presently-used SPC [17] potential, in that ST2 overestimates tetrahedral structure in liquid water with a temperature of maximum density (TMD) of around 320 K and a freezing point of ∼300 K [26]. Ricci et al. [27] have indicated that librational modes above ∼400 cm−1 may possibly be associated with fast-sound waves. In any event, the temperature-independent ∼3−3.5 ps period of system-density fluctuations (and associated ∼10−11 cm−1 modes) observed in the present work with SPC, which is an under-structured model, exhibiting its freezing point and TMD, and other properties at lower temperatures than real water [18], along with the clear influence of temperature-dependent librational modes at 600−800 cm−1 , mirrors only qualitatively the two collective modes identified by IXS structure-factor measurements. The difference between ST2 and SPC models, with their respective over- and under- structuring propensities, renders direct quantitative accuracy more problematic in terms of direct comparison with experimental data. Certainly, the speed of sound at ambient temperature in understructured SPC water of 1635 ms−1 is not in as good accord with the experimental value of 1497 ms−1 [28] as would be desired, ideally. The prominent low-frequency mode identified in the present work is not reflective of the speed of sound, either apparent or real, directly, in that it is clearly independent of system size and temperature, and so cannot be regarded as a standing, or stationary, wave per se. Still, this collective-system-density mode at ∼10−11 cm−1 appears to reflect directly local-density ∼3−3.5 ps relaxation.

122

C.J. Waldron, N.J. English / Chemical Physics Letters 649 (2016) 119–122

4. Conclusions System-mass-density fluctuation analysis displays a prominent vibrational mode at molecular-librational frequencies of ∼600−800 cm−1 (with pronounced temperature dependence). There is also a marked system-size dependence, disappearing for larger system sizes, as individual molecular librations, arising from the hydrogen-bonding network, no longer affect system-density fluctuations. Its marked temperature dependence is not unexpected, given the intimate connection with local hydrogen-bonding behaviour. The persistent system-density fluctuations were clearly evident at circa 10−11 cm−1 for all system sizes were essentially independent of temperature, with lower-amplitude, double, triple and quadruple ‘overtones’ evident at up to ∼40 cm−1 only in larger systems as the manifold of available collective-vibrational modes becomes more extensive. It was conjectured that the ‘fundamental’ mode at 10−11 cm−1 reflects the timescale of ∼3 ps observed in earlier MD studies of for the dissipation of local (as opposed to system) density fluctuations in liquid water in this 200−300 K temperature range [2]: collective vibrations appear to affect both the local [2] and system fluctuations. References [1] V. Holten, C.E. Bertrand, M.A. Anisimov, J.V. Sengers, J. Chem. Phys. 136 (2012) 094507.

[2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28]

N.J. English, J.S. Tse, Phys. Rev. Lett. 106 (2011) 037801. F.H. Stillinger, Science 209 (1980) 451. A. Luzar, D. Chandler, Nature 379 (1996) 55. N.J. English, Mol. Phys. 103 (2005) 1945. F.W. Starr, J.K. Nielsen, H.E. Stanley, Phys. Rev. Lett. 82 (1999) 2294. D. Laage, J.T. Hynes, Science 311 (2006) 832. J.T. Titantah, M. Karttunen, Sci. Rep. 3 (2013) 2991. S.N. Chakraborty, N.J. English, J. Chem. Phys. 143 (2015) 154504. M. Yasutomi, Front. Phys. 3 (2015) 8. N.J. English, J.S. Tse, Phys. Rev. Lett. 103 (2009) 015901. N.J. English, J.S. Tse, Phys. Rev. B 80 (2009) 134306. N.J. English, Mol. Phys. 106 (2008) 1887. N.J. English, M. Lauricella, S. Meloni, J. Chem. Phys. 140 (2014) 204714. N.J. English, J. Chem. Phys. 141 (2014) 234501. N.J. English, J.S. Tse, Phys. Rev. E 92 (2015) 032132. H.J.C. Berendsen, J.P.M. Postma, W.F. van Gunsteren, J. Hermans, in: B.P. Reidel (Ed.), In Intermolecular Forces, Dordrecht, 1981, p. 331. C. Vega, E. Sanz, J.L. Abascal, J. Chem. Phys. 122 (2005) 114507. D.J. Evans, G.P. Morriss, Statistical Mechanics of Nonequilibrium Liquids, Academic Press, London, 1990. U. Essmann, L. Perera, M.L. Berkowitz, T. Darden, H. Lee, L.G. Pedersen, J. Chem. Phys. 103 (1995) 8577. D. Bertolini, A. Tani, Mol. Phys. 75 (1992) 1047. I.M. Svishchev, P.G. Kusalik, J. Chem. Soc. Faraday Trans. 90 (1994) 1405. B. Santra, R.A. DiStasio Jr., F. Martelli, R. Car, Mol. Phys. 113 (2015) 2829. G. Ruocco, F. Sette, J. Phys.: Condens. Matter 11 (1999) R259. A. Rahman, F.H. Stillinger, Phys. Rev. A 1 (1974) 368. P.H. Poole, R.K. Bowles, I. Saika-Voivod, F. Sciortino, J. Chem. Phys. 138 (2013) 034505. M.A. Ricci, M.D. Rocca, G. Ruocco, R. Vallauri, Phys. Rev. A 40 (1989) 7266. I. Shvab, R.J. Sadus, J. Chem. Phys. 139 (2013) 194505.