Properties of Liquid Water: Origin of the Density Anomalies - The

Mar 1, 1994 - Mary Vedamuthu, Surjit Singh, G. Wilse Robinson ... Stanley , Hajime Tanaka , Carlos Vega , Limei Xu , and Lars Gunnar Moody Pettersson...
2 downloads 0 Views 2MB Size
2222

J . Phys. Chem. 1994,98, 2222-2230

FEATURE ARTICLE Properties of Liquid Water: Origin of the Density Anomalies Mary Vedamuthu, Surjit Sin& and G. Wilse Robinson’ Subpicosecond and Quantum Radiation Laboratory, Departments of Physics and Chemistry, Box 41 061, Texas Tech University, Lubbock, Texas 79409-1 061 Received: September 30, 1993; In Final Form: December 10, 1993”

-

-

Strong support for a “mixture” model of liquid water can be found from an analysis of the accurate experimental density data (H20) over the range t -30 “C in the supercooled regime to t +70 OC. Published density data can be fit to this mixture model with six- to seven-decimal-point precision. Remarkably, the output parameters from these fits indicate the presence of capacious intermolecular bonding with a density extremely close to that of ordinary ice-Ih, intermixed with compactly bonded regions having a density near that of the common dense forms of ice, in particular ice-11. Densities a t higher temperatures could also be fit to good precision with such a model, though the model must clearly become less valid as the temperature rises and more varied bonding forms contribute. The fitting procedure also shows that both the capacious and dense components have positive thermal expansion coefficients that are similar in magnitude to those of their respective ice forms. As T approaches the vicinity of 225 K in the deep supercooled regime, the structure of the liquid approaches disordered ice-I-type bonding, with no contribution from the densely bonded component. Combined with the differential X-ray scattering data of Bosio, Chen, and Teixeira on liquid water, and structural data on the ice polymorphs from Kamb’s work, it can be concluded that the bonding differences between the dense and capacious structures are not a t the nearest-neighbor level but occur instead in the outlying non-hydrogen-bonded nextnearest-neighbor O--O structure. Because of the long-range structural implications of this conclusion, uncertainties arise in molecular dynamics modeling of the liquid and on the usefulness of attempts to learn about the liquid from the study of small gas-phase clusters.

I. Introduction

A little over 6 years ago, a semiquantitative explanation of the various properties of liquid water was presented in this journal by Bassez, Lee, and one of the present authors.’ That paper attempted to reveal the origins of the temperature and pressure dependence of the viscosity and of various rate processes, such as spin-lattice and dielectric relaxation. Additionally dealt with were the relationships connecting these dynamic properties of water with “anomalous” characteristics such as the density and the isothermal compressibility. Apparent in this analysis was the presence of a “feedback” of cooperativity phenomenon that enters the problem in the same way that is does for rotations in solids:* when water molecules begin to “rotate”, the local fields become more isotropic, and the activation barrier for the rotation of other molecules in the vicinity becomes lower, creating still more rotating molecules. In spite of the importance of and the general interest in liquid water, the successes of the original paper have not been extended very far by ourselves, nor have they been exploited by others. Though the earlier results were suggestive, one reason for the lack of exploitation is that the ideas on which that paper was based had already been discredited by Kauzmann3 and otherse4 In fact, Kauzmann in 1975 wrote,3 “Since the time of Rontgen (our ref 5 ) several of the unusual properties of water have been qualitatively interpreted in terms of mixture models involving an equilibrium between ‘bulky’ and ‘dense’ components, or their equivalent.” Kauzmann then went on to say, “There are, however, a number of aspects of the thermodynamic behavior of water, particularly its response to pressure, which are not easy to explain by these models.” Kauzmann then concluded, “There is considerable reason to question the basic validity of these theories.” 0

Abstract published in Advance ACS Absrracrs, February 1, 1994.

It develops, however, that essential data in thesupercooled region were not available in 1975. Nevertheless, there remains today a general feeling that a complete understanding of the properties of liquid water is unfathomable and in any case cannot be found through “mixture model” concepts. The purpose of the present paper, and future papers in this series, is to try to get the venerableS mixture model for liquid water back onto its feet. As more experimental data accumulate, it becomes more and more difficult to ignore the accurate quantitative assessments of the properties of this substance that can be obtained from this simple model. A premiere set of experimental data on water, the set quantitatively most amenable to a precision description because of its high level of accuracy and extensive thermodynamic range, is the liquid density, which is considered an “anomalous” property. The density of HzO will thus constitute the main focus of the present paper. Future papers will deal with pressure effects, including the isothermal compressibility, isotope effects, and finally the dynamical properties of water. Preliminary work186 in these areas indicates that all of these features of liquid water can be quantitatively understood using the model outlined in the present paper. 11. Qualitative Descriptions

An explanation of the anomalous density maximum in liquid water was the primary motivation for the initial formulation of “mixture models”.5 Secondary motivations,s of which workers were aware in 1891, were “the minimum in the isothermal compressibility near 50 “C, the fact that the density anomaly moderates and the maximum moves to lower temperatures as the pressure is increased, and the decrease in viscosity with rising pressure near 18 O C ” . In the mixture model for liquid water described in this paper, there are supposed to be two or more general types of intermo-

0022-3654/94/2098-2222%04.50~0 0 1994 American Chemical Society

Feature Article lecular bonding configurations, a “bulky” or “capacious” bonding form, with a low density, such as occurs in ice-Ih, plus a dense bonding form, such as occurs in the most thermodynamically stable dense forms of ice, e.g., ice-11, -111, -V, and -VI? As the temperature T (or pressure P) is raised, the capacious structure, held intact by open tetrahedral hydrogen bonding, becomes relatively less stable and begins to break down, creating more of the dense structure. At the beginning, it is important to stress that this model is not conceived to be a “mixture of ices” having even modest longrange order, but rather is a rapidly fluctuating mixture of intermolecular bonding types found in the most stable polymorphs of ice. The “long-range” order, except at very low temperatures of supercooling, is conceived to extend generally no farther than about two or three van der Waals diameters (C-9 A). Anomalous effects accompanying changes of the temperature and pressure in the properties of any substance will be evident if gross structural and quantum-state characteristics of that substance substantially vary with T and P. The breakdown of the capacious structure in liquid water is a good example. Occurring then would be deviations from the situation commonly met in normal systems, such as an ensemble of gaseous diatomic molecules, where thermal population of a fixed array of energy levels takes place. As a matter of fact, it is well established experimentally that intermolecular librational frequencies for liquid water strongly vary with temperature.E-’O The observed decrease of these frequencies with increasing temperature is a corollary to the idea, a version of which was first proposed by Segr6,” that the intermolecular potential surfaces in liquid water cooperatively soften and flatten out with increasing T (and P), an effect that is very likely related to the Pauling-Fowler mode12J2J3 for cooperative librational melting in solids. Such temperature (pressure)-dependent potential surface changes are surely responsiblelJ4 for the anomalous heat capacity in liquid water15 (see ref 2 for some solid-state examples). The potential surface flattening, which accompanies the breakdown of the capacious structure, gives rise to a greater fluidity for the dense structures, contrary to ordinary intuition. The strong isotope dependence1V6in the behavior of many properties of liquid water is caused by the fact that the proportions of the capacious and dense bonding structures at any given temperature depend on frequencies and zero-point effects in the associated temperature-dependent intermolecular potentials, particularly those related to molecular rotational librations, which have the largest H, D isotope effects. As the temperature is raised from the supercooled regime, the capacious bonding weakens, and the structure fluctuates more and more between the diverse, but energetically similar (see Table 3.6 of ref 16), types of intermolecular bonding. The rotational cooperativity effect, and the concomitant potential surface flattening, gives rise to non-Arrhenius temperaturedependent activation energies, a well-known property of all dynamical effects in liquid water. See, for example, refs 1 and 17-20. This cooperativity effect, of course, dies out at higher temperatures (or higher pressures), where more of the higher density, less viscous type of bonding is already present, whence the dynamic properties of water begin to behave more normally with variations in T o r P. Increasing the pressure a t a constant temperature diminishes the proportion of the bulky ice-I-like regions and thus reduces the probability for the formation of homogeneous nuclei. This allows pressurized water to be supercooled to much lower temperatures than normal ~ a t e r . 1This ~ ~pressure ~~ effect also flattens intermolecular librational potential surfaces in water, mollifying the temperature effect, thus causing the temperature anomalies to tend to disappear at elevated pressures.5J5J7,21 We have always maintained14J2 that the breakdown of the open tetrahedral structure in liquid water is not caused by a deformation of nearest-neighbor hydrogen bonding. In fact, this

