# Dynamics of simulated water under pressure

## I Introduction

The “slow dynamics” and glass transition of both simple and molecular liquids has been a topic of significant interest in recent years. The initial slowing down of liquids at temperatures down to $`T_c1.2T_g`$, where relaxation times approach $`1`$ ns, has been well described by the mode-coupling theory (MCT) . MCT has been successfully applied to a wide variety of real and model systems , including hard spheres , Ni<sub>80</sub>P<sub>20</sub> , SiO<sub>2</sub> , and polymer melts . However, there has not been an extensive test of the validity of the MCT predictions for a model system over a wide range of pressures and along different thermodynamic paths.

At low pressure, it was shown previously that the power-law behavior of dynamic properties in the SPC/E model  can be explained using MCT . Furthermore, the possible relationship between the experimentally-observed power-law behavior and the predictions of MCT has been discussed . The experimentally-observed locus of apparent power-law singularities of dynamic and thermodynamic properties \[Fig. 1\] is of particular interest , and has catalysed the development of three scenarios to explain the anomalous properties of water: (i) the existence of a spinodal bounding the stability of the liquid in the superheated, stretched, and supercooled states ; (ii) the existence of a liquid-liquid phase transition line separating two liquid phases differing in density ; (iii) a singularity-free scenario in which the thermodynamic anomalies are related to the presence of low-density and low-entropy structural heterogeneities . The predictions of MCT are of interest since MCT might account for the apparent power-law behavior of dynamic properties on cooling, thereby removing the need for a thermodynamic explanation of the dynamic properties of water.

In this article we focus on two related issues: (i) the possibility of using MCT to explain the slow dynamics of water under pressure, and (ii) a test of the validity of the MCT predictions over an extremely wide pressure range in a system with dramatic structural changes. We find that MCT provides a good account of the slow dynamics of the SPC/E model for water at all pressures, with the structure evolving continuously from an open tetrahedral network to a densely packed fluid, similar to a Lennard-Jones type liquid. By examining the wave-vector dependence of collective dynamics, we are able to discover how these structural changes are reflected in the dynamic behavior of the liquid. We are also able to test the validity of the relationship predicted by MCT for the diffusivity exponent $`\gamma `$ and the von Schweidler exponent $`b`$ over a wide range of values $`\gamma `$ and $`b`$ \[Fig. 2\]. Our results support the predicted relationship of these exponents. A brief report of a subset of the present results for the SPC/E potential has recently appeared . The dynamic properties of the ST2 potential  at one pressure, and of the TIP4P  potential at several pressures, have also recently been discussed.

## II Mode Coupling Theory

We will focus our discussion on the idealized form of MCT, originally formulated to describe spherically-symmetric potentials. Recent extensions have been made to account for the rotational motion present in non-spherical molecular systems , such as water. The idealized version of MCT has been shown to provide a good account for the center-of-mass motion for the SPC/E model . We provide only a brief account of the MCT predictions relevant to the results of this article, and we refer the reader to extensive reviews for more information .

MCT assumes that localization, or “caging” of molecules due to the slow rearrangement of neighboring molecules, is the source of the dramatic increase of relaxation times on cooling, leading to a strong coupling between single particle motion and the density fluctuations of the liquid. Indeed, according to MCT, the static density fluctuations, measured by the structure factor $`S(q)`$, entirely determine the long time dynamic behavior. MCT accounts for the loss of correlation by the interaction of density mode fluctuations, ignoring other possible mechanisms for relaxation. MCT predicts the asymptotic power-law divergence of correlation times, and power-law vanishing of the diffusion constant

$$DD_0(T/T_c1)^\gamma $$
(1)   
at a critical temperature $`T_c=T_c(P)`$, where we refer to $`\gamma =\gamma (P)`$ as the diffusivity exponent. In real systems, “freezing” of the system dynamics is avoided at $`T_c`$, as

relaxation mechanisms not accounted for by MCT become significant. However, $`T_c`$ can still be interpreted as a “crossover temperature” where the dynamics change from being dominated by density fluctuations to being controlled by “activated” processes. Some recent work has also demonstrated the significance of $`T_c`$ as a crossover temperature where relaxation occurs primarily through basin hopping , in the energy landscape view of liquid dynamics .

MCT predicts that the Fourier transform of the density-density correlation function  or intermediate scattering function

$$F(q,t)\frac{1}{S(q)}\underset{j,k=1}{\overset{N}{}}e^{i𝐪[𝐫_k(t)𝐫_j(0)]}$$
(2)   
decays via a two-step process. In the first relaxation step, $`F(q,t)`$ approaches a plateau value $`F_{\text{plateau}}(q)`$ which is described, to leading order in time, by a power law with exponent $`a=a(P)`$,

$$F(q,t)F_{\text{plateau}}(q)t^a,$$
(3)   
At larger times, $`F(q,t)`$ decreases from $`F_{\text{plateau}}(q)`$ and MCT predicts the decay obeys the von Schweidler power law to leading order in time

$$F_{\text{plateau}}(q)F(q,t)t^b,$$
(4)   
where $`b=b(P)`$ is known as the von Schweidler exponent. The region of validity of Eqs. (3) and (4) can be quite limited.

The slow relaxation of $`F(q,t)`$ has a characteristic relaxation time $`\tau `$ that is also predicted to have asymptotic power law dependence on temperature,

$$\tau \tau _0(T/T_c1)^\gamma $$
(5)   
with the same value of the exponent $`\gamma `$ as for the diffusion constant. Hence, Eqs. (1) and (5) predict that the product $`D\tau `$ is not singular as $`TT_c`$, hence we take the product to be constant over the range that Eqs. (1) and (5) are valid (neglecting corrections to scaling).

MCT predicts that the scaling exponents $`a`$, $`b`$, and $`\gamma `$ are not independent; $`a`$ and $`b`$ are related by the exponent parameter $`\lambda `$ using the relationship

$$\lambda =\frac{[\mathrm{\Gamma }(1a)]^2}{\mathrm{\Gamma }(12a)}=\frac{[\mathrm{\Gamma }(1b)]^2}{\mathrm{\Gamma }(1+2b)}$$
(6)   
where $`\mathrm{\Gamma }(x)`$ is the gamma function. MCT also relates $`\gamma `$ to $`a`$ and $`b`$ via

$$\gamma =\frac{1}{2a}+\frac{1}{2b}.$$
(7)   
Because of Eqs. (6) and (7), only one exponent value is needed to determine all others, so calculation of two exponents determines if the dynamics of a system are consistent with the predictions of MCT. Furthermore, these exponents are expected to depend on the path along which $`T_c`$ is approached.

After $`F(q,t)`$ departs from the plateau, $`F(q,t)`$ is well-described by a Kohlrausch-Williams-Watts stretched exponential

$$F(q,t)=A(q)\mathrm{exp}\left[\left(\frac{t}{\tau (q)}\right)^{\beta (q)}\right],$$
(8)   
where $`\tau (q)`$ is the relevant relaxation time. Moreover, it has been shown that the exponent $`\beta =\beta (q)`$ is related to the von Schweidler exponent 

$$\underset{q\mathrm{}}{lim}\beta (q)=b.$$
(9)   
This relation facilitates evaluation of $`b`$, since the region of validity of Eq. (4) is difficult to identify in practice.

## III Simulations

We perform MD simulations of 216 water molecules interacting via the SPC/E pair potential . The SPC/E model treats water as a rigid molecule consisting of three point charges located at the atomic centers of the oxygen and hydrogen, which have an OH distance of 1.0 Å and HOH angle of 109.47, the tetrahedral angle. Each hydrogen has charge $`q_H=0.4238e`$, where $`e`$ is the fundamental unit of charge, and the oxygen has charge $`q_O=2q_H`$. In addition, the oxygen atoms of separate molecules interact via a Lennard-Jones potential with parameters $`\sigma =3.166`$ Å and $`ϵ=0.6502`$ kJ/mol.

Our simulation results are summarized in Table I. For $`T300`$ K, we simulate two independent systems to improve statistics, as the long relaxation time makes time averaging more difficult. We equilibrate all simulated state points to constant $`T`$ and $`\rho `$ by monitoring the pressure and internal energy. We control the temperature using the Berendsen method of rescaling the velocities , while the reaction field technique with a cutoff of 0.79 nm  accounts for the long-range Coulombic interactions. The equations of motion evolve using the SHAKE algorithm  with a time step of 1 fs, except at $`T=190`$, where a time step of 2 fs is used due to the extremely slow motion of the molecules. Equilibration times at high temperatures are relatively small. At low $`T`$, extremely long equilibration times are needed. The structural and thermodynamics properties may be obtained after relatively short equilibration times. However, dynamic properties show significant aging effects (i.e. dependence of measured properties on the chosen starting time) if great care is not taken in equilibration.

For production runs, it is desirable to make measurements in the isoenergetic/isochoric ensemble (NVE). However, a small energy drift is unavoidable for the long runs presented here, so we again employ the heat bath of Berendsen, using a relaxation time of 200ps . The large relaxation time prevents an energy drift but achieves results that are very close to those that would be found if it were possible to perform a simulation in the NVE ensemble.

Since we perform long runs for many state points, we store the molecular trajectories $`\{𝐫_i,𝐩_i\}`$ at logarithmic intervals to avoid storage problems that linear sampling presents. Specifically, we sample configurations at times growing in powers of $`2`$ up a maximum time $`t_{\text{max}}`$. We begin a new sampling cycle each time $`t_{\text{max}}`$ (relative to the cycle starting time) is reached. This sampling method allows for calculation of dynamic properties on time scales spanning eight orders of magnitude (from 1 fs to 100 ns) using a relatively small amount of disk space. Still, more then 2 GB of storage was required for storing configurations at $`T=210K`$. Our simulations have a speed of approximately 200 $`\mu `$s per update per molecule on a MIPS R10000 processor, representing a total calculation of approximately 8.4 years of CPU time, including the systems of 1728 molecules discussed in the appendix. For the larger systems, we utilize a parallelized version of our simulation code on eight processors to improve performance.

## IV Static Structure Factor

We first summarize the structural properties of our simulations in order to better understand the relationship between the changes in structure with the changes in dynamic behavior, which we will detail in Sec. VIII. Other studies have considered the structural and thermodynamic properties of SPC/E in a large region of the $`(P.T)`$ plane , so the present discussion is brief.

The MCT theory requires as input the static density-density correlation functions. In the case of water, the structure of the system is very sensitive to the value of the external control parameter $`(P,T)`$. Hence, for all state points simulated, we calculate the oxygen-oxygen partial structure factor 

$$S(q)\frac{1}{N}\left|\underset{j=1}{\overset{N}{}}e^{i𝐪𝐫_j}\right|^2.$$
(10)   
Several studies have carefully calculated the structure of simulated water, and found surprisingly good agreement with experiments . At $`T=210`$ K, we show the structural changes from low to high density \[Fig. 3(a)\]. The structure at low density/pressure is similar to that observed for low-density amorphous (LDA) solid water, consisting of an open tetrahedral network. At high density/pressure, water is very similar to high density amorphous (HDA) solid water, where core-repulsion dominates, similar to simple liquids under pressure.

We show the evolution of $`S(q)`$ as a function of $`T`$ along the $`\rho =1.0`$ g/cm<sup>3</sup> isochore in Fig. 3(b). We note that in the temperature range from 190 K to 300 K, where the dynamics show the most dramatic change in behavior, $`S(q)`$ shows only small changes in the first two peaks. Also, the location of the first maximum $`q_0`$ in $`S(q)`$, the wave vector at which $`F(q,t)`$ typically shows the slowest relaxation, does not appear to change significantly. All other densities and temperatures show a relative smooth interpolation of Figs. 3(a) and (b).

## V Mean-Squared Displacement and Diffusion

The mean-squared displacement $`r^2(t)|𝐫(t)𝐫(0)|^2`$ is shown in Fig. 4. All the curves show $`t^2`$ dependence at small time, as expected in the “ballistic” regime. For low T (e.g. $`T=210`$ K \[Fig. 4(a)\]), $`r^2(t)`$ shows relatively flat behavior over 3-4 decades in time. This is the “cage” region, in which a molecule is trapped by its neighbors and cannot diffuse, and is only vibrating within its “cage.” At low $`P`$, the cage consists of hydrogen-bonded neighbors in a tetrahedral configuration. This cage is relatively strong, compared to simple liquids, because of the H bonds. The size of the cage may be estimated by the value $`r^2(t)`$ at the plateau, as shown in the inset of \[Fig. 4(a)\]. Surprisingly, the size of the cage is not monotonic with density, and has a maximum at $`\rho 1.1`$ g/cm<sup>3</sup>. We shall see that this corresponds roughly to the $`\rho `$ at which $`D`$ also has a maximum. We observe a small bump in $`r^2(t)`$ at $`t0.35`$ ps, as observed in Ref. . A system size study indicates that this may be attributed to finite size effects .

For long times, all the curves show linear $`t`$ dependence, indicating that our simulations are in the diffusive regime. We extract the diffusion constant $`D`$ using the asymptotic relation $`r^2(t)=6Dt`$. We plot the density dependence of $`D`$ in Fig. 5 and find that the SPC/E potential, like water, shows an anomalous increase in $`D`$ on increasing density. We also point out the feature that $`D`$ shows a slight increase at very low density; namely, at $`\rho =0.90`$ g/cm<sup>3</sup> and $`T=210`$ K. This can be attributed to the fact that the liquid is extremely stretched at this density, causing an increase in the defects of the bond network, and thus increased diffusivity. We show $`D`$ as a function of pressure along several isotherms to compare with experimental measurements \[Fig. 6. The anomalous increase in $`D`$ is qualitatively reproduced by our calculations for the SPC/E model, but the quantitative increase of $`D`$ is significantly larger than that

observed experimentally. This discrepancy may arise from the fact that the SPC/E potential is under-structured relative to water , so applying pressure allows for more bond breaking and thus greater diffusivity than observed experimentally. We also find that the pressure where $`D`$ begins to decrease with pressure — normal behavior for a liquid — is larger than that observed experimentally . This comparison of $`D`$ with experiment leads us to expect that while the qualitative dynamic features we observe in the SPC/E potential may aid in the understanding of the dynamics of water under pressure, they will likely not be quantitatively accurate.

We estimate $`D`$ along the isobars $`P=80`$ MPa, 0 MPa, 100 MPa, 200 MPa, 300 MPa, and 400 MPa from the isochoric data. We confirm that along the -80 MPa isobar, our estimates agree with the $`80`$ MPa calculations of Ref. , which employs the same truncation of the potential used here (see Sec. III). Along the 0 MPa

isobar, our estimates of $`D`$ are smaller than those calculated for SPC/E in Ref. , perhaps because Ref.  chooses a different truncation of the electrostatic terms—highlighting the extreme sensitivity of the dynamics to changes in the potential.

We fit $`D`$ by the power law of Eq. (1) along both isochores and the estimated isobars for $`T300`$ K \[Figs. 7 and 8\]. The values of the two fit parameters $`T_c`$ and $`\gamma `$ are given in Table II . We also include the data from Ref.  along the $`\rho =1.0`$ g/cm<sup>3</sup> isochore and from Ref.  along the $`P=80`$ MPa isobar to improve the quality of the fits. At $`\rho =1.40`$ g/cm<sup>3</sup>, we exclude $`T=210`$ K when fitting $`D`$ and obtain $`T_c=209.3`$, since we expect the power law of Eq. (1) to fail for $`TT_c+5`$ K, because activated processes – such as “hopping” not accounted for in the idealized MCT – become significant and aid diffusion. To demonstrate the presence of hopping at $`\rho =1.40`$ and $`T=210`$ K, we plot the “self” part of the van Hove correlation function $`G_s(r,t)`$, which measures the distribution of particle displacements $`r`$ at time $`t`$, for several densities at $`T=210`$ K \[Fig. 9\]. For these densities where a power law adequately describes $`D`$, there is a single peak. At $`\rho =1.40`$ g/cm<sup>3</sup>, we see a “shoulder” in $`G_s(r,t)`$ at $`r0.2`$ nm in addition to a well-defined peak at $`r0.05`$ nm, indicating that particle hopping is significant.

Along the $`\rho =1.00`$ g/cm<sup>3</sup> we have simulated to significantly lower $`T`$, allowing us to study the temperature dependence of $`D`$ for $`TT_c`$. Fig. 10 shows that the lowest temperatures are consistent with the Arrhenius form

$$D=D_{\mathrm{}}\mathrm{exp}(E/k_BT).$$
(11)   
Arrhenius temperature dependence of $`D`$ is not surprising in this region, since for $`T<T_c`$ it is expected that the energy barriers the system must overcome to

rearrange exceed the thermal energy . Hence the motion of the system is dominated by activated jumps over the energy barriers, as described by Goldstein . We obtain an activation energy of $`E65`$ kJ/mol and extrapolate a glass transition temperature $`T_g125`$ K , surprisingly close to the experimental value of 136 K. The extrapolated value of $`T_g`$ is similar to that estimated in Ref.  which studied hydrogen bond dynamics. Moreover, our results are consistent with a crossover from “fragile” behavior (the behavior described by MCT) for $`TT_c`$, to “strong” behavior (Arrhenius behavior with

$`Ek_BT_g/2540`$ kJ/mol for our estimate of $`T_g`$) for $`TT_c`$. The possibility of a “fragile-to-strong” crossover in water has been discussed recently based on experimental findings , but lower temperatures are required to test this possibility in the SPC/E model.

## VI Isochrones of $`D`$ and the Locus of $`T_c(P)`$

To construct isochrones of $`D`$ (lines of constant $`D`$), we first estimate $`T(D)`$ at values of $`D=10^5`$cm<sup>2</sup>/s, $`10^{5.5}`$ cm<sup>2</sup>/s, $`10^6`$cm<sup>2</sup>/s, and $`10^7`$cm<sup>2</sup>/s, using the fits of Figs. 7 and 8. Along the isobaric paths, we know already $`P`$ for these points, and along isochores we may estimate the value of $`P`$ using the results presented in Table I. We plot the isochrones in Fig. 11.

We also show the the loci of $`T_c(P)`$ in Fig. 11(a), obtained from the fits in the previous section. We know $`P`$ at $`T_c`$ along the isobaric paths, and we estimate the $`P`$ at $`T_c`$ along isochores by extrapolating $`P`$ in Table I to $`T_c`$.

Using the experimental diffusion data of Ref. , we also construct the behavior of the experimental isochrones following the same technique \[Fig. 11(b)\]. The shape of the locus of $`T_c(P)`$ compares well with that observed experimentally , and changes slope at roughly the same pressure \[Fig. 11\]. Therefore, an explanation of the SPC/E dynamics using the MCT would support using the MCT framework as an interpretation of the experimentally found locus of $`T_c(P)`$. We find, however, that $`\gamma `$ decreases with $`P`$ for the SPC/E model, while $`\gamma `$ increases with $`P`$ \[Fig. 12\]. This disagreement underscores the need to improve the dynamic properties of water models, most of which already provide an adequate account of static properties .

## VII Intermediate Scattering Function

We plot the intermediate scattering function $`F(q_0,t)`$ in Fig. 13(a) for all $`T`$ along the $`\rho =1.00`$ g/cm<sup>3</sup> isochore, where $`q_0=18.55`$ nm<sup>-1</sup>, the approximate value of the first peak of $`S(q)`$ where the relaxation of $`F(q,t)`$ is slowest. We define the relaxation time $`\tau `$ by $`F(t=\tau )=e^1`$. We show $`\tau `$ along isotherms in Fig. 5(b), from which it is obvious that $`\tau `$ has very similar behavior to $`D^1`$. Indeed, MCT predicts that the product $`D\tau `$ is constant along isochores, which we test in Fig. 5(c). We find that $`D\tau `$ increases slightly on cooling, but remains relatively constant along each isochore. The weak residual

$`T`$-dependence in $`D\tau `$ should be subjected to a deeper scrutiny to find out if it is related to a $`q`$vector dependent correction to scaling (since $`D`$ is a $`q=0`$ quantity) or to the progressive breakdown of the validity of the ideal MCT on approaching $`T_c`$.

The study of the time dependence of $`F(q_0,t)`$ allows us to test the predicted relation between the exponents $`b`$ and $`\gamma `$ (see Eqs. (6) and (7)). Since the value of $`b`$ is completely determined by the value of $`\gamma `$ , calculation of these exponents for SPC/E determines if MCT is consistent with our results. The range of validity of the van Schweidler power law \[Eq. (4)\] is strongly $`q`$-dependent , making unambiguous calculation of $`b`$ difficult.

Fortunately, according to MCT , at large $`q`$-vectors,

the stretching exponent $`\beta (q)`$, which characterizes the the long-time behavior of $`F(q,t)`$ (see Eq. (8)), is controlled by the same exponent $`b`$ at large $`q`$. Fits of $`F(q_0,t)`$ according to Eq. (8) are shown for many $`q`$ values at $`T=210`$ K and $`\rho =1.00`$ g/cm<sup>3</sup>. The same fit quality is observed for all other low $`T`$ state points. The $`q`$-dependence of $`\beta (q)`$ for $`\rho 1.30`$ g/cm<sup>3</sup> and $`T=210`$ K is shown in Fig. 14  for $`F(q,t)`$. In addition, we show the expected value of $`b`$ according to MCT, using the values of $`\gamma `$ extrapolated from Fig. 12. The large-$`q`$ limit of $`\beta `$ appears to approach the value predicted by MCT. Hence we conclude that the dynamic behavior of the SPC/E potential in the pressure range we study is consistent with slowing down as described by MCT \[Fig. 2\]. We also checked that the values of $`b`$ calculated from Eq. (9) are consistent with the von Schweidler power law Eq. (4), but that corrections to scaling in $`t^{2b}`$ are relevant at several $`q`$ vectors, as discussed in Ref. .

## VIII Relationship of Structure to Dynamics

The results shown in Fig. 14, and the observed power-law dependence of diffusivity, suggest that MCT is able to predict the dynamical behavior of SPC/E water in a wide range of $`P`$ and $`T`$. As discussed above, the structure of the liquid changes significantly under increased pressure. To highlight the effect of structural changes on dynamic properties, we consider an approximately isochronic path – along which $`D`$ remains nearly constant – such that the changes in dynamic properties we observe on increasing $`P`$ are confined to their $`q`$-vector dependence. We select 5 state points with $`D=(0.30\pm 0.09)\times 10^{}6`$ cm<sup>2</sup>/s: (i) $`T=220`$ K, $`\rho =1.00`$ g/cm<sup>3</sup>, (ii) $`T=210`$ K, $`\rho =1.05`$ g/cm<sup>3</sup>, (iii) $`T=210`$ K, $`\rho =1.10`$ g/cm<sup>3</sup>, (i) $`T=210`$ K, $`\rho =1.20`$ g/cm<sup>3</sup>, (iv) $`T=220`$ K, $`\rho =1.30`$ g/cm<sup>3</sup>, and (v) $`T=240`$ K, $`\rho =1.40`$ g/cm<sup>3</sup>. We show in Fig. 15 the $`q`$-dependence of the $`\alpha `$-relaxation time $`\tau (q)`$ extracted from the fit of $`F(q,t)`$ to the stretched exponential of Eq. (8). For all state points, the $`q`$-dependence of $`\tau `$ follows the $`q`$-dependence of $`S(q)`$, as commonly observed in supercooled liquids and in solutions of the full $`q`$-vector dependent mode coupling equations. We also note that $`\tau (q)`$ is well described by the relation

$$\tau (q)S(q)/q^2$$
(12)   
(the de Gennes narrowing relation), as shown in the same figure. The MCT prediction for the $`q`$-dependence of $`\tau `$ is often very close to the relation (12).

## IX Discussion

We have presented extensive simulations that provide evidence for interpreting the dynamics of the SPC/E potential in the framework of MCT. Our calculations also provide a necessary test of the relation predicted between the diffusivity exponent $`\gamma `$ and the von Schweidler exponent $`b`$ for a wide range of values $`\gamma `$ and $`b`$. Our results support interpretation of the experimental locus $`T_c(P)`$ as the locus of MCT transitions.

We found that on increasing pressure, the values of the exponents become closer to those for hard-sphere ($`\gamma =2.58`$ and $`b=0.545`$) and Lennard-Jones ($`\gamma =2.37`$ and $`b=0.617`$) systems , thereby confirming that the hydrogen-bond network is destroyed under pressure and that the water dynamics become closer to that of normal liquids, where core repulsion dominates. A significant result of our analysis is the demonstration that MCT is able to rationalize the dynamic behavior of the SPC/E model of water at all pressures. In doing so, MCT encompasses both the behavior at low pressures, where the mobility is essentially controlled by the presence of strong energetic cages of hydrogen bonds, and at high pressures, where the dynamics are dominated by excluded volume effects. We also showed how these structural changes are reflected in the $`q`$-dependence of $`F(q,t)`$.

Our results underscore the need to improve the dynamic properties of potentials for realistic simulations of

water and other other materials. Of the many potentials available for studying water, only the SPC/E potential is known to display the power-law dependence of dynamic properties, but even SPC/E fails to reproduce the power law quantitatively. A recent study of the ST2 potential  found that the $`T`$-dependence of $`D`$ is consistent with an Arrhenius $`T`$-dependence for $`T300`$ K, crossing over to a another region of Arrhenius behavior for $`T275`$ K , in contrast to the non-Arrhenius behavior observed in real water and to our interpretation based on MCT for $`TT_c`$. The presence of a low-$`T`$ Arrhenius regime in the ST2 potential might be due to activated processes, that are expected to dominate the dynamics of fragile liquids below $`T_c`$, as we observed for the SPC/E potential. Hence the ST2 potential may provide an excellent opportunity to study these activated processes on a smaller time scale than is typically observed for most fragile liquids.

Finally, we stress that a full comparison between theory and simulation data requires a complete solution of the recently proposed molecular-MCT (the extension of MCT to molecules of arbitrary shape) . A detailed solution of the complicated molecular-MCT equations in such large region of $`T`$ and $`P`$ values would requires computational effort beyond the present possibilities, but a detailed comparison between molecular-MCT and MD data for one selected isobar is underway .

## X Acknowledgments

We thank C.A. Angell, A. Geiger, E. La Nave, A. Rinaldi, S. Sastry, A. Scala and R.J. Speedy for enlightening discussions and comments on the manuscript. We especially thank S. Harrington for his contributions to the early stages of this work. We thank the Boston University Center for Computational Science for access to the 192 processor SGI/Cray Origin supercomputer. F.S. is supported in part by MURST (PRIN 98). The Center for Polymer Studies is supported by the NSF Grant No. CH9728854.

## A Finite-Size Effects

Some recent work  indicates that significant finite-size effects can affect results at temperatures close to the MCT $`T_c`$. Most of our systems are farther than 5% from $`T_c`$ (i.e., $`T/T_c10.05`$) so that the relatively small size of our system should not affect our results. We simulated two independent systems of 1728 molecules at both $`T=200`$ K and $`190K`$ and $`\rho =1.00`$ g/cm<sup>3</sup> to check if significant finite size effects appear at low $`T`$. The results are shown in Table. III. We observe no significant deviations from the system of 216 molecules at $`T=200`$ K. Hence we believe that no strong finite-size effects are present in the 216 molecule system for $`T200`$ K. However at $`T=190`$ K, the potential energy of the larger system appears to be significantly smaller than that in the smaller system. We lack adequate computer resources to make a reliable estimate the diffusivity in the larger system, but simulations are continuing in order to check the possible finite-size effects at this temperature.