Waterlike Structural and Excess Entropy Anomalies in Liquid Beryllium

13294

J. Phys. Chem. B 2007, 111, 13294-13300

Waterlike Structural and Excess Entropy Anomalies in Liquid Beryllium Fluoride Manish Agarwal and Charusita Chakravarty* Department of Chemistry, Indian Institute of Technology-Delhi, New Delhi 110016, India ReceiVed: July 9, 2007; In Final Form: September 7, 2007

The relationship between structural order metrics and the excess entropy is studied using the transferable rigid ion model (TRIM) of beryllium fluoride melt, which is known to display waterlike thermodynamic anomalies. The order map for liquid BeF2, plotted between translational and tetrahedral order metrics, shows a structurally anomalous regime, similar to that seen in water and silica melt, corresponding to a band of state points for which average tetrahedral (qtet) and translational (τ) order are strongly correlated. The tetrahedral order parameter distributions further substantiate the analogous structural properties of BeF2, SiO2, and H2O. A region of excess entropy anomaly can be defined within which the pair correlation contribution to the excess entropy (S2) shows an anomalous rise with isothermal compression. Within this region of anomalous entropy behavior, qtet and S2 display a strong negative correlation, indicating the connection between the thermodynamic and the structural anomalies. The existence of this region of excess entropy anomaly must play an important role in determining the existence of diffusional and mobility anomalies, given the excess entropy scaling of transport properties observed in many liquids.

1. Introduction Water displays a number of thermodynamic and kinetic anomalies when compared with simple liquids, particularly at low temperatures and pressures. The density anomaly is the best known of these unusual properties of water and corresponds to a density-temperature (FT) regime for which the isothermal expansion coefficient R is negative. The existence of such a density anomaly is associated with the existence of anomalies in thermodynamic response functions, such as the isobaric heat capacity (Cp) and the isothermal compressibility (κT). The kinetic anomalies of water are associated with an increase in various measures of molecular mobility, such as diffusivity and viscosity, with increasing density. The presence of these liquid-state anomalies is believed to be due to anisotropic hydrogen-bonding intermolecular interactions which lead to the formation of a three-dimensional random fluctuating network with local tetrahedral order.1-3 Such anomalous behavior is, however, not unique to water but has been identified in a diverse range of liquids. A number of inorganic melts possess such networks and are structural analogues of water despite the difference in the underlying interatomic interactions and, consequently, in the melting temperatures. Examples of such network-forming liquids are molten SiO2, BeF2, silicon, germanium, bismuth, and GeO2.4-12 Both thermodynamic and kinetic anomalies have been demonstrated, using either experiment or simulations, in the case of SiO2 and BeF2. More recently, waterlike anomalies have been demonstrated in model liquids with isotropic core-softened pair interactions, such as the Gaussian-core and one- and two-scale ramp liquids.13-20 Experimental realizations of such coresoftened fluids are provided by certain categories of colloidal and star polymer solutions.21 A necessary condition for the existence of waterlike anomalies is the existence of an entropy anomaly, corresponding to a temperature regime for which the entropy of the liquid increases * Corresponding author. Phone: (+) 91 11 2659 1510. Fax: (+) 91 11 2686 2122. E-mail: [email protected].

on isothermal compression, instead of decreasing as in the case of simple liquids.1,2 The region of the entropy anomaly must be bounded by points at which the total entropy S reaches an extrema with respect to the density F such that (∂S/∂F)T ) (-V/ F)R/κT ) 0. Since R must be negative in the region of the entropy anomaly, an entropy anomaly implies the existence of a waterlike density anomaly. Along an isobar, R must be zero at the temperature of maximum density (TMD). On thermodynamic grounds, it can be shown that a negatively sloped locus of TMD points on the TP plane should lead to anomalous behavior of the compressibility and isobaric heat capacity.22,23 Thus, the presence of the entropy anomaly will typically signal the presence of a set of connected thermodynamic anomalies. The total entropy, S, of a classical liquid can be partitioned into ideal, Sid, and excess, Se, contributions. At a given temperature, Sid must decrease with compression. Therefore, any entropy anomaly must originate from the behavior of the positional or excess contributions. Rosenfeld-type scaling relationships of the form

X ) A exp(-RSe)

(1)

where X is a suitably scaled transport property have been demonstrated in a wide range of liquids where the coefficients R and A depend on the nature of the interparticle interactions as well as the type of transport property, for example, diffusivity, viscosity, or thermal conductivity.24-28 Since such scaling relationships have been demonstrated for core-softened fluids as well as binary ionic melts, the existence of an entropy anomaly implies the likelihood of diffusional and related kinetic anomalies, in addition to the thermodynamic anomalies.29-32 It must be kept in mind, however, that the deviations from Rosenfeld scaling will occur in the supercooled regime, especially in the neighborhood of the glass transition, and this must be taken into account when considering the relationship between mobility and entropy anomalies.

10.1021/jp0753272 CCC: $37.00 © 2007 American Chemical Society Published on Web 10/27/2007

Anomalies in Liquid Beryllium Fluoride

J. Phys. Chem. B, Vol. 111, No. 46, 2007 13295

It is useful to compare the microscopic behavior underlying the excess entropy anomaly in the case of core-softened and network-forming liquids.29-37 In the case of core-softened liquids, the length scales associated with pair correlations in the high-density and low-density limits are different while in the intermediate density regime; there are two competing length scales. As the system is compressed along an isotherm starting from a low density, the additional length scale comes into play, and the resulting increase in the number of possible packing arrangements results in an increase in entropy. In the case of tetrahedral, network-forming liquids, there is a competition between two different local order metrics, rather than between two length scales. The preferred local symmetry at low densities is tetrahedral but shifts to being icosahedral at high densities when the isotropic, short-range repulsions become dominant. Because of this density-driven change in local symmetry, the local structural order reaches a minimum, and the excess entropy reaches a corresponding maximum at some intermediate density along an isotherm. While this change in local orientational order will be reflected in the pair correlation function and therefore can be regarded as a change in effective length scales, the extent to which an effective pair interaction can completely represent a range of thermodynamic properties for a complex fluid is subject to some debate.38,39 The diffusional mechanisms are likely to be qualitatively different in the two categories of fluids. In a core-softened fluid, diffusion is likely to proceed by a combination of collisions and cage relaxation, similar to those in simple liquids while in a network-forming liquid, at low temperatures, an additional diffusional mechanism will exist because of coupling of local vibrational modes with network reorganizations. The importance of the interplay between two different order metrics in determining the anomalous properties of water was first demonstrated by Errington and Debenedetti by introducing order metrics to measure the degree of orientational as well as translational correlations.3 The orientational order can be gauged using a local tetrahedral order parameter, qtet, associated with an atom i and is defined as

qtet ) 1 -

3

3

4

∑ ∑ (cos ψjk + 1/3)2

8 j)1 k)j+1

(2)

where ψjk is the angle between the bond vectors rij and rik where j and k label the four nearest neighbor atoms of the same type.3 The translational order parameter, τ, measures the extent of pair correlations present between tetrahedral centers (e.g., O in H2O and Be in BeF2) and is defined as

τ ) (1/ξc)

∫0ξ |g(ξ) - 1| dξ c

(3)

where ξ ) rF1/3, r is the pair separation, g(ξ) is the pair correlation function and ξc is a suitably chosen cutoff distance.40 Since τ increases as the random close-packing limit is approached, it can be regarded as measuring the degree of density ordering. At a given temperature, qtet will show a maximum and τ will show a minimum as a function of density; the loci of these extrema in the order define a structurally anomalous region in the density-temperature (FT) plane within which qtet and τ are strongly correlated. The region of the density anomaly, where (∂F/∂T)P > 0, is bound by the structurally anomalous region. The diffusionally anomalous region ((∂D/∂F)T > 0) closely follows the boundaries of the structurally anomalous region, specially at low temperatures. In water, the structurally anomalous region encloses the region of anomalous diffusivity.

TABLE 1: TRIM Potential Parameters for BeF242-44 σ+ (Å)

σ- (Å)

s (Å)

z+

z-

n+

n-

b (kJ mol-1)

0.93

1.24

0.29

2

-1

2

8

34.33

Subsequently, a similar nesting of regions corresponding to structural, diffusional, and density anomalies was shown in the case of silica, except that the order of the structural and diffusional anomalies was reversed.7 In this paper, we study the relationship between the local order metrics and the existence of waterlike thermodynamic anomalies, using beryllium fluoride melt as a convenient and interesting example of a waterlike, network-forming liquid. On the basis of stoichiometry and radius ratio rules, BeF2 is known to be a close structural analogue of SiO2.41 Molecular dynamics simulations of BeF2 using a rigid-ion model potential suggest the presence of density and mobility anomalies similar to those seen in silica and water.42-45 Interest in this system from the point of view of potential applications in nuclear technology has stimulated some recent experimental and computational studies, including the development of accurate, many-body potentials parametrized using data from ab initio electronic structure simulations.46-48 The behavior of local order metrics in this system is of special interest since it has been argued that the constraint of local tetrahedrality should increase in importance in the following sequence: SiO2 < BeF2 < H2O < Si, resulting in a more pronounced anomalous behavior along the series.11 The present study focuses on the structural anomaly in molten BeF2 and complements another recent work in which we have studied the thermodynamic properties of liquid beryllium difluoride (BeF2).51 Since our focus is to understand the relationship between structural order, entropy, and waterlike anomalies, it is necessary to cover a range of number densities and temperatures comparable to those previously used for SPC/E water and BKS silica. For an exploratory study of this type, the well-studied transferable rigid ion model (TRIM) potential for BeF242-44 was felt to be more appropriate than a computationally expensive many-body parametric potential. The TRIM potential is an adaptation of the Born-Mayer-Huggins or Tosi-Fumi potential and has been shown to be reasonable for a number of molten salts.49,50 In the specific case of BeF2 where polarization effects are small, recent comparisons with first principles results suggest that the TRIM potential overestimates pressures but is in reasonable agreement with the equation of state with regard to behavior along the 1.96 g cm-3 isochore which corresponds to the experimental density at 1 atm of pressure.47 We summarize briefly the results from our study of the thermodynamic properties of the TRIM model of BeF2 here51 (Table 1) to place the present study on structural anomalies in context. We have studied the density range from 1.5 to 3.0 g cm-3 and a temperature range from 1250 to 3000 K in our simulations. Hemmati et al. studied a temperature range from 500 to 5000 K and a density range from 1.5 to 2.3 g cm-3. Problems with equilibration occur at temperatures of 1250 K and below, especially at low densities. In the present study on structural order, we focus therefore on temperatures of 1500 K and above. Our density range has been extended, in comparison to the previous studies, to study the effects of compresssion. The pressures along the highest density isochore lie between 10 and 25 GPa for which it is, in principle, accessible in highpressure experimental studies.52,53 Using the equation of state information from the simulations, we mapped out the locus of the temperature of maximum density (TMD) in the temperature-pressure and density-temperature planes. By using a pair

13296 J. Phys. Chem. B, Vol. 111, No. 46, 2007

Agarwal and Chakravarty

correlation approximation to the excess entropy, an anomalous rise in excess entropy on isothermal compression was demonstrated. The anomalous behavior of the excess entropy was found to be largely connected with the behavior of the Be-F pair correlation function. The internal energy was shown to have T3/5 temperature dependence, similar to that seen in SiO2. The pair correlation entropy showed a T-2/5 temperature dependence only at high densities and temperatures. The correlation plots between internal energy and the pair correlation entropy for isothermal compression showed the characteristic features expected of network-forming liquids with waterlike anomalies.29 The tagged particle potential energy distributions were shown to have a multimodal form at low temperatures and densities similar to those seen in other liquids with three-dimensional tetrahedral networks, such as water and silica.33,34,37 This paper is organized as follows. The details of the molecular dynamics simulations are identical to those given in ref 51, but for completeness, we summarize them in section II. Sections III and IV contain the results and conclusions, respectively. 2. Computational Details The transferable rigid ion model (TRIM) potential is used to model interatomic interactions since it has been extensively used to study liquid state anomalies and the glass transition in ionic melts.42-45 The TRIM potential is pair-additive and models electrostatic long-range interactions with Coulomb interactions and short-range repulsions using an exponential repulsion term. The pair interaction between atoms i and j is given by:

φTRIM(rij) )

(

) (

)

zizje2 z i zj σi + σj - rij + 1+ + b exp 4πorij ni nj s

(4) where rij is the distance between atoms i and j. The parameters associated with an atom of type l are the charge zl, the number of valence-shell electrons nl, and the ionic size parameter σl. The repulsion parameter b and the softness parameter s are assumed to be the same for all three types of pair interactions. The parameters for the TRIM potential used in this work are identical to those used in our previous study. Molecular dynamics simulations of a system of 150 Be and 300 F ions were carried out in the canonical (NVT) ensemble, using the DL_POLY software package.54,55 We used cubic periodic boundary conditions in conjunction with Ewald summation to account for long-range contributions from electrostatic interactions.56,57 A Berendsen thermostat, with the time constant τB ) 200 ps, was used to maintain the desired temperature for the production run. The Verlet leapfrog algorithm with a time step of 1 fs was used to integrate the equations of motion. The system was simulated at 8 temperatures in the 1250-3000 K range along 8-10 isochores in the density range from 1.5 to 3.0 g cm-3. After the equilibration procedure described in ref 51 was completed, fairly long production runs of 6-8 ns were performed with a goal of obtaining reliable structural data. The thermodynamic excess entropy was estimated using the pair correlation approximation.29,32,58,59 For a binary system consisting of a total number of N particles enclosed in a volume V at temperature T, the pair correlation entropy will be given by60

S2/NkB ) -2πF

xRxβ ∫0 ∑ R,β

∞

{gRβ(r) ln gRβ(r) [gRβ(r) - 1]}r2 dr (5)

where gRβ(r) is the pair correlation function between spherically symmetric particles of type R and type β and xR is the mole fraction of particles of type R. Considering an MX2 ionic melt as a binary mixture of M and X ions interacting via isotropic pair interactions, we can apply the above formula to obtain the pair correlation entropy, S2. The reliability of this approximation for binary ionic melts has been discussed previously.29,51 We are currently further testing the accuracy of this approximation for ionic melts by comparing the pair correlation entropy with the thermodynamic excess entropy for BKS silica and TRIM BeF2. 3. Results 3.1. Structural Anomaly. In order to identify a structurally anomalous regime in a tetrahedral, network-forming liquid, it is necessary to consider the relationship between the average local tetrahedral order, as defined in eq 2, with the translational order, as defined in eq 3. Figure 1 shows qtet as a function of the density, F. At low temperatures, the variation of local tetrahedral order with density is marked, and there is a maximum at approximately 1.6 g cm-3. The locus of these maxima in qtet forms the low-density boundary of the region of structural anomaly. As temperature increases and the network connectivity is attenuated, the variation in qtet with density is reduced, and the maxima become less prominent and shift slightly to higher densities. At 3000 K, the maximum in qtet occurs at approximately 1.8 g cm-3. Figure 2 shows the behavior of the translational order parameter, τ, as a function of density. Normal liquid behavior, characterized by a strong positive correlation between τ and F is seen at all temperatures for densities greater than 2.4 g cm-3. At low temperatures and densities, τ decreases with F, and a clear minimum in translational order can be observed between 2.2 and 2.4 g cm-3 for isotherms below 2250 K. For isotherms at 2500 K and above, τ rises only very slowly with density until approximately 2.2 g cm-3, after which it shows the characteristic rise with F as increasing density enhances pair correlations. The locus of minima in the τ(F) curves forms the high-density boundary of the structurally anomalous region. Therefore, for temperatures of 2250 K or less, the existence of a structurally anomalous region, bounded by maxima on qtet(F) and minima in τ, is well-defined. Figure 3 shows that the order map plotted between translational and orientational order for liquid BeF2 is very similar to that seen for water and liquid silica. The structurally anomalous region where qtet and τ are strongly correlated can be identified. The existence of a “forbidden” region in the order map is obvious in that, for any value of τ, there is a maximum observable value of qtet. The order map of BeF2 resembles that of silica in that state points in the anomalous region lie in a band in the qtet-τ plane rather than along a line in the case of water. A more careful consideration of the curvature around the minimum in the τ versus qtet plots suggests that the correlation between qtet and τ is stronger in BeF2 than in SiO2 but is not as strong as in the case of SPC/E water where there is a rigid molecular entity. 3.2. Tetrahedral Order Parameter Distributions. We examine the tetrahedral order parameter distributions for BeF2 to further establish the analogy with waterlike anomalies in silica. Figure 4a shows the P(qtet) distribution for the 1.6 g cm-3

Anomalies in Liquid Beryllium Fluoride

J. Phys. Chem. B, Vol. 111, No. 46, 2007 13297

Figure 1. Variation of tetrahedral order parameter (qtet) for Be2+ BeF2 melt with density (F) along several isotherms.

Figure 2. Variation of translational order parameter (τ) for Be2+ in BeF2 with density (F) along several isotherms.

Figure 3. Order map showing correlation between translational (τ) and tetrahedral (qtet) order.

isochore. At relatively low temperatures of around 1500 K, the qtet distribution has a peak at approximately 0.9 indicating that there is a high degree of local tetrahedral order. As temperature increases, the peak shifts to lower values of qtet and becomes bimodal at between 2000 and 2250 K. Further increase in temperature results, enhancing the second peak at qtet ≈ 0.5. Beyond 2750 K, the distributions are unimodal with a peak at ≈0.5. Thus, the transition in local symmetry preference from tetrahedral to simple liquid behavior is obvious. The region of the broad bimodal distribution corresponds to the disappearance of the minima in the pair correlation entropy. Figure 4b shows the P(qtet) distribution at different densities along the 2000 K isotherm. The transition from a skewed unimodal distribution with high local tetrahedral order at low densities to a bimodal distribution at intermediate densities in the anomalous region to a broad, nearly symmetric distribution at high densities is evident. To further characterize the P(qtet) distributions, we also examined the density dependence of the first four semi-invariant moments of the tetrahedral order distribution along different isotherms. The third (I3) and fourth (I4) moments measure of deviation from Gaussian behavior in the tetrahedral order

Figure 4. Tetrahedral order parameter distributions, P(qtet), along (a) at different temperatures along 1.6 g cm-3 isochore and (b) different densities along the 2000 K isotherm. Shift in peak positions with temperature (part a) or density (part b) are indicated by arrows. Distributions at 2250 K and 2500 K in (a) and 2.0 g cm-3 in (b) are shown with lines for clarity.

distribution. Along a given isotherm, I3 approaches zero as density increases while I4 approaches a nearly constant nonzero value. 3.3. Relationship between Structural and Entropy Anomalies. Figure 5a shows the relationship between qtet and S2. At temperatures of 2000 K or less, qtet shows a negative correlation with S2 for densities below 2.5 g cm-3. In this regime, decreasing tetrahedral order implies an increase in positional disorder and therefore an increase in S2. At higher densities, qtet is decorrelated from S2. Isotherms at higher temperatures show a region where qtet is positively correlated with S2, implying that an increase in qtet implies an increase in S2. These are, however, state points for which the magnitude of qtet is small, implying that the tetrahedral network is essentially destroyed; moreover, qtet decreases with increasing density. Since qtet reflects fourcoordinate, tetrahedral order, a decrease in qtet with increasing density for these isotherms probably reflects an increase in coordination number to values well above 4. Such an increase in coordination number will be consistent with a correlated decrease in S2 and qtet in the regime of simple liquid behavior. Figure 5b shows the correlation between τ and S2. In the anomalous regime, τ shows a strong negative correlation with S2 while, in the normal regime, τ shows a very weak negative correlation with S2. SiO2 shows a similar change in slope of the τ versus S2 curves at the boundary between the normal and the anomalous regimes. The negative correlation between S2 and τ is obvious from the definition of the two quantities in terms of g(r) in the case of monatomic fluids. In the case of SiO2 and H2O, τ has been typically defined in the literature in terms of the pair correlation function connecting tetrahedral centers. In the present case, this corresponds to using the gBe-Be(r) pair correlation function when defining τ. In contrast, S2 contains contributions from the three possible pair correlation functions in a binary mixture. Clearly, in the anomalous regime, the Be-Be pair correlation function is very sensitive to local

13298 J. Phys. Chem. B, Vol. 111, No. 46, 2007

Agarwal and Chakravarty

Figure 7. Regions in the F-T plane corresponding to anomalies in total entropy (S), excess entropy (S2), and density (F) anomaly. The boundaries of the region of density anomaly are defined by the locus of temperatures of maximum density (TMD), as derived from the equation of state. The locus of extrema in S ≈ Sid + S2 are nearly co-incident with the TMD locus.

Figure 5. Correlation of with pair correlation entropy, S2, with (a) tetrahedral order metric, qtet, and (b) translational order metric, τ. As a guide to the eye, the highest and lowest densities along some isotherms have been marked by vertical and horizontal arrows, respectively. The density (in grams per cubic centimeter) at which the τ(S2) curve shows a sharp change in slope has also been marked on the 1500, 1750, and 2000 K isotherms.

to somewhat higher densities than the structural anomaly. This is a consequence of the fact that τ is defined only in terms of the gBe-Be(r) while S2 contains contributions from the three possible pair correlation functions. Given the difference in behavior of S2 and τ, it is important to consider which of the two quantities provides a physically more reasonable measure of local order because of pair correlations. In the case of a rigid molecular fluid like water, where the tetrahedral center and the center-of-mass essentially coincide, defining τ in terms of the pair correlation between tetrahedral centers of mass is appropriate. In the case of ionic melts, such as SiO2 and BeF2, it would appear to be more appropriate to use S2 as an order metric for quantifying pair correlation order. In this work, we have considered only the pair correlation contribution to the excess entropy. We are currently in the process of evaluating the thermodynamic excess entropy for this system. With regard to the present work, we provide a measure of the accuracy of this approximation by considering the behavior of the total entropy, S, as approximated by Sid + S2. Here, Sid is defined as:

Sid/NkB ) 1 +

Figure 6. Boundaries of excess entropy and structural anomaly for BeF2 in the in the FT plane. The locus of extrema in the S2(F) curves defines the boundaries of the region excess entropy anomaly. The loci of maxima in qtet and of minima in τ enclose the structurally anomalous region.

tetrahedral order and therefore shows a strong negative correlation with S2. In the normal regime, tetrahedral order is irrelevant, and the Be-Be interaction, which is not a nearest neighbor interaction, plays a much less critical role in determining the excess entropy leading to a weak negative correlation between S2 and τ. Figure 6 shows the boundaries of the structurally anomalous regime in the density-temperature (FT) plane as defined by the loci of maxima in qtet and minima in τ. One can also define a region of anomalous excess entropy behavior by considering loci of minima as well as maxima in Se. The region corresponding to this excess entropy anomaly, determined by locating extrema in S2, is also demarcated in Figure 6. Since qtet is negatively correlated with S2 in the low-density anomalous regime, it is not surprising that the locus of maxima in qtet is nearly co-incident with the locus of minima in S2. The locus of maxima in S2, however, does not closely follow the locus of minima in τ, and the region of excess entropy anomaly extends

∑R xR [0.5VR - ln(FRΛR3)]

(6)

The boundary of the region of thermodynamic or density anomaly must correspond to points where (∂S/∂F)T ) 0, that is, to the locus of extrema in the total entropy. One can also determine this boundary directly from the equation of state by determining the state points for which (∂P/∂T)V ) 0 without any free energy or entropy estimation. Figure 7 shows that the locus of TMD points determined from the equation of state as well as from extrema in S ≈ Sid + S2. This gives one a reasonable degree of confidence that the pair correlation approximation captures much of the essential physics. Figure 7 suggests an alternative view of the cascade of anomalies picture originally suggested by Errington and Debenedetti and discussed in the introduction. One can define a region of anomalous excess entropy behavior as well as of anomalous entropy behavior. The latter region clearly corresponds to the thermodynamically anomalous region. Given a Rosenfeld type of excess entropy scaling of diffusivities, one would expect that the loci of extrema in the excess entropy would be critical in determining the boundaries of the diffusionally anomalous region. Moreover, the monotonic behavior of Sid with density suggests that the region of the excess entropy anomaly will always enclose the region of total entropy anomaly. This is consistent with the nesting of thermodynamically anomalous regions within diffusionally anomalous regions for

Anomalies in Liquid Beryllium Fluoride all systems with waterlike anomalies that have been studied so far. A close relationship between the structural and diffusional (or excess entropy) anomalies is expected because the order metrics and the excess entropy are both measures of positional disorder in the system. However, the exact definition of the order metrics will change the boundaries of the structurally anomalous region relative to that of the diffusionally anomalous region. 4. Conclusions The order map for liquid BeF2, plotted between translational and tetrahedral order metrics, is very similar to that seen in water and liquid silica. The anomalous region in the order map corresponds to a band of state points for which tetrahedral (qtet) and translational (τ) order are strongly correlated. In water, this correlation is very strong, and the band reduces to almost a straight line while, in SiO2, the correlation in the anomalous regime is somewhat weaker than in BeF2. This provides an interesting corroboration of an earlier suggestion that the constraint of local tetrahedrality should increase in tetrahedral, network-forming liquids in the following sequence: SiO2 < BeF2 < H2O < Si.11 The normal regime and the inaccessible or forbidden region in the qtet-τ plane are well-defined. The tetrahedral order parameter distributions of BeF2 melt have been studied in order to further substantiate the analogous structural properties of BeF2, SiO2, and H2O. In each case, when the system enters the anomalous regime, the degree of bimodality increases, and correspondingly, the variance rises.7,61 The normal liquid regime at high densities, high temperatures, or both is characterized by unimodal distributions with very small third- and higher-order moments. To understand the connection between the structural, diffusional, and density anomalies, we compute the boundaries of a region of the excess entropy anomaly within which S2 increases with density. The low-density boundary of this excess entropy regime is very close to the locus of maxima in the tetrahedral order. The region of excess entropy anomaly, however, extends to greater densities than the region of anomalous structural behavior. Given the Rosenfeld-type excess entropy scaling of diffusivities and viscosities of fluids, the excess entropy and the mobility anomalies are expected to be closely connected. We also compute the boundaries of the region of total entropy anomaly, corresponding to the locus of extrema in Stot ≈ S2 + Sid. As a measure of the reliability of the pair correlation approximation in capturing the basic physics of waterlike anomalies, we show that the boundaries of the total entropy and density anomaly coincide even when the total entropy is represented as a sum of ideal and pair correlation contributions. The total entropy anomaly is nested within the region of excess entropy anomaly and provides an interesting perspective on the nesting of regions of diffusional, structural, and density anomalies. The idea of nested regions corresponding to anomalies in excess and total entropy should prove to be a useful connection to explore when studying other fluids with waterlike anomalies. Acknowledgment. C.C. would like to thank the Department of Science and Technology, New Delhi for financial support under the Swarnajayanti Fellowship scheme. M.A. thanks Indian Institute of Technology, New Delhi, for the award of a Junior Research Fellowship. Computational resources and technical support provided by Computer Science Centre, IIT-Delhi is gratefully acknowledged. References and Notes (1) Mishima, O.; Stanley, H. E. Nature 1998, 396, 329. (2) Debenedetti, P. G. J. Phys.: Condens. Matter 2003, 15, R1669.

J. Phys. Chem. B, Vol. 111, No. 46, 2007 13299 (3) Errington, J. R.; Debenedetti, P. G. Nature 2001, 409, 318. (4) Madden, P. A.; Wilson, M. J. Phys.: Condens. Matter 2000, 12, A95. (5) Scala, A.; Starr, F. W.; La Nave, E.; Sciortino, F.; Stanley, H. E. Nature 2000, 406, 166. (6) Saika-Voivod, I.; Poole, P. H.; Sciortino, F. Nature 2001, 412, 514. (7) Shell, M. S.; Debenedetti, P. G.; Panagiotopoulos, A. Z. Phys. ReV. E. 2002, 66, 011202. (8) Poole, P. H.; Hemmati, M.; Angell, C. A. Phys. ReV. Lett. 1997, 79, 2281. (9) Saika-Voivod, I.; Sciortino, F.; Poole, P. H. Phys. ReV. E 2000, 63, 011202. (10) Micoulaut, M.; Cormier, L.; Henderson, G. S. J. Phys.: Condens. Matter 2006, 18, R753. (11) Angell, C. A.; Bressel, R. D.; Hemmati, M.; Sare, E. J.; Tucker, J. C. Phys. Chem. Chem. Phys. 2000, 2, 1559. (12) Tanaka, H. Phys. ReV. B 2002, 66, 064202. (13) Jagla, E. A. J. Chem. Phys. 1999, 111, 8980. (14) Jagla, E. A. Phys. ReV. E 2001, 63, 061509. (15) Yan, Z.; Buldyrev, S. V.; Giovambattista, N.; Stanley, H. E. Phys. ReV. Lett. 2005, 95, 130604. (16) Kumar, P.; Buldyrev, S. V.; Sciortino, F.; Zaccarelli, E.; Stanley, H. E. Phys. ReV. E 2005, 72, 021501. (17) Yan, Z.; Buldyrev, S. V.; Giovambattista, N.; Debenedetti, P. G.; Stanley, H. E. Phys. ReV. E 2006, 73, 051204. (18) de Oliveira, A. B.; Netz, P. A.; Colla, T.; Barbosa, M. C. J. Chem. Phys. 2006, 124, 084505. (19) Wilding, N. B.; Magee, J. E. Phys. ReV. E 2002, 66, 031509. (20) Gibson, H. M.; Wilding, N. B. Phys. ReV. E 2006, 73, 061507. (21) Likos, C. N.; Hoffmann, N.; Lo¨wen, H.; Louis, A. A. J. Phys.: Condens. Matter 2002, 14, 7681. (22) Sastry, S.; Debenedetti, P. G.; Sciortino, F.; Stanley, H. E. Phys. ReV. E 1996, 53, 6144. (23) Rebelo, L. P. N.; Debenedetti, P. G.; Sastry, S. J. Chem. Phys. 1998, 109, 626. (24) Hoyt, J. J.; Asta, M.; Sadigh, B. Phys. ReV. Lett. 2000, 85, 594. (25) Rosenfeld, Y. Phys. ReV. A 1977, 15, 2545; Chem. Phys. Lett. 1977, 48, 467. (26) Rosenfeld, Y. J. Phys.: Condens. Matter 1999, 11, 5415. (27) Rosenfeld, Y. Phys. ReV. E 2000, 75, 7524. (28) Dzugutov, M. Nature 1996, 381, 137. (29) Sharma, R.; Chakraborty, S. N.; Chakravarty, C. J. Chem. Phys. 2006, 125, 204502. (30) Mittal, J.; Errington, J. R.; Truskett, T. M. J. Chem. Phys. 2006, 125, 076102. (31) Mittal, J.; Errington, J. R.; Truskett, T. M. J. Phys. Chem. B 2006, 110, 18147. (32) Errington, J. R.; Truskett, T. M.; Mittal, J. J. Chem. Phys. 2006, 125, 244502. (33) Mudi, A.; Chakravarty, C. J. Phys. Chem. B 2004, 108, 19607; 2006, 110, 4502. (34) Mudi, A.; Chakravarty, C.; Ramaswamy, R. J. Chem. Phys. 2005, 122, 104507; 2006, 124, 069902. (35) Mudi, A.; Chakravarty, C J. Phys. Chem. B 2006, 110, 8422. (36) Mudi, A.; Chakravarty, C.; Milotti, E. J. Chem. Phys. 2006, 125, 074508. (37) Sharma, R.; Mudi, A.; C. Chakravarty J. Chem. Phys. 2006, 125, 044705. (38) Johnson, M. E.; Head-Gordon, T.; Louis, A. A. J. Chem. Phys. 2007, 126, 144509. (39) Louis, A. A. J. Phys.: Condens. Matter 2002, 14, 9187. (40) Truskett, T. M.; Torquato, S.; Debendetti, P. G. Phys. ReV. E 2000, 62, 993. (41) Greenwood, N. N.; Earnshaw, A. Chemistry of the Elements; Butterworth-Heinemann: Oxford, 1984. (42) Woodcock, L. V.; Angell, C. A.; Cheeseman, P. J. Chem. Phys. 1976, 65, 1565. (43) Hemmati, M.; Moyinihan, C. T.; Angell, C. A. J. Chem. Phys. 2001, 115, 6663. (44) Hemmati, M.; Angell, C. A. J. Non-Cryst. Solids 1997, 217, 236. (45) Busing, W. R. J. Chem. Phys. 1972, 57, 3008. (46) van der Meer, J. P. M.; Konings, R. J. M. J. Nucl. Mater. 2007, 360, 16. (47) Heaton, R. J.; Brookes, R.; Madden, P. A.; Salanne, M.; Simon, C.; Turq, P. J. Phys. Chem. B 2006, 110, 11454. (48) Salanne, M.; Simon, C.; Turq, P.; Heaton, R. J.; Madden, P. A. J. Phys. Chem. B 2006, 110, 11461. (49) Tosi, M. P.; Fumi, F. G. J. Phys. Chem. Solids 1964, 25, 45. (50) Tosi, M. P.; Price, D. L.; Saboungi, M.-L. Ann. ReV. Phys. Chem. 1993, 44, 173. (51) Agarwal, M.; Sharma, R.; Chakravarty, C. J. Chem. Phys. 2007, in press.

13300 J. Phys. Chem. B, Vol. 111, No. 46, 2007 (52) Hemley, R. J.; Mao, H. K.; Shen, G. Y.; Badro, J.; Gillet, P.; Hanfland, M.; Hausermann, D. Science 1997, 276, 1242. (53) Ahart, M.; Yarger, J. L.; Lantzky, K. M.; Nakano, S.; Mao, H. K.; Hemley, R. J. J. Chem. Phys. 2006, 124, 014502. (54) Smith, W.; Forester, T. R. J. Mol. Graphics 1996, 14, 136. (55) Smith, W.; Yong, C. W.; Rodger, P. M. Mol. Simul. 2002, 28, 385; http://www.cse.clrc.ac.uk/msi/software/DL_POLY/. (56) Allen, M. D.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1986.

Agarwal and Chakravarty (57) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic Press: London, 2002. (58) Baranyai, A.; Evans, D. J. Phys. ReV. A 1989, 40, 3817. (59) Chakraborty, S. N.; Chakravarty, C. Phys. ReV. E 2007, 76, 011201. (60) Laird, B. B.; Haymet, A. D. J. J. Chem. Phys. 1992, 97, 2153. (61) Mudi, A. M., Ph.D. Thesis, Indian Institute of Technology Delhi.