The Journal of Physical Chemistry, Vol. 98, No. 9, 1994 2223 same view had already been clearly described by Kamb7 in a much earlier paper. “We may expect that the increase in density (of the liquid with increasing temperature) is accomplished primarily by a structural change similar to the change that allows the closer molecular packing of the denser solid forms. -.The increase in density is thus achieved not by an increase in the coordination number of bonded neighbors, but instead by an increase in non-nearest-neighbor coordination by a bending of the hydrogen bonds.” In this connection, it is important to note from ref 7 (see also Table 3.4 of ref 16) that, in the four most stable-P I 103 MPa (10 kbar), t I 0 OC-high-density polymorphs, ice-11, -111, -V, and -VI, every 0 atom has four nearest neighbors a t a distance of about 2.75 A, as in ice-Ih (or ice-IC), but in these dense structures there are non-hydrogenbonded O.-O neighbors lying in the range 3.24-3.51 A. This is a distance not represented in either ice-Ih or ice-IC, so it acts as a “fingerprint” for the presence of the dense structures. The squeezing together, by bond bending, of these non-nearestneighbor 0 atoms (from -4.5 %I in ice-Ih) is the major cause of the >25% increase in density of these polymorphs. Also, in iceI11 there are five-membered rings and in ice-V and -VI, fourmembered rings. Mixing these structures would provide a plethora of local bonding configurations in the liquid state as the temperature is raised. Such a multiplicity of bonding could also play a central role in the problem of interfacial or surface-pertubed water, including ionic solutions. Depending on the surface, the water molecules would have to decide which of all these possible bonding types give the optimum free energy: “structure making” surfaces promoting the open tetrahedral bonding, but “structure breaking” surfaces preferring the denser, less viscous bonding of ice-11, etc. As a matter of fact, it would be quite amazing if bonding characteristics in crystalline phases did not carry over to the liquid state, and Kamb cites other examples in his paper.7 While close to our own views, Kamb’s picture is a t variance with those presented in most other descriptions of liquid water, which depict structural changes in terms of nearest-neighbor hydrogen bond disfigurations. See, for example, the discussions in Chapter 5 of ref 16, which cover most prototype water models. In thenext two sections,a quantitativedescription of the density anomaly in liquid water will be presented. The resulting picture will be seen to be indeed a “mixture model”, but certainly not a “two-state model”3 because of the probable presence of more than a single high-density bonding type and the dependence of the states themselves on T and P.

-..

--.

..-

111. The Density Anomaly

In the mixture model picture to be described here, the expression for the pressure/temperature-dependent specific volume of liquid water is elementary. It is

where the parametersf represent the mass fractions as a function of P and T of the capacious and dense components, which from now on are labeled I and 11, respectively. For simplicity, it will be assumed that there are no components of intermediate density, so thatfIl = 1 -fI. This assumption does not rule out the presence of a variety of structures as long as such structures can be grouped according to their local densities into one of the two categories described above. In this connection, there is no known ice polymorph that has a density between p = 0.92 g cm-3 (ice-Ih) and those of the most stable high-density forms:’ 1.18 (ice-11), 1.15 (ice-HI), 1.26 (ice-V), and 1.34 g cm-3 (ice-VI). According to the previous section, fir might be expected to approach zero as some low temperature is approached and to increase with increasing Tor P . For eq 1 to be convincing, both forms I and I1 must have normal positive thermal expansion coefficients, expected to be in the neighborhood of 10-3-10-4 K-1, as is the case for ice-Ih.16 Thus, the specific volumes of both the

2224 The Journal of Physical Chemistry, Vol. 98, No. 9, 1994

Vedamuthu et al.

TABLE 1: Liquid Water Densities (g cm") at 1 atm from Kell's Equation Compared with Experimental Supercooled Densities and Fit i Densities t, O C KelP Zheleznyib Schufle HSd fit i -34 -32 -30 -25 -20 -1 5 -10 -9 -8 -7 -6

-5 -4

-3 -2 -1 0 +4 +10 +15 +20 +25 +30 +40 a

0.98385 0.98959 0.993 55 0.996 28 0.998 12 0.99840 0.99865 0.99887 0.999 08 0.99926 0.999 41 0.999 55 0.999 67 0.999 76 0.999 840 0.999972 0.999700 0.999 100 0.998204 0.997045 0.995647 0.992216

0.9751 0.9793 0.9829 0.9895 0.9936 0.9964 0.9983 0.9986 0.9989 0.9991 0.9993 0.9995 0.9996 0.9997 0.9998 0.9998 0.9999 1 .om

0.993 36 0.996 16 0.99796 0.998 25 0.9985 1 0.99875 0.998 98 0.999 23 0.99938 0.999 53 0.99967 0.99977 0.999 85 1 .oooo

0.9775 0.9809 0.9839 0.9895 0.9935 0.9963 0.9982 0.9984 0.9987 0.9989 0.9991 0.9993 0.9994 0.9996 0.9997 0.9998 0.9999

0.97736 0.98087 0.983 85 0.989 59 0.993 55 0.996 28 0.998 12 0.998 40 0.99865 0.99887 0.999 08 0.999 26 0.99941 0.999 55 0.999 67 0.999 76 0.999 840 0.999972 0.999 700 0.999 100 0.998 206 0.997050 0.995659 0.992255

Reference 24 . Reference 27. Reference 26. Reference 28.

imantZ5 containing seven adjustable parameters. The output I and I1 bonding forms in eq 1 increase as the temperature is values of this equation for 133 temperatures between -30 and raised. Why then does the specific volume V(P,T)of the mixture 150 "C are listed in his Table 111. These tabulated data provide concomitantly decrease; Le., why does the density increase with probably one to two decimal point precision beyond the actual increasing temperature in the range 238-277 K, giving a density experimental accuracy. For this reason, in the presentation of maximum (H2O) near the latter temperature? the remainder of our paper, the final decimal point in the Kell To help answer the above question, it can be seen from eq 1 that ( a v / d T ) ~ =[(VII- V,) a f r t / e T + f I a V , / a T + f I t d V ~ ~ / e T l ~ .data will usually be rounded, though it was retained in all our fitting procedures. The first term in the square brackets is negative, since VI > VII, and the fractionfrl of dense structure increases with increasing For the supercooled regime, Kell states that his equation does T. The other two terms have positive signs because of normal not extrapolate to low-temperatures well. Then, in the footnote thermal expansion. The magnitude of the first term dominates to his Table 111, he further declares that, "Values below 0 O C the magnitudes of the other terms at low temperatures because were obtained by extrapolation, and no claim is made for their of the steep increase of fr1 with increasing T . The first term accuracy." Kell's densities24in this important temperature range therefore gives a dominating negative contribution to the volume are compared in Table 1 with published experimental data.2628 change a t these temperatures. Simply speaking then, it is the A number of observations can be made. strongly increasing "proportionality" of I1 that gives rise to the If there are no remaining systematic errors in the most recent, initial decrease of volume with increasing temperature and thus and probably best, experimental density measurements for to the volume minimum (density maximum). The behavior of supercooled water,28 Kell's extrapolated data are apparently much the pressure derivatives of eq 1 at each temperature T, noting better than he realized. Systematic errors in the very difficult that f r t increases with increasing pressure but less so a t higher temperature range below -20 OC can arise from sample size effects. temperatures, will be seen in later work, without calling on any Such effects are now well established.28-31 The density, as well singular phenomenon,23 to be responsible for the minimum in the as other measured properties of water, does vary with the sample isothermal compressibility, -V-l(aV/aP), at 46.5 O C in liquid size. In particular, for capillary or emulsion measurements, the water.24 usual techniques for supercooled water, the measured apparent Thus, the density of liquid water a t deep supercooled temdensity increases with decreasing sample size below about 50 pm. peratures, because of the "proportionality effect", anomalously The early measurements of Mohler3* to -13 "C employed increases as the temperature is raised. Eventually, however, as capillaries whose diameters ranged from 188 to 380 pm, so those the proportion of I1 becomes larger, the increase off1 "saturates", data should present minimal size effects. In fact, Mohler's data the magnitude of the negative (VII- VI) (afIt/a7'),p term agree within about (1-2) X 10-4 g cm-3 with the single data run diminishes, the positive terms begin to dominate, and the density reported by Schufle26using a 133-pm-diameter capillary. Schudecreases with increasing temperature, as is common for normal fle's density data extended Mohler's to -21 OC, but compared substances. Applying pressure at a constant temperature gives with the other data in Table 1, they are seen to be a bit too small more of component 11, lessening the density anomaly, an effect a t the lower temperatures. that is observed experimentally for liquid water15s21*24 and was More recently, Hare and Sorensen reported size effects, when found to emerge from our initial analysis.' the capillary diameter was changed from 253l to 300 pm:28 of the order of a few parts in 104 for density measurements down to -30 IV. Density Data O C and greater deviations down to -34 O C . The often-cited densities of Zheleznyi,27 using small capillaries, seem too large The aim here is to fit closely the smoothed experimental density between -2 and -10 O C and then become too small by >1 X 10-3 data of KelP4 a t atmospheric pressure. The number of decimal g cm-3 compared with the H a r d o r e n s e n densities at temperplaces provided by this author is six from -30 to -1 OC and from atures below about -30 OC. This is in the wrong direction for + 101to +150 OC and seven decimal places for the data between a capillary size effect. Interestingly, the H a r d o r e n s e n data28 0 and +lo0 "C. Kell defines, as we do here, 100 OC as precisely from their 300-pm capillary agree within (1-2) X 10-4 g cm-3 373.15 K. Kell's density smoothing function is a Pad6 approx-

The Journal of Physical Chemistry, Vol. 98, No. 9, 1994 2225

Feature Article with Kell's extrapolated data.24 For this reason, we select the K e l P smoothed extrapolated data and the Hare-Sorensen 300pm experimental data28 as the most reliable throughout the supercooled regime.

V. Density Analysis Armed with the data and the concepts in sections I1 and 111, precise quantitative fits of the temperature-dependent density data for HzO at atmospheric pressure can now be carried out. To accomplish this, analytical forms for the pressure/temperaturedependent quantities VI(P,T),&I(P,T), andfrl(P,T) in terms of a set of temperature-independent fitting parameters at fixed atmospheric pressure are required. Standard forms for VI and &I are employed for use in eq 1,

Seven adjustable parameters are hereby introduced: VI(P,TO), &I(P,To), M'), aldP), PdP), BII(P),and To(P). Because of the properties expected from the discussions in sections I1 and I11 of the proportionality factor AI, some such form as the following

fII(P,T)= D tanh

[

+

A ( T - To) B ( T - T0)Z 1 +C(T-To)

can be selected. It is seen that this mathematical form has the desired temperature behavior, approaching zero A , B, and Cwhen T + TOand reaching a constant value, D 1 1, as the temperature is raised. This form is certainly not unique and in fact is somewhat arbitrary. It is the simplest of other successful forms that gave the six- to seven-decimal-point precision required for the density fits. The pressure-dependent parameters A , B, and C in themselves have little or no physical significance. All the physical significance resides in the quantity f,,(P,T). Equation 4 thus introduces four more adjustable parameters, giving in all 11 adjustable parameters for each pressure P,four more than in the Kell smoothing equation. Clearly, this number of adjustable parameters, even considering the six- to seven-decimal-pointaccuracy of the smoothed density data, is too large. Any fitting procedure involving such a large set of parameters would be subject to errors, particularly in the smaller a,8, and B parameters, because of mutual compensation effects, as the fits attempt to match the meaningless zeros beyond the reported number of significant figures in the experimental data. In many preliminary and test fits, we have indeed found this to be the case. Thus, it will be necessary to find ways of reducing the number of adjustable parameters without in any way forcing a fit toward the proposed mixture model. The parameter D reflects the "degeneracy" of dense component configurations. If there is only onecapacious component I, having open tetrahedral bonding, one dense component 11, and no other structural forms, D would be equal to l/2, since the populations of the dense and capacious components would become equal at large T. However, this is very likely not the case for water, since there are expected to be at least four distinct (Le., having different densities) higher-density bonding forms in the liquid, corresponding to the bonding configurations of the most stable ice polymorphs. Probably then, D should be closer to 0.8. However, as the temperature rises, all the specific bonding forms must give way more and more to fully disordered structure. As a matter of fact, the best preliminary fits up to -40 OC always showed D to be very close to unity, which, taken literally, would mean that the open tetrahedral bonding tends to disappear entirely at

high temperatures in lieu of more densely packed and disordered structures. It was therefore elected to set D equal to unity in the fits. However, we found that, for all practical purposes, it really does not matter too much whether D has a value of ' / 2 , 1 , or any value in between, since the A , B, and C parameters were found to compensate,so that the resulting physically significantfI1values for a given temperature always turned out to be about the same. This fact further illustrates that the form o f h i is immaterial to the fitting procedure, providing of course that the chosen form has the right mathematical properties. It is the value offil itself as a function of temperature (and pressure) that really matters, and eq 4 is just one of many forms that could have been chosen to express the temperature-dependent fI1 in terms of a set of temperature-independent fitting parameters. Thus, excluding D,there still remain a total of 10 adjustable parameters for the density fits. Through a large number of preliminary data fits, it was found that the parameter TOusually lay between 224 and 230 K. In the present mixture model description, this temperature is the temperature at which all intermolecular bonding in liquid water resembles that in disordered or amorphous ice-Ih or -1c.33 Very interestingly, these To's from the density fits are not far from the "singular temperature", 228 K, found by Angel1 and co-workers15 from a variety of independent considerations. Of further interest are the output values of the two volumes, VI(P,TO) and VII(P,TO), together with their thermal expansion coefficients. The values offII(P,T), through the parameters A, B, C, and TO,are also interesting, since this quantity should have connections' with modern experimental data in the supercooled regime: radial distribution functions,34small-angle X-ray scatterir~g,~s and activation energie~l~.~O.~'j-~* for dynamical processes as a function of P and T. Thus, it is seen that the output Vs, a's, and Ps, in addition tofI1 and TO,but not the particular mathematical form offiI, will have physical significance for the density discussed in this paper and for other properties of liquid water to be discussed in future papers. The pressure dependence of these parameters, for example, will be required in future discussions of the isothermal compressibilities. Since only atmospheric pressure is being considered here, the pressure dependence of all the variables will be suppressed in the remainder of the paper. The primary analyses consider Kell's data only over the temperature range -30 to +70 OC. Densities for temperatures above 70 OC were alsoexamined. The high-temperature mixture model fits were found to give adequate agreement with the density data, but the physically significant quantities clearly begin to vary as more high-temperature data are included in the fits. The greatest sphere of reliability of the model, and thus its main interest, should reside in the lower-temperature regions anyway, where not too many intermolecular bonding forms play a role and where fully disordered structure is at a minimum. In spite of such problems, the major parameters -VI(TO), VII(TO),TO, andfII-did not vary more than about 6% in a variety of fits tried for temperatures up to 100 OC,so a remnant of the mixture model must exist even at fairly high temperatures. Among the approximately 60 fits employed, only four are reported here in detail: (i) an exact low-temperature (14 OC) algebraic fit of eight selected density data points to determine the eight parameters A , B, C, TO,VI(TO),VII(TO),a],and ~ I I with , D = 1 and the (I's set equal to zero for this very low-temperature regime; (ii) a least-squares fit of 39 data points from -30 to +40 OC (the number of data points being equally balanced above and below +4 "C) to attempt to refine, for this higher temperature regime, the values of the four basic volume parameters-VI( To), VI](TO),a],and cqr-when the A , B, C, and TOparameters, and thus the values offII(T), are fixed at their fit i output values; (iii) a least-squares fit similar to fit ii, but with the A , B, C, and TO parameters, not from fit i, but deduced from temperaturedependent activation energies for dynamic processes as outlined in ref 1; and (iv) a least-squares fit, with the same fixed fit ii

2226

The Journal of Physical Chemistry, Vol. 98, No. 9, 1994

TABLE 2

Vedamuthu et al.

Fitting Parameters from Fits i-iva

A

B C

To WTo) PI(TO) VIdTo) PII(T0) a1 a11

fit i

fit ii

4.54866E-2b 3.6522E-4 8.691 96E-2 225.334 1.0876 1 0.9195 0.84632 1.1816 4.78005E-4 1.29473E-3

4.54866E-2‘ 3.6522 E-4‘ 8.69 196E-2‘ 225.334‘ 1.08739 0.9196 0.84745 1.1800 4.57558E-4 1.29374E-3

fit iii 4.529 12E-2‘~~ 4.3486E-4‘~~ 9.56974E-2‘sd 225.334‘ 1.09180 0.9159 0.8 1340 1.2294 8.37014E-4 1.49776E-3

fit iv 4.54866E-2‘ 3.6522E-4‘ 8.69 196E-2‘ 225.334‘ 1.08717 0.9198 0.83369 1.1995 1.09085E-3 3.39868E-3 7.4987E-6 3.2862E-6

0‘ 0‘ oe 0‘ oc 011 oe a The units are as given in the text. Fit i: exact algebraic fit of 8 Kell densities between -30 and +4 OC. Fit ii: least-squares fit of 39 Kell densities between -30 and +40 OC. Fit iii: least-squares fit, with fixed A, B, and C from dynamic data, of 39 Kell densities between -30 and +40 O C . Fit iv: least-squares fit of 39 Kell densities between -30 and +70 OC. Read as 4.54866 X 1 P 2 . Parameter fixed. Parameter from dynamic data. PI

parameters, of 39 data points between -30 and +70 OC to determine best values of all six volume parameters (Vs,a’s, and Ps) when a higher temperature regime is considered. A briefer description of some of the other fits examined is included in section VIII.

and the fact that TOis close to but somewhat lower than the “singular temperature”15 at 228 K. Recalculated densities from these fit i parameters for 24 temperatures between -34 and +40 O C arecompared with theexperimental datazGz8 and thesmoothed KellZ4 data in Table 1.

VI. Exact Algebraic Fit

VII. Least-Squares Fits

Fit i. Of the fits mentioned at the end of the last section, fit i is perhaps the most interesting, since it provides unequivocal connections between the density data and the mixture model of eq 1. As already mentioned, the model is expected to be the most accurate at lower temperatures because the multiplicity of intermolecular bonding structures is minimized. Such an algebraic fit cannot in any way be “forced”, since eight density points are to be fit exactly by equations containing eight unknowns, with no freedom for manipulation or choice of the output parameters. Of course, if the experimental data to be fit have errors of a few percent, it would be expected that parameter variations of at least that magnitude can occur without greatly modifying the “exact” algebraic solution. Thus, the approach to a unique set of correct output parameters is only possible from “noiseless” experimental density data, such as Kell’s smoothed data,providing there is at least one-decimal-point accuracy in the densities beyond that desired for the parameters. Actual experimental data in the supercooled regime that could be used to calibrate Kell’s extrapolated data to the required level of accuracy have only recently become available.28 For fit i, we use only data from the supercooled regime up to +4 OC. This temperature range is one that provides the most rapidly increasing value offrI( 7‘). Taking eight uniformly spaced data points in this temperature range to create the eight required equations, we simply solve algebraically for the eight parameters in eqs 1-4 (01= 011= 0, D = l), accepting only the physically relevant solution from the nonlinear equation set. Specifically, the K e l P density data at -30, -25, -20, -15, -10, -5,0, and +4 OC were chosen. Since it is to be remembered that Kell’s data are from a smoothing equation and only round-off “noise” in the final decimal place (one beyond those given in Table 1) is present, the chosen procedure should give accurate low-temperature liquid parameters for comparison with the density and the thermal expansion coefficients of the ice polymorphs, and, in future papers, with other properties of water and ice. The resulting highprecision parameters effectively just “drop out” of such a fit, which of course matches exactly the six- to seven-decimal-point density precision given by Kell for these eight data points in this temperature range. Furthermore, the agreement remains excellent for other temperatures in and somewhat beyond this range as well. See Table 1, observing that the deviations at higher temperatures are in the expected direction for the missing terms. The resulting parameters for fit i are given in Table 2. Note here the near agreement of the output p’s with ice-Ih and ice-I1 densities

When more experimental data points than the number of parameters are considered, a least-squares fitting procedure is employed. See the Appendix. A large number of such fits have been performed. One concern is that Kell’s data set above about 25 OC contains many more density points above the temperature of maximum density than below it. Therefore, the dominating high-temperature data could distort the more interesting limiting low-temperature parameters. A more serious problem with any fitting procedure of this type is the “mutual compensation” effect mentioned earlier. With all the fitting parameters variable, even the six- to seven-decimal-point agreement with the Kell data is an insufficient criterion for obtaining unique values for the parameters, particularly for the minor parameters, a’s, Ps, and B. For example, in a number of “blind fits”, where all eight parameters are allowed to vary, VI(To) was found to range from 1.05 to 1.10 cm3 g-I, while fir( TO)ranged from 0.80 to 0.85 cm3 g-I, while still reproducing good agreement with the density data. The thermal expansion coefficients (YI and a11 changed more radically, usually between about 5 X 10-3 and 1 X 10-5 K-1, but were sometimes even negative, as these two small quantities attempted to compensate for variations in the two volumes. See next section. In these “blind” least-square fits, it was also found that, even though the A, B, and Cparametersvarieda good deal, theresulting values offII( T) were pretty stable and agreed within better than 10% with the values obtained from the fit i parameters. Thus, it was decided that, in order to reduce the variable parameters to a manageable number, A, B, C, and TOcould be fixed in the least-square fits. Fit ii. This least-squares fit of Kell’s density data for temperatures 1+40 O C employs just four adjustable parameters, setting the @’s to zero, D to unity, and fixing the A, B, C, and TO parameters at their fit i output values. The analysis utilizes 39 of Kell’s density data points in the range -30 to +40 O C . Only every fifth data point above +20 O C is selected in order to give an approximately equal number of data points above and below the temperature of maximum density. The output parameters from fit ii are listed in Table 2. Fit iii. In ref 1, it was found that both the densities and dynamic data could be fit to rather high temperatures with the current mixture model, though the precision demanded there for the densities was much less (a few percent) than we wanted to achieve in the present paper. In any case, we felt it still might be worthwhile to follow the lead provided in ref 1 to see whether a procedure similar to the one there could also provide highly

The Journal of Physical Chemistry, Vol. 98, No. 9, 1994 2221

Feature Article

TABLE 3: Recalculated Densities for Least-Squares Fits ii-iv Compared with Those from Ref 24* t , oc Kellz4 fit ii fit iii fit iv -30 -25 -20 -15 -10 0 +4 +10 +20 +30 +40 +60 +80

+loo

0.983 85 0.989 59 0.993 55 0.996 28 0.998 12 0.999840 0.999 972 0.999700 0.998 204 0.995647 0.992216 0.983 199 0.971 798 0.958364

+1 -1 -1 0 0

-1 +I +1 +1 0

0

+1

-4

+2 0 -3

+1

+4

-1 0 0 0 -1 0 -1 0

-4

+1

+1

-7 +5 +174 +705

0 -1 +8 +7 1

+7 +135 +547 +1457

-2

Shown in columns 3-5 are density differences in the final tabulated decimal place between the Kell data and the recalculated densities using the parameters from the three fits. For example, for t < 0 O C , +2 means that the recalculatedvalueisgreaterthan Kell’s by twoin the fifthdecimal place; for t 1 0, the difference is in the sixth decimal place. Of course, differences are affected by round-off to & 1. 0

accurate densities. A measure of self-consistency with the dynamic data would be retained in this way, and some idea about the effect of fixing the four internalfil parameters, as was done in fit ii, can be gained. Therefore, for this set of least-squares fits, the parameters A, B, and Care obtained from dynamic data according to ref 1, where it was assumed39 that f i ( T ) = AH+(T)/AH+(To),the ratio of activation energies at T and TO. To carry out such a fit, the temperature-dependent dynamic rates from ref 1 were smoothed, and AH(T0) was chosen to be 14 820 J mol-’. A preliminary least-squares fit to eq 4 of the dynamicfi( T ) = 1 -fi~( T ) values in the temperature range 238423 K was first performed to obtain TO,A , B, and C. The TO value found was 224.8 K, again somewhat lower than the Angelli5 singular temperature of 228 K but matching fairly closely the TO value obtained from fit i. Since these values are so close to one another, it wasdecided tofix Toat 225.334 Kfor thedetermination of the A , B, and C parameters from eq 4 using least squares. These dynamic A , B, and C parameters along with the fixed TO = 225.334 K (see Table 2) were then used in the density fits, with only V,( TO),VI[(TO),a[,and a11being varied. This fit is called fit iii, and its output parameters are listed in Table 2. Fit iv. This fit extends the density analysis to +70 OC. For both of the +40 OC fits, the 0 values were set to zero (but see the next section). This simplification undoubtedly distorts the values of the a parameters and could affect the values of VI( TO) and VII(TO)to some extent as well. A fit of the density data to +70 O C certainly requires the inclusion of in eqs 2 and 3. Thus, for fit iv, we fix A , B, C, D,and To,as was done for fit ii, but allow the other six parameters to vary. The Kell data points are again balanced above the below 4 OC, with a total of 39 data points being used. The resulting output parameters from this fit are listed in Table 2. A comparison with the Kell data of the recalculated densities in the temperature range-34 to +IO0 O C from the least-squares fits is provided in Table 3, and Table 4 compares somefi(T) values obtained from the dynamical data with those derived directly from the densities.

VIII. Other Fits Of the more than 60 different fits that were considered, only four have been reported in detail in the tables. A briefer description of some of the other fits will be given in this section. From all the various fits performed, it was clear, in accordance with some of the discussion in the last two sections, that good fits to the density data could be achieved for ranges of the parameters that preserved the general tenor of the proposed mixture model yet did not necessarily point to a mixture of ice-Ih/ice-I1 bonding types.

TABLE 4 Temperature-Dependent Fractions of Ice-I-Type Bonding from Fit i, Compared with Those Used in Fit iii from Temperature-Dependent Dynamic Data t, O C fit i fit iii -34 -30 -20 -10 0 +4 +20 +40 +60 +80

+loo

0.6930 0.6517 0.5757 0.5199 0.4746 0.4585 0.4014 0.3422 0.2925 0.2501 0.2137

0.7039 0.6642 0.5891 0.5336 0.4871 0.4703 0.4104 0.3475 0.2945 0.2494 0.2108

One of these fits fixed PI( TO)a t the ice-Ih value near 225 K, i.e., 0.9237. (From now on the dimensions of the density and volume parameters will be omitted and understood to be as given earlier in the text.) This fit gave good results, but the higher density compared with the fit i value was compensated for by changes in the other parameters, a1 and a11respectively being 1.31 X lo4 and 1.41 X lk3. This a1 value is roughly the same as it is in ice-Ih (Table 3.10 of ref 16), leaving open the question of whether or not, in the liquid where the long-range crystal forces are missing, the Vand a volume parameters are somewhat larger than in ice-Ih. However, the value of TO= 225.815 K was about the same as that obtained for fit i, as was p11( TO)= 1.188. Another fit tried was an algebraic fit of the eight parameters of fit i, but extended to +40 OC. The output parameters in this case were A = 4.11 X le2, B = 2.8 X lV, C = 8.83 X 10-2, TO = 225.34 K, p~(To)= 0.919, p[[(To)= 1.215, 2.61 X 10-4, and a11= 1.SO X 10-3. The only significant difference in these parameters and those obtained from fit i is in the P I [ ( TO)value, which shows a modest increase toward the density if ice-V and -VI as more high-temperature density data are included in the analysis. Fit ii has also been extended to include the 0 parameters, so that compensation among the V, a,and 0 parameters could be evaluated. In this case, the six output parameters have values as follows: PI(TO) = 0.9197, PII(TO) = 1.1953, a1 = 9.34 X 10-4, a11=6.06X 1V,BI=5.34X 10-6,p1~=2.42X10-6. Itisnoticed that the p ~ ( T o )value hardly changed at all, still accurately matching the density found for fits i and ii. The density P I [ ( TO) has moved only very slightly toward ice-V and -VI densities, thus fairly well preserving the ice-Ih/ice-I1 mixture model picture. However, a1 has doubled, while a11has diminished by a factor of 2 to compensate for the presence of the nonzero 0 values. The exact reason for these latter changes is not known but could have to do only with the increased uncertainty in the minor parameters when a greater number of adjustable parameters is employed. One may well ask a t this point how far from the ice-Ih or ice-I1 densities could the related VI(To) and VII(TO)parameters stray and still provide good fits to the density data. The answer to this question is fairly obvious if one thinks in terms of the fit i algebraic procedure. By fixing one parameter, seven are adjustable and thus seven data points can be fit exactly, no matter what value the fixed parameter was given. Furthermore, with this number of parameters, the same number used by Kell,24least-squares fits might also give excellent agreement with the experimental data. So why do we say in the Abstract that the experimental density data provide “strong support” for a mixture model? Certainly the fits presented in Tables 1-3 constitute a necessary condition for the validity of this type of model. The question now being posed is whether those fits provide a sufficient condition. Clearly, VI( TO)and V11(TO)cannot take on absurdly different values than we have derived. Otherwise, nonphysical values of the other parameters would arise: negative a’s, but much more seriously, negative V, negativefi orfrI, or negative TO.To evaluate this point, we carried out four algebraic seven-parameter data fits for the same temperature range -30 to +4 OC,as fit i, where

2228 The Journal of Physical Chemistry, Vol. 98, No. 9, 1994

Vedamuthu et al.

one of the V( To) values was fixed, D = 1, and the /3 values were set equal to zero. The fixed Vvalue in each of these fits was set either 10%higher or 10%lower than the corresponding value in fit i. Briefly, the results in the usual units are as follows: VI(To)[fixed] = 1.20289 ( p = 0.831), Vll(To)= 0.49206 (p

= 2.032), To = 2 19.96 K, aI= 5.56 X 1 O-',

all=

4.84 X lo-', frl(4 "C)= 0.36

VI(To)[fixed] = 0.98418 ( p = 1.016), VIl(To)= 0.93096 ( p = 1.074), other parameters unphysical because of closeness of Vs

VI(To) = 1.05442 ( p = 0.948), VI/(T,)[fixed] = 0.94067 (p

= 1.063), To = 229.67 K, ar = -4.72 X lo4, all=

"c)

8.25 X lo4, frr(4 = 0.58 VI(To)= 1.04572 ( p = 0.956), VI,(To)[fixed] = 0.76964 (p

= 1.299), To = 231.51 K, al = -9.98 X lo4, cyI, = 6.84 X lo-', fr1(4"c)= 0.16

These results are seen not to be very pleasing. The first set corresponds to an unexpectedly high pl1( TO), trying tocompensate for the fixed low pl(7'0). The second set, because of the trend toward similar values of VI(TO)and V11(TO),remained totally unphysical, with negative TO,and fr(7') values always close to zero, which lead to an unstable QII. The third set is not so bad except that a1 is negative, as it also is for the fourth set. In addition to the less than pleasing fitting parameters from these four fits, the densities above the fitting range are not so well represented as for fit i. For example, they are respectively for the three sets of physically acceptable parameters at +40 OC (36 O C above the fitting range) 0.992 13, 0.992 36, and 0.992 08, compared with the fit 1 density at this temperature of 0.992 26 and the Kell density of 0.992 22.

IX. Structural Relationships Looking at a three-dimensional model of ice-11, or comparing illustrations of ice-Ih and ice-II,a which we have reproduced here in Figure 1, one sees that in i c e 4 half the hexagonal channels in the ice-Ih structure have been filled in with collapsed structure having roughly twice the density of the open hexagonal structure. It is as if a new hexagonal tunnel structure has been created in ice-11, butwith thewallsofthetunnelsthickened with thecollapsed structure. This is the origin of the 28%loss of molar volume in the ice-Ih to ice-I1 transformation. An interesting description of the structural relationships between nineof the ice polymorphs, including ice-Ih, -11, -111, -V, and -VI,in terms of "ordering and displacive" mechanisms from a common parent has recently been presented in Physical Review Letters by Dmitriev, Rochal, and Toledano.4' One can, in fact, picture a fluctuational network of ice-I-type bonding interlaced with ice-11-type bonding (plus ice-111, -V, or -VI types at elevated temperatures or near certain interfaces). In the denser forms of bonding, probably some stabilizationenergy is gained because of the near-favorable O-*Ovan der Waals distance of -3.3 A in these structures. Of course, in the liquid state, compared with the ice-I1 crystal, the free energy of component I1 is probably lowered, since the H atoms making up the hydrogen bonds are very likely more disordered as far as to which of the two 0 atoms they belong. It is to be remembered that the entropy of proton disorder"'2in the iceIh crystal provides an important contribution to that polymorph's stability. (Compare the AG values from the

Figure 1. Illustrations of (a, top) the ice-Ih crystal structure looking down the c axis, clearly showing the hexagonal tunnels, and (b, bottom) the i c e 1 1 crystal structure along its hexagonal c axis, showing that half the tunnels in ice-Ih have been replaced with condensed structure, which surrounds each of the ice-I1 tunnels. Reproduced with permission from ref 40. Copyright 1964 W. H. Freeman and Co.

data provided in Table 3.6 of ref 16 with and without the entropy of proton disorder.) Near the lowest experimentally available temperature of supercooling, -34 OC (the homogeneousnucleation temperature), Table 4 indicates that the fractionj; of ice-I-type bonding from the fit i analysis is -0.69. At +4 "C this fraction drops to -0.46. The intermixture of thesevery large proportions of material having such great structural and density disparity should surely be detectable experimentally. In the next section, one such experiment on the differential X-ray scattering of the liq~id,3~ will be

Feature Article seen to lend important support to the presence of two temperaturedependent bonding forms, of the types expected from the density analysis. Is there any additional evidence for the discrete regions of differential density in liquid water proposed by the model? Two things7 may make experiments difficult: (1) rapid fluctuations involving the two bonding forms would tend to give an average for many types of experimental data; (2) the spatial extent of the density variations is not very large, generally not much greater than about twice the molecular diameter of a water molecule, with the collapsed structure, say of ice-11, perhaps forming highdensity striations through the open structure, as can be imagined from Figure la,b. The non-nearest-neighbor aspect of the structural reorganization in liquid water ruins usual intuition about chemical or intermolecular bonding. It also raises serious questions about the realism of all prior computer molecular dynamics (MD) simulations of liquid water,43 particularly under variations of temperature and pressure. The “best” simulations are those that use highly complex intermolecular potential functions, inclusive of intramolecular bond flexibility and rapidly responsive “electronic polarization”,u but for just this reason they are the ones limited to a small sample size, usually no more than 256 molecules. Can expected non-nearest-neighbor effects be well represented in such a small sample volume? It would seem that any MD model of the liquid would have first of all to reproduce well the ice polymorph structures under appropriate thermodynamic conditions. Seldom4s has this been a goal of such calculations. Other concerns with this type of structure in the liquid would be the loss of simple connections between water clusters and the liquid state. Considering the subtle non-nearest-neighbor effects expected from the liquid model proposed here, it would be most difficult to extract much useful information about liquid water from cluster studies, particularly a t the dimer or trimer le~e1.~6

X. Independent Support for the Model The most direct independent support for the model follows from the isochoric differential X-ray diffraction data of Bosio, Chen, and T e i ~ e i r a . 3Clearly ~ seen from these data is a peak at about 3.3 A, which grows in with increasing temperature differential. This distance, as mentioned earlier, is one that does not occur at all in ice-I structures but is very close to prominent features in ice-11, -111, -V, and -VI.’ At the same time, in the liquid, a robust decrease in the 4.5-A feature, a non-nearestneighbor O-.Odistance not occurring in the dense ice forms,7 marks the disappearance of ice-I-type structure with increasing AT. The data of Bosio et al. also show a concomitant, though relatively smaller, decrease in the 2.8-A peak over their ATrange. This is close to the hydrogen-bonded nearest-neighbor 0-0 distance that occurs in all these ice f o r m ~ , ice-Ih, ~ J ~ -11, -111, -V, and -VI. The decrease of this peak is very likely a direct measure of molecules that are in some sort of a mixed disordered structure or an intermediate rotational configuration between the structures. The results of Bosio et al. then seem a clear indication that a dynamic, temperature-dependent mixture of ice-I-, -11-, -111-, -V-, and/or -VI-type bonding is present in the liquid in the manner expected for the mixture model described in the present paper. Kamb (see fig 10 of ref 7) in his paper also pointed out the relationship between the radial distribution function of the liquid and one created by a superposition of ice-Ih, -11, and -111 radial distribution functions. What about the internal vibrational frequencies? It is known that the Raman spectrum of liquid water can be approximately resolved into two component~,g+47~~ one of which has been termed “openn or “hydrogen-bonded” and the other as “closed“ or “nonhydrogen-bonded”. In fact, for temperatures above 0 “C, a clear isosbestic point has been fo~nd,*,4~ possibly indicating the presence of two distinct species. In the low-temperature work, particularly that carried out on 0-D vibrations in D20 or HOD in H20,49

The Journal of Physical Chemistry, Vol. 98, No. 9, 1994 2229 it is also clear that the lower frequency component in the supercooled liquid tends toward the ice-Ih or amorphous ice-I Raman spectrum as the temperature is reduced. These aspects seem to support the mixture model proposed here. However, it remains to decipher the entire temperature-dependent infrared and Raman spectrum of liquid water at low temperatures in terms of the ice-Ih/ice-11 model. This chore is made difficult by the probable presence of temperature-weighted ice-I-type and ice11-typeintra- and intermolecular modes, including phonon modes, whose temperature-dependent frequencies are coupled together,49 with the coupling parameters themselves probably being temperature dependent. Whether or not the model is consistent with the very recent low-angle x-ray scattering data35remains to be seen. One aspect of our present density studies does seem worth mentioning in this connection. At the beginning of our analysis, a variable power law parameter n was incorporated into the fitting functionfir of eq 4 so that this function would approach zero as ( T - To)”. However, preliminary fits immediately showed that n was extremely close to unity, so, subsequently, n was set to one. Such a power law dependence shows that the functionfil is not singular; Le., none of its derivatives diverge as T approaches TO.Correspondingly, the density and other related thermodynamic properties should exhibit analytic behavior in the vicinity of TO.This was at first a surprise to us, since other published work15123 indicated otherwise. However, ref 35, finding no strong increase in correlation length in the vicinity of TO,also concludes that there is no singular behavior near this temperature. It should also be remarked here that the cooperative melting phenomenon2 discussed in the earlier work1 does not demand a singular function for the fraction of “rotating” molecules. The important feature for the water problem is the effect that these rotating molecules have on softening and flattening the intermolecular potential functions, somewhat reminiscent of phonon “mode softening” in solid-state physics.39,So

XI. Concluding Remarks It should be kept in mind that the main purpose of the present paper was not just to fit density data accurately. This has already been done by K e l F and others. Rather, the purpose was to check the consistency of the derived mixture model parameters with solid-state bonding parameters, and eventually with other types of experimental results, and thus to create an alternate way of thinking about liquid water when future experiments are being designed and analyzed. The exact agreement in Table 1 between calculated mixturemodel densities using fit i parameters and actual H2O densities in the temperature range -34O to nearly +40 OC, in addition to the fact that the output parameters in Table 2 have close relationships with bonding characteristics of the open and dense polymorphs of ice, indicates that a necessary condition for the validity of the mixture model has been achieved. The near uniqueness of the fitting parameters, discovered from the examination of a large number of fits, goes partway in providing a sufficient condition for such validity. We hope in future work, through consistency checks for other types of experimental datal96951 on liquid water, to create further convincing arguments supporting the validity of the model. The somewhat weaker ability of the equations to reproduce the density data into the higher temperature regimes is not caused, as can be inferred from Table 4, by the form offII, but rather is caused by the expected breakdown of eqs 1-3 a t high temperatures, where a plethora of bonding types can occur. This bonding multiplicity has been one of the main objections in the past to “two-state models” and “mixture models” of liquid water. Because of the results presented in this paper, it is now felt that such objections lose much of their validity at sufficiently low temperatures. Water, after all, is a completely different type of liquid than liquid argon, on which most theories of the liquid

2230 The Journal of Physical Chemistry, Vol. 98, No. 9, 1994

state have been based and for which the objections would have more validity. Therefore, the mixture model concept, which is certainly simpler than other proposed 'liquid theory" modds, should be able to form a sound low-temperature base for extending theunderstanding of liquid water to other thermodynamic regimes. For most substances, the density, though a necessary and useful concept, is a barren, even pedestrian, quantity. Yet, for liquid water, the density has long been s~spected59~.5*,53 to hold the key to this substance's secrets. In fact, the year 1991 marked the 100th anniversary of such ideas.$ Thanks to modern experimental technology, the visions of these past workers may be coming close at last to being vindicated. Acknowledgments. Discussions with Walter Kauzmann and George Walrafen helped the authors formulate a better approach to the mixture model concept for liquid water than in the earlier paper.' A very fruitful telephone conservation with C. M. Sorensen, as well as his sending us a preprint of ref 35, are gratefully acknowledged. A number of discussions with Barclay Kamb were extremely helpful, in particular his point by point criticisms of an earlier version of the paper. It was he who pointed out to us the existence of the drawing of ice-I1 in ref 40, which we have reproduced here with thekind permission of LinusPauling and W. H. Freeman and Co. Financial support has been shared by the National Science Foundation (CHE-9112002), the Robert A. Welch Foundation (D-0005and D-1094), and the state of Texas Advanced Research Projects program (003644-004).

Appendix: Least-Squares Procedure Here, we give a brief summary of our least-squares fitting methods. The procedure employed is the nonlinear technique of Levenberg and Marquardt as described in the book by Bevington.54 The Fortran routines used are taken from Press et aLS5 Most of the routines were run in quadruple precision, so as not to lose important 5th to 7th decimal place digits in the final densities because of round-off a t each step in the iterative procedure. The need for this increased precision lengthened the computational times, but these were not excessively long (typically no more than an hour) even on our rather ancient in-house VAX 11/730. The Levenberg-Marquardt fitting routine was found to be remarkably stable and reliable even for most of the six- to eightparameter fits applied to the full reported Kell density values (6-7 decimal places). Even starting 'blindly" with initial values far from the final values, the routine, as monitored by the running x2 values, converged reasonably rapidly.

References and Notes (1) Bassez, M.-P.; Lee, J.; Robinson, G. W. J. Phys. Chem. 1987, 91, 5818. (2) Fowler, R. H. Statistical Mechanics, 2nd ed.; Cambridge University Press: Cambridge, 1966; pp 810-816. (3) Kauzmann, W. L'Eau Syst. Biol., Colloq. Int. C.N.R.S. 1975,246, 63. (4) See, for example: Endo, H. J . Chem. Phys. 1980, 72,4324. (5) Rdntgen, W. C. Ann. Phys. Chem. 1892, 45, 91. (6) Robinson, G. W.; Vedamuthu, M. Unpublished results. (7) Kamb, B. InStruciural Chemistry andMolecular Biology; Rich, A., Davidson, N., Eds.; W. H. Freeman and Co.: San Francisco, 1968; pp 507542. (8) Larsson, K. E.;Dahlborg, U. React. Sci. Technol. 1962, 16, 81. (9) Walrafen, G. E.; Fisher, M. R.; Hokmabadi, M. S.;Yang, W.-H. J . Chem. Phys. 1986, 85, 6970. (10) Chen, S.-H.; Toukan, K.; Loong, C. K.; Price, D. L.; Teixeira, J. Phys. Rev. Lett. 1984, 53, 1360. (11) SegrC, E. Rend. Lincei 1931, 13, 929. (12) Pauling, L. Phys. Rev. 1930, 36, 430. (13) Fowler, R. H. Proc. R . SOC.London 1935, A149, 1. (14) Robinson, G. W.; Zhu, S.-B. In Reaction Dynamics in Clusters and Condensed Phases; Jortner, J., Levine, R. D., Pullman, B., Eds.; Kluwer Academic Publishers: Dordrecht, 1994; pp 423-440. (15) Angell, C. A. Annu. Rev. Phys. Chem. 1983,34, 593. (16) Eisenberg,D.; Kauzmann, W. TheStructureandProperties of Water; Oxford University Press: New York, 1969. (17) Lang, E. W.; Liidemann, H.-D. Ber. Burrsen-Ges. Phys. Chem. 1981, 85, 603. (18) Tables 4.5 and 4.8 and accompanying discussion in ref 16. (19) Pruppacher, H. R. J. Chem. Phys. 1972,56, 101.

Vedamuthu et al. (20) Bertolini. D.; Cassettari, M.; Salvetti, G. J. Chem. Phys. 1982, 76, 3285. (21) Kanno, H.; Speedy, R. J.; Angell, C. A. Science 1975, 189, 880. (22) See ref 1, section 5. (23) Speedy, R. J.; Angell, C. A. J. Chem. Phys. 1976,65,851. (24) Kell, G. S. 1. Chem. Eng. Dara 1975,20,97; in particular, see Table 111 of this reference. (25) Baker, G. A., Jr. Essentials ofpad6 Approximants; Academic Press: New York, 1975. (26) Schufle, J. A. Chem. Ind. 1965, 16, 690. (27) Zheleznyi, 8. V. Russ. J. Phys. Chem. (Engl. Transl.) 1969, 43, 1311. (28) Hare, D. E.; Sorensen, C. M. J. Chem. Phys. 1987,87,4840. (29) Schufle, J. A.; Venugopalan, M. J . Geophys. Res. 1967, 72, 3271. (30) Leyendekkers, J. V.; Hunter, R. J. J. Chem. Phys. 1985,82, 1440; 1447. (31) Hare, D. E.; Sorensen, C. M. J. Chem. Phys. 1986, 84, 5085. (32) Mohler, J. F. Phys. Rev. 1912, 35, 236. (33) Sceats, M. G.; Rice, S. A. In Water-A Comprehensive Treatise; Franks, F., Ed.; Plenum: New York, 1982; Vol. 7. (34) Bosio, L.; Chen, S.-H.; Teixeira, J. Phys. Rev. 1983, A27, 1468. Chen, S.-H.; Teixeira, J. Adv. Chem. Phys. 1986,64, 1. See Figure 8 of this paper. (35) Xie, Y.; Ludwig, K. F., Jr.; Morales, G.; Hare, D. E.; Sorenscn, C. M. Phys. Rev. Lett. 1993, 71, 2050. (36) Lang, E.; Liidemann, H.-D. J. Chem. Phys. 1977,67, 718. (37) Lang, E.; Liidemann, H.-D. Ber. Bunsen-Ges. Phys. Chem. 1980.84, 462. (38) Lang, E. W.; Girlich, D.; Lildemann, H.-D.; Piculell, L.; Mfiller, D. J. Chem. Phys. 1990,93,4796. Lang,E. W.;Liidemann,H.-D. N M R Basic Princ. Prog. 1990, 24, 129. (39) A rough justification for this form can be found in: Kittel, C. Introduction to Solid State Physics, 5th ed.; John Wiley & Sons: New York, 1976;p 422, noting that for a sinusoidal rotational potential function, the peak is proportional to the square of the frequency. to valley energy (m) (40) Pauling, L.; Hayward, R. The Architecture of Molecules; W. H. Freeman and Co.: San Francisco, 1964; plates 41 and 42. (41) Dmitriev, V. P.; Rochal, S. B.; Toledano, P. Phys. Rev. Lett. 1993, 71, 553. (42) Pauling, L. J. Am. Chem. SOC.1935,57, 2680. (43) For a recent comprehensive overview of computational studies of liquid water, see: Zhu, S.-B.; Singh, Surjit; Robinson, G. W. Adu. Chem. PhyS. 1993,85 (3). 627-731. (44) Zhu, S.-B.; Robinson, G. W. J. Chem. Phys. 1991, 95, 2791. (45) As examples, see: Morse, M. D.; Rice, S.A. J . Chem. Phys. 1982, 76, 650. Kuwajima, S.;Warshel, A. J. Phys. Chem. 1990, 94, 460. (46) Saykally, R. J. Multidimensional Tunneling Dynamics in Hydrogen Bonded Clusiers: The Ammonia Dimer and Water Trimer;paper PHYS 94; 206th American Chemical Society National Meeting, Chicago, IL, Aug 2227, 1993. (47) d'hrigo, G.; Maisano, G.; Mallamace, F.; Migliardo, P.; Wanderlingh, F. J. Chem. Phys. 1981, 75,4264. (48) Hare, D. E.; Sorensen, C. M. J. Chem. Phys. 1990, 93, 25. (49) Hare, D. E.; Sorensen, C. M. J. Chem. Phys. 1990,93,6954; Chem. Phys. Lett. 1992, 190, 605. (50) Perry, C . H.; Agrawal, D. K. Solid State Commun. 1970, 8, 225. (51) Forexample,insomeoftheunpublished workcitedin ref6,preliminary fits of DzO density data, of the same types as fit i for H20, yielded p~(To)[fixed] = 1.020 g cm-3, pn(T0) = 1.301 g cm4, a1 = 4.01 X l W K-I, a11= 1.30X lo-3K-1, T~=231.83K,andfi(11.18)=0.46,where11.18°Cisthe temperature of maximum density in liquid DzO. Note the agreement between thefivalues for the twoisotopicmodificationsattheir temperaturesofmaximum density and the near agreement in the I- and 11-component density ratios (1.109 and 1.101) with the D20/H20 ratio of 1.I 102 0.0002 in ice-Ih. See: Lonsdale, K.Proc. R . Soc. London 1958, A247,424. (52) Tammann, G. Z . Anorg. Chem. 1926, 158, 1. (53) Be", J. D.; Fowler, R. H. J. Chem. Phys. 1933, I , 515. (54) Bevington, P.R. Dara Reduction and Error Analysisfor the Physical Sciences; McGraw-Hill: New York. 1969; Chapter 11. (55) Press, W. H.; Flanncry, B. P.; Teukolsky, S. A.; Vetterling, W. T. Numerical Recipes in C: The Art of Scientific Computing; Cambridge University Press: Cambridge, 1988; Chapter 14.

*