1554
J. Phys. Chem. 1995, 99, 1554-1565
Computer Simulation Study of an sN2 Reaction in Supercritical Water Perla B. Balbuena,*J Keith P. Johnston,*jt and Peter J. Rossky*9$ Department of Chemical Engineering and Department of Chemistry and Biochemistry, The University of Texas at Austin, Austin, Texas 78712 Received: June 10, 1994; In Final Form: November 10, 1994@
Molecular dynamics computer simulation is used to examine solvation of species of varying polarity along the reaction coordinate for the S N reaction ~ of C1- and CH3C1 in supercritical water. At reduced densities down to 0.5, the solute-solvent interactions are sufficiently strong to preserve solvent configurations with solute coordination numbers nearly as high as in ambient water. Hydrogen-bonding interactions are converted to less specific polar interactions. This loss of hydrogen bonds from ambient water to supercritical water decreases as the hydrogen bond strength increases for the series of solutes from the equivalent Cl's in the transition state complex, to C1- in the ion dipole complex, to free C1- in the reactant state. Two-dimensional angle averaged cylindrical and (1-D) spherically averaged pair distribution functions, along with energy distribution functions, provide a clear explanation of the origins of AA, AE,and TAS along the reaction coordinate.
Introduction Supercritical water (SCW) exhibits a unique blend of nonaqueous as well as aqueous characteristics. At the critical point, (Tc = 374 "C; P, = 221 bar; ec = 0.322 g/cm3) the volume is 3 times larger than at ambient conditions, and the dielectric constant E is only 5.3. As a result, SCW can dissolve most organics and gases such as oxygen; thus it is a promising solvent for oxidation of organics in waste Furthermore, nearcritical water and SCW can replace certain organic solvents for chemical synthesis. To develop these applications, a better understanding is needed of the molecular interactions in SCW solutions and how they influence chemical reactions. For example, as the density increases from gaslike to liquidlike values, the dominant reaction chemistry can change from free radical (homolytic) to ionic (heterolytic) mechanism^.^,^ As in all supercritical fluids, pressure and temperature have a large effect on the density and 6 , which can lead to pronounced values for partial molar volumes and enthalpies of solute^.'.^^^ Likewise, temperature and pressure can have a large influence on dissociation reactions of inorganic acids, bases, and salts* and on redox reaction^.^ Recently, the equilibrium constant for an organic reaction, the proton transfer between 2-naphthol and OH- has been measured spectroscopically in SCW up to 400 "C and 7000 psia.1° As the temperature increases from 25 to 400 "C, the behavior changes from exothermic to endothermic due to the larger electrostriction of water about OH- than the naphtholate anion. Whereas a number of studies have reported macroscopic thermodynamic properties of SCW solutions, experimental studies focusing on molecular properties are relatively rare. Hydrogen bonding for dilute HDO in H20 has been characterized qualitatively in terms of frequency shifts for OD vibrations." The radial distribution function for SCW has been measured by neutron scattering.I2 Ion pairing has been characterized in nitrate solutions with Raman ~pectroscopy'~ and in 2-naphtholate solutions with UV-vis spectroscopy.1° Furthermore, solvatochromic measurements indicate that the local densities of water molecules defined in the first coordination shell about benzophenone and acetone are augmented above Department of Chemical Engineering.
* Department of Chemistry and Biochemistry.
@Abstractpublished in Advance ACS Abstracts, January 1, 1995.
0022-365419512099- 1554$09.OOlO
the bulk densities.14 The so-called solute-solvent clustering has been measured in other supercritical fluids such as C0215 and has been described by simu1ation.l6 In both CO2I7 and water,14 solute-solvent interactions have been characterized in three density regions: gas-phase solute-solvent clustering, clustering in the supercritical region, and "liquid-like" solvation. Clustering usually persists as the density is decreased from well above ecto well below the critical density and for some systems is prevalent even down to er = @lec= 0.1.16-zo However, clustering is a short-range phenomenon that need not track the bulk thermodynamic solvent compressibility, which diverges at the critical point. Solvent-solute clustering has been included in thermodynamic models for the free energy of an ion in S C W . ~2122 ~, Given the lack of measurements of molecular interactions in supercritical water, due to significant experimental challenges at these severe conditions, computer simulation is an especially important tool. Gibbs ensemble Monte Carlo simulation23with the simple point charge (SPC) potentialz4predicts Tc = 587 K, and ec = 0.270 glcm3,which is reasonably accurate considering that the model parameters were optimized at ambient conditions. The degree of hydrogen bonding in SCW solutions is somewhat controversial. While the energy of a hydrogen bond in pure waterz5 is 10 kT at 298 K ( k is the Boltzmann constant), it is still about 4 kT at the critical temperature. According to computer simulation with the TIP4P model, hydrogen bonding persists up to the critical temperature for densities as low as 0.1 m g h L Z 6 A similar conclusion was reached from spectroscopic measurements for hydrogen bonding between water and acetone.l4 Furthermore, good agreement was observed between the simulation data for pure water and predictions of the lattice fluid hydrogen bonding However, hydrogen bonding was found to be considerably weaker in the above neutron scattering studyI2 than in these simulations. Clearly, further study is needed to understand hydrogen bond donor and acceptor strengths in SCW solutions. Recently, the structure of SCW about various ionic and nonionic solutes has begun to be characterized with MD computer ~ i m u l a t i o n . ~ *For - ~ ~sodium and chloride ions, the first solvation shell in SCW is observed to be similar to that in ambient water, indicating significant clustering of solvent about the solute, above the bulk water density.28 For the interaction 0 1995 American Chemical Society
sN2 Reaction in Supercritical Water
J. Phys. Chem., Vol. 99, No. 5, 1995 1555
between the sodium and chloride ions,3oentropic forces drive ion association which dominate over enthalpic forces that favor dissociation. The opposite occurs in ambient solution. For a chemical reaction, computer simulation can provide a detailed picture of solvent configuration along the reaction coordinate. Although the influence of solvent-solute clustering has been invoked to explain experimentally measured reaction rate constants, only approximate information about the structure of the solvent has been a ~ a i l a b l e ,from ~ ~ . ~either ~ spectroscopy or molecular thermodynamic models.9b The same is true for clustering involving a hydrogen bond donor and acceptor in an inert ~ o l v e n t or ~ ~for , ~ two ~ species undergoing a chemical r e a c t i ~ n . ~ Therefore, ~ , ~ ~ . ~ ~computer simulation studies can provide great insight into the understanding of how the solvation structure and forces influence the free energy barrier for a chemical reaction. In a previous c o m m u n i ~ a t i o nwe ~ ~ presented the first molecular simulation studies of the effect of a supercritical fluid on the rate of a chemical reaction, in particular SCW. The free energy surface for the sN2 reaction of c1- and CH3C1, and its energetic and entropic contributions, were compared with those in ambient water and in an organic solvent, dimethylformamide (DMF). For el = e /ec= 1.5, and Tr = T/T, = 1, the behavior in SCW lies between that in ambient water and DMF, despite the low dielectric constant of 9.75.40We concluded that the details of solvation on a molecular scale play an important role for the reaction at these conditions. The present objective is to relate the free energy surface to the microscopic structure based on comprehensive two-dimensional angle averaged cylindrical and (one dimensional) spherically averaged distribution functions, as well as energy distributions, coordination numbers, numbers of hydrogen bonds, and dielectric properties. The structure of water molecules about the C1 nucleophile, C1 leaving group, and CH3 group will be examined in detail along the reaction coordinate. Computer simulation offers a new opportunity to discover how hydrogen bonding and clustering vary as a function of the charges on the various groups, as these charges change from the reactant state to the ion-dipole complex to the transition state. Furthermore, the full three dimensional spatial variation in microscopic solvation about the reacting species may be resolved for the first time, providing significant advancement in the understanding of solvent effects on reactions in supercritical fluids. Method and Model
When a chemical reaction takes place in the gas phase, the free energy barrier for the reaction is completely determined by the interactions among the reactants. In a condensed phase, it is also influenced by the forces exerted by the solvent. The total change in free energy must include the change in the free energy of the solvent as well. The resulting free energy for the reaction as a function of the reaction coordinate gives a free energy surface also known as the potential of mean force (PMF). The reaction we have chosen to study, C1- CH3Cl- ClCH, C1-, is a prototype sN2 substitution reaction. The free energetics have been investigated in water at ambient conditions by Monte Carlo simulation41and by integral equation theories.42 The gas-phase reaction profile has been characterized in previous ab-initio calculation^^^^^^ which give the geometries, charge distributions, and energies along the reaction path. These results also provide the parameters for the intermolecular potential functions needed for solution-phase simulations as fully described e l s e ~ h e r e . ~The ~ . ~solute-water ~ and water-water interactions are assumed to be painvise additive and described through Coulomb and Lennard-Jones terms, acting between sites
+
+
TABLE 1: Electrostatic Charges of Solute Atoms (Unitsof e)
structure reactants
QCH?
9Cl’
0.178
IDC
0.220
TSC
0.492
-0.178 -0.247 -0.746
9c1
-1.0
-0.973 -0.746
on the monomers. The methyl group is represented as a single site. The SPC was used for water. An asymmetric reaction coordinate was defined as rA = rccr - reel-, the difference between the distances from the methyl carbon to each of the chlorine atoms (denoted C1 and Cl’). This choice of reaction coordinate reflects the symmetry of the reaction. Therefore, lrAl ranges from infinity (separated reactants or products) to zero (transition state complex). All parameters vary smoothly along the reaction coordinate through a functional form that was fitted to ab-initio results (see Table 3, ref 41, and Table 2, ref 42). Statistical mechanical perturbation theory43is used in order to obtain the incremental change in free energy in the solution at each step along the reaction coordinate:
Here dA is the change in free energy when the system is perturbed from rA to rA -k drA, Uws is the water-solute configurational energy (which depends parametrically on rA) and the brackets indicate an average taken for the solution at fixed rA. A critical feature of this system in the present context is the redistribution of charge during the reaction, with the net negative charge localized on one C1 in the reactant state and distributed over the complex in the transition state. The atomic charge parameters of the model solute for the three geometries of special interest are summarized in Table 1. The complete model is specified in earlier For each point along the reaction coordinate, molecular dynamics simulations were performed for about 20 ps to equilibrate the system. Then the system was perturbed (for 1525 ps) by drA = 0.125 8, to obtain the average in eq 1. The time step was in all cases 1 fs. A total of 484 water molecules was used to solvate this reaction. The size of the simulation box was L x L x 1.55. At each state point, the bulk density was obtained by changing the volume of the box, keeping constant the number of water molecules. Thus, L was the appropriate value in angstroms to obtain the desired bulk density, for example, L = 28.84 A for Tr = 1 and er = As is shown later when we discuss the variation of the local densities for different conditions of bulk density, a proper check was done in order to assure that the correct bulk density was obtained in the regions far from the solute. Periodic boundary conditions were applied, and the forces were calculated using the minimum image convention. For each value of rA, the solutes were constrained to be collinear with the longest axis of the simulation box. The PMF was calculated as
where AA(rA) is the change in free energy in the solution at rA, which results from the sum of all the contributions dA(rA) when the solutes are brought together in small steps drA from infinite separation to the distance rA. The second term, Aug,(rA), is the corresponding change in the gas phase which is additive in the model!* In our calculations, we have specified the reference value W (IrAl = 8 A) = 0 for a given T and 4. That this choice is adequate for present purposes is justified on the basis that
Balbuena et al.
1556 J. Pkys. Ckem., Vol. 99,No. 5, 1995
In AW, the minima corresponding to an ion-dipole complex essentially disappear, and the activation barrier grows to 26.3 kcdmol as given by Monte Carlo ~imulation.~'This well disappears since the energy of ion-dipole attraction is offset by the net weakening of the nucleophile-water hydrogen bond interaction!* The free energy profile obtained in SCW39closely resembles that in ambient water, despite the reduction in solvent density of about two, and in the dielectric constant E from 78 to 9.75. However, a slight minimum is resent corresponding to the ion-dipole complex (IrAl * 1.5 ). As shown below, there are fewer hydrogen bonds between the nucleophile and water in SCW compared with ambient water. Thus, the smaller degree of desolvation (reduction in coordination number) in SCW does not completely compensate for the stabilization upon ion-dipole complex formation. In SCW, the activation barrier between the ion-dipole complex and the transition state is calculated to be 25.3 kcallmol, 1 kcdmol higher than that between the reactants and the transition state. Consequently, there is an equilibrium between the reactants and the IDC. To calibrate the results against those in a polar organic solvent, we consider dimethylformamide ( E = 37). The activation barrier according to Monte Carlo simulation4*is 19.3 kcdmol, and the well for the ion-dipole complex is -2.7 k c d mol. For each of these properties, the results in SCW lie between those in ambient water and DMF, even though E is only 9.75.40 Three-Dimensional Solvation Analysis. In order to understand the role of microscopic solvation (defined by the molecular solvent configuration near the solute) in SCW, we present the structural properties for er = 1.5 and Tr = 1. These results are for three values of the reaction coordinate: the separated reactants (IrAl = 8 A), the IDC ( 1 1 ~ 1= 1.441 A), and the TSC ( r A = 0). The distribution of electrostatic charges on each solute atom is given in Table 1. The reactants are located along the long axis of the simulation cell, with cylindrical symmetry (for the model in which CH3 is represented as an extended atom) about this axis. Therefore, we can calculate angle averaged cylindrical resolved distribution functions for the water molecules near the solutes at the above three geometries along the reaction coordinate. Plots of these solute-water functions gso( r ~are ) shown in Figures 1-3 , where r and z are the radial and axial cylindrical coordinates, respectively. Here g s o ( r , z ) is defined as the relative local density of oxygen (water) atoms with coordinates I and z with respect to the solute, compared to the bulk density. The variation of the mean local density of water is represented by the color scale, increasing from violet to red according to the colors of the rainbow. In each of these figures, only water density is represented; the locations of the solutes are given by "holes", that is, regions of zero water density. Clearly, the figures demonstrate that solvation is dominated by electrostatic forces in this case. For example, there is a pronounced accumulation of water molecules near the chloride ion (Figure 1) in the case of the separated reactants or products and only a slight increase in water about the dipole. As the ion approaches the dipole to form the IDC complex, water is squeezed out from between the two. Again (Figure 2), there is a clear preference for the water molecules to solvate the ion Results versus the dipole. The water concentrations on the left and right Potential of Mean Force. In our previous comm~nication,~~ sides of the ionic end (top in Figure 2) of the complex increase. The peaks observed on the axis of the reaction coordinate ( I = we have reported the PMF as a function of the reaction 0) are due to statistical uncertainties, and they should decrease coordinate in SCW and the comparison with published data41*42 in magnitude for longer simulation times. The IDC is as well for ambient water (AW) and the gas phase. In the gas solvated as the reactants, if not more so, as will also be seen in the double-well energy surface indicates the formation of an ion-dipole complex and a symmetric transition state complex. the next section. In contrast, the TSC (Figure 3) is highly
the value at ambient conditions is less than 1 kcal/mol!* The complete free energy calculation required about 55 h on a Cray Y-MP at a single T and e. The energetic contribution A E ( r A ) to the free energy was calculated from the simulations by a method based on finite difference derivatives with respect to temperaturea for all values of r A ; except at the transition state, r A = 0 where it was obtained directly from a 150 ps simulation. The entropy change, - T A s ( r A ) is given by AA - AE. For accurate configurational analysis, specific fixed geometries were simulated for 150 ps, which takes about 9 h on a Cray Y-MP. The solutes were studied at fixed geometries corresponding to the separated reactants, at lrAl = 8 A; iondipole complex (IDC), at ( I A ~= 1.441 A, and the transition state complex (TSC) at lrAl = 0, at TI = 1.0 and eI = 1.5. In addition, we have analyzed the geometry of the IDC for the following state points: er = 0.5, 0.9, 1.5, and 2.0 at TI =1.0 and eI = 1.5 at Tr = 1.1. Some limitations of the present approach are inherent to the use of semiempirical Hamiltonians. Ab-initio molecular dynamics calculation^^^ have recently been applied to SCW. However, even if accurate, the large amount of computer time required by these techniques for a single simulation makes them unfeasible for free energy calculations with the present computational availability. Within the context of classical molecular dynamics, the applications of free energy calculations to several chemical and biochemical systems have been reviewed re~ e n t l y .It~ ~is well-known that perturbations dominated by electrostatic forces (the case of this reaction) converge more quickly than van der Waals perturbations. Standard tests46have been performed to check the accuracy of the results at each value of the reaction coordinate. The average statistical error for the incremental free energy in our simulations was 0.1 kcall mol. The uncertainties are larger for AE and TAS, as discussed p r e v i ~ u s l y .In ~ ~our case, one point of concern is the absence of counterions in the solution. In ambient conditions, the activation barrier calculated bv molecular simulations with the same approximation, for this sN2 reaction4' agrees well with experimental results. However, under supercritical conditions, the presence of counterions leads to ion pairing at lower densities, thus making the reaction less favorable. Work is in progress 'to evaluate this effect in our calculations. Another important point is the use of a rigid model for water at elevated temperatures, with polarization effects treated in a mean field way by using effective charges on the atoms. In a recent study, Zhu and Wong41 have presented a polarizable model for water and analyzed the sensitivities of the structural and thermodynamic properties to variations of the parameters in the intermolecular potentials. Among other findings, they concluded that the use of a point-charge model (as opposed to dispersed charge) is in fact a good approximation to describe the electrostatic interactions among water molecules. We note also that the comparisons to AW presented here are based on a different model (TP4P) for water!l However, considering that these rigid point-charge models are similar and both well tested, we assume in this comparison that the effect on the free energy calculations is small.
-
1
sN2 Reaction in Supercritical Water
J. Phys. Chem., Vol. 99, No. 5, 1995 1557
r Figure 1. Cylindrical 2-D oxygen-solute pair distribution function gso(r,z)corresponding to the separated reactants at Qr = 1.5 and Tr = 1.O. The chloride nucleophile approaches the CHjCl from the top of the figure. The range is from 5.2 to < 1.3.
desolvated, due to the dispersal of charge (see Table l), with relatively little augmentation in water density over the bulk. In particular, the water concentrations at the sides of the complex are much lower than about the Cl- nucleophile in the reactant state. Spherical Pair Distribution Functions at Qr = 1.5 and Tr = 1. The character of the figuresjust discussed provides spatial information about the distribution of the molecules within each shell that is not shown in the usual fully angle-averaged representation of the spherical pair distribution functions. However, we have calculated these radial distribution functions (rdf s), denoted g(r), for several reasons. First, these rdf s may be compared with published values in ambient Also, the r d f s are used to determine the orientation of the water molecules in the first solvation shell of the solutes, by examining both the solute-oxygen and solute-hydrogen rdf s. Finally, the g(r)’s are somewhat quantitative, i.e., the value of r for each peak may be easily identified. In the following, the functions go(r) are the conventional correlation functions, giving the relative local density of atomic species j at a radial distance r from an atom of species i,
compared to the bulk density. Figure 4 shows the pair function gc~o(r) where C1 is the chloride nucleophile and 0 is the oxygen in the water molecule. The expected decrease in the number of molecules in the first solvation shell of the chloride nucleophile is clearly seen in going from the reactants to the TSC, due to the loss in charge. The first peak is centered at 3.4 8, for the reactants and the IDC and at 3.7 8, for the TSC. This peak corresponds to a shell of water molecules, in which many of the water molecules are hydrogen bonded to the chlorine. Below we will find the actual number of water molecules that are hydrogen bonded to C1- from the chlorinehydrogen distribution functions and from the chlorine-water pair energy distributions. In Figure 4, a second solvation shell is weak but clearly seen in the case of the reactants, and it is still present but less clearly defined in the IDC. It is replaced by a much broader shell about 5 8, wide for the TSC. From the corresponding distribution functions for ambient water:’ one finds that the positions of the first peaks are coincident for AW and SCW. Other aspects of the distribution of water molecules are different. For the TSC, the relative heights for the first peak are slightly higher in AW. Due to the
Balbuena et al.
1558 J. Phys. Chem., Vol. 99, No. 5, 1995
“I
r
Figure 2. gso(r,z)corresponding to the ion-dipole complex at Qr = 1.5 and T,= 1.0.The ionic part of the complex is located toward the top of the figure.
ratio in bulk densities of 0.45,
(3) According to this equation, Figure 4 indicates that the number of water molecules surrounding the chloride ion at the distance corresponding to the maximum of the first peak is about 2 times larger in AW4’ than in SCW at these conditions for the case of the reactants and about 2.7 times larger for the TSC. However, due to the broader peak in gclo(r) at SCW conditions, the total number of molecules in the first solvation shell is similar to that at ambient conditions. The coordination numbers will be discussed in detail below. In AW, the second peak for the TSC was noted4’ to be distinctly structured; this was attributed to the presence of the second chloride at the other extreme of the TSC. At SCW conditions, the distribution functions are less structured at distances farther than the first solvation shell. This is a
consequence of the large fluctuations in water density at near critical conditions. Near a critical point, density fluctuations have a range that far exceeds molecular dimensi0ns.4~ Correlation functions still exhibit a first maximum, although broader than at lower temperatures, provided that the potential energy well is comparable in size with RT,, with Tc the critical temperature. At distances further than the first solvation shell, the pair correlation function develops a relatively broad “tail” instead of the well-defined oscillations observed at higher density and lower temperatures. The excess number of water molecules at a distance R compared to a random distribution may be expressed as MX(R)= 4n@bhR(g(33(r)-
>r2dr
(4)
where @b is the bulk density. As shown in Figure 5 , the sharp increase in Nex (R) for the reactants and IDC resulting from the first solvation shell (between 3 and 4 A) contrasts the weak increase corresponding to the TSC. The increase is due to
J. Phys. Chem., Vol. 99, No. 5, 1995 1559
sN2 Reaction in Supercritical Water
’I
r
Figure 3. gso(r,z) corresponding to the transition-state complex at e, = 1.5 and T, = 1.0. The distribution of charges on the solute is symmetric (see Table 1).
electrostriction about the ionic C1-. However, after the first shell (r > 5 A) the slope for the TSC is more pronounced due to the larger tail for gao. This larger slope is due to the attraction of water to the second chloride at the other extreme of the TSC, which has a charge of -0.746 e. The radial distribution function g c ~ ~ (for r ) the leaving chloride is shown in Figure 6. For the case of the reactants and IDC in SCW, the C1’ is part of a neutral molecule, so that gc1to(r)has the features of a weakly attractive similar to that of methanol in SCW.29*50Here weakly attractive means that the first peak in the solute-water pair correlation function is lower than the one corresponding to the water-water interaction and that at distances farther than the first shell the solute-water g(r) approaches 1 from above.5o For the TSC, a significant first peak is present due to the development of more charge on Cl’. The g c 1 ~(r) of the TSC should be the same as the gc&) shown in Figure 4, due to the symmetry of the complex; differences are indicative of the quantitative accuracy of the results. The corresponding ambient water distribution
functions$2 also indicate a relative lack of structure for the water near the leaving chloride for the reactant geometry, although a first peak is present, in contrast to SCW. However, a greater change in gcyo than in gclo can be noted in going from AW at 25 “C to SCW. As expected, the increase in kT by a factor of 1.98 reduces the number of interactions between C1’ and H20 more than that of the stronger interactions between C1- and H20.
The distribution of water near the CH3 group is shown in Figures 7a,b. A large broad first peak from 4 to 6 8, is observed at supercritical conditions for the TSC. In ambient water,’” first and second peaks for the TSC are clearly defined and separated by a deep valley. The valley completely disappears in SCW. Also in contrast with AW, a peak is not present at 4.5 8, in gCH30 for the reactant state for the same reasons as discussed in Figure 6 for gC1’0. The valleys between each of the solvation shells in SCW are much shallower than in AW. gCH@(r) is compared for the reactants and the IDC in Figure 7b. There is a broad first peak for the IDC, centered about 6.5 A. This peak,
Balbuena et al.
1560 J. Phys. Chem., Vol. 99, No. 5, 1995
2.0 -Reactants
-Reactants
4
n
L
1.5 n
3 -
bl
W
0
v
6
pc)
1.0
b w
2-
0.5 1 -
0.0
1
3
5
7
9
r(A) Figure 4. Spherical 1-D pair distribution function reactants, IDC and TSC, at er = 1.5 and T,= 1.0.
I
’”[ n Lc W
l2
I
I
I
I
I
I
for the
gC&)
I
I
-Reactants
t
A
#
also similar. In SCW the f i s t peak is centered at about 2.5 A, with a slightly lower height than in AW. However, there is a higher local density at the minimum after the first peak, and the second peak is much broader, as noted already for gclo(r) (Figure 4). Coordination Numbers. Having discussed the distribution functions in detail, we now present a brief summary of the results in terms of the number of nearest neighbors in the first solvation shell, ns. This quantity is defined by the relationship
(5)
t
4t 1
Figure 6. Spherical pair distribution functiongC,,&) for the reactants, IDC and TSC,at er = 1.5 and TI = 1.0.
3
5
r
7
9
(A,
Figure 5. Excess number of water molecules (with respect to bulk) surrounding the chloride nucleophile at er= 1.5 and T, = 1.0.
along with that for the leaving chloride (Figure 6) and for the nucleophile (Figure 4), shows a very well developed structure surrounding the ion-dipole complex. We will refer to this fact when discussing the entropic contribution to the free energy of solvation. The solute-hydrogen pair correlation functions, gSH(Y), are shown in Figures 8 and 9 for the free nucleophile and the leaving chloride, respectively. For the nucleophile (Figure 8), two clear peaks are observed for the reactants, which are displaced to slightly larger radii with respect to AW.41 The relative local density at the second peak in SCW is higher than in AW, indicating a greater tendency of one of the H atoms in the water molecules in the first solvation shell to point outward than in AW. The corresponding curve for the TSC in SCW shows far less structure than in AW;41e.g., the first peak is much lower. It is evident that the number of hydrogen atoms associated with the chloride ion in both geometries decreases with respect to ambient conditions. For the reactants in SCW (Figure 9) g&r) resembles that in AW.41 For the TSC, the structures are
where Rm is the distance corresponding to the first minimum after the first peak in the solute atom(s)-water oxygen gso(r), i.e., enclosing the first solvation shell. For the case of the chloride nucleophile, ne1 is similar at each geometry and does not vary significantly from AW4* to SCW as shown in Table 2. For a single chloride ion, the same conclusion has been observed previously.28 The larger value for the TSC results from five water molecules being coordinated to each equivalent C1. The coordination to the CH3C1’ reactant is not important since these weak interactions contribute little to the differential solvation energy from the reactant to TSC, as seen earlier for AW.41 Because of the low density and high compressibility in SCW, attractive forces move molecules into energetically favorable locations more readily. To further understand solvation, it is thus useful to determine the extent of hydrogen bonding. One measure of the number of hydrogen bonds nm can be calculated by a geometric criterion using the expression analogous to Eq 5 but substituting the solute-hydrogen radial distribution functions and integrating to the minimum after the first peak. Unlike ncl, the value obtained from this integral, ~ H decreases B substantially from the reactants to the TSC (see Table 2 ) . A second method to determine nm is to examine the water-solute pair energy distribution. The number of water molecules interacting with the chloride nucleophile for a given energy is shown in Figure 10 for the reactant, IDC, and TSC geometries. The main peak centered near 0 kcallmol is characteristic of the large number of water molecules that are located far away from the solute. The peaks with a large negative energy (-- 12 kcal/ mol for the reactant and IDC) represent the H bond between
S N Reaction ~ in Supercritical Water 2.0
a
I
I
I
I
I
J
J. Phys. Chem., Vol. 99,No. 5, 1995 1561 I
I
I
I
I
4
I
I
l
I
l
I
l
-Reactants
1
-Reactants
I
I-TSC
l
1.5
3
n
h
bl
bl
W
a 1.0
3U
8 W
2
W
I
0.5
0.0
1
3
1
il
5
7
0
9 1 1 1 3
1 3 5 7 9 r(h
2.0 b
-Reactants 1-IDC
I
Figure 8. Spherical pair distribution function gc&) and TSC, at er = 1.5 and T, = 1.0.
2.0
I
I
I
for the reactants
/
-Reactants
1.5 1.5
n
bl
8 1.0
W
n
h
U
W
Ee 1.0
W
a
W
0.5 0.5
0.0
1
3
5
7
9
11
13
4)
Figure 7. Spherical pair distribution function gCHjO(r) at e, = 1.5 and T, = 1.0: (a) reactants and TSC; (b) reactants and IDC.
C1- and H20. In the TSC this peak is replaced by a shoulder with an average energy of only -7 kcaymol. These energy values are in agreement with those for AW4' although the distributions of these energies are different. The number of H bonds calculated by integrating these energy distribution peaks from valley to valley is shown in Table 2. The values from the energetic criterion are in good agreement with those determined geometrically. The new results for SCW along with earlier results of Jorgensen et al4I for AW are reported in Table 2. The ratio of the number of H bonds (nm) to the coordination number (ncl) for the nucleophile decreases from 0.9 in AW to 0.68 in SCW for the case of the reactants. Even though the coordination number remains high in SCW, the number of hydrogen bonds decreases. This decrease becomes more pronounced as the strength of the hydrogen bond decreases from the reactant state to TSC. The larger disruption for weaker hydrogen bonds by any perturbation is to be expected. Because the coordination number is essentially constant, hydrogen bond interactions must be converted to less specific polar interactions by thermal energy and the reduction
0.0
3
5
7
9
d) Figure 9. Spherical pair distribution function gclq(r) for the reactants and TSC, at e, = 1.5 and T, = 1.0. TABLE 2: Coordination Numbers and Number of H Bonds for Several Geometries in SCW and AW nm
nmm
geometric energetic energetic Rm- nC1 criterion criteriono criterion geometry reactants, SCW 1.4000~ 7.7 5.0 5.3 0.68 IDC, SCW 4.1 4.4 1.5~00 7.5 0.58 TSC, SCW 1.5000 10 3.5 3.8 0.38 reactants, AW4' 1.4000 7 6.3 0.90 TSC, AW4' 1.4000 8 6.6 0.83 For the reactants and IDC we use the upper energetic limit of -7.17 and -4.6 kcaUmol for the TSC. 000 is the Lennard-Jones diameter for the oxygen-oxygen interaction [3.166 8, for the SPC model (this work) and 3.153 8, for l"4P (ref 41)]. in density. Thus, the differences in solvation in AW and SCW are somewhat greater than would be expected based on Its alone. Discussion of AA, AS,and AE in Terms of the Microscopic Solvation. Previously, we have reported39 the energetic and
Balbuena et al.
1562 J. Phys. Chem., Vol. 99, No. 5, 1995
TABLE 3: Solvent Reorganization Energies in SCW and AW
2.0
~~~
geometry reactant, SCW IDC, SCW TSC, SCW reactant, AW41 TSC, AW4I
1.5
( E,, )IN,"
( E,, )*IN, (bulk)
-4.95 -4.92 -4.98 -9.73 -9.79
-5.05
A&,
48.4 58.1 29.0 72.5 57.5
-5.05
-5.05 -10.02 -10.02
The energies are in kcaVmol of solute, and N, is the number of water molecules in the simulation sample.
1.o
MSCW= 19.3kcal/mol AESCW = 22.8 kcal/mol
0.5
TASSCW = 3.5 kcallmol
[ Reactants ] SCW
0.0 -15
"SCW / AW
-10
0
-5
10
5
Energy, kcal/mol
CI
CI
T
f
1.1 "2
7'8
Figure 10. Distribution of pair interaction energies for water molecules with the nucleophile for the reactants, IDC and TSC at el = 1.5 and T,
0.58
= 1.0.
'
[ Reactants ] AW
25
[TSClsCW
T
I1
ITSCIAW
MAW= 22.7 kcallmol
20
AEAW
m
15
TASAW = -0.7 kcalhol
v
y"
10
h
5
= 22.0 kcal/mol
Figure 12. Thermodynamic cycle for the energetic changes. vertical paths, AW to SCW; horizontal paths, reactants to TSC.
9 E o W -5 0
2
4
r
A
6
8
(A>
Figure 11. Calculated changes in free energy in solution and its energetic and entropic contributions as a function of the reaction coordinate, rA,
entropic contributions to the free energy at the three geometries (Figure 11). Our present objective is to rationalize these results in terms of the structural features of solvation. The properties AA(IA),~ ( I A )and , -TAS(~A)refer to the differences between the value of each property at rA and the corresponding value at lrAl = 8 8, (taken as the separated reactants). The energy curve AE(rA) increases monotonically from the reactants to the TSC ( r A = 0), reflecting the weakening of water-solute interactions concomitant with charge dispersal. The entropic curve is also informative. As Ir~ldecreases from 8 to 1.441 8, at the IDC, AS decreases. This ordering of solvent does not occur at AW conditions, as indicated by integral equation calculation^.^^ As is decreased from the IDC to the TSC, AS increases and becomes positive due to loss of electrostriction from solute charge delocalization. However, since the net value of this entropic contribution to the free energy at the TSC is small compared to the energetic term, the free energy remains large and positive and the activation barrier remains high. The decrease in AS for the IDC in SCW may be understood with the new structural information. The analysis has to consider not only the structure surrounding the C1-, given by g c d r ) , but also the other two rdf's, gcl,o(r) and gCH,O(r). A comparison of the IDC and reactants shows that the water structure near the nucleophile is mostly preserved (see Figure
4), whereas water near the leaving chloride (Figure 6) and the methyl group (Figure 7b) becomes more structured due to the small increase in charge for C1' in the model. Here the shift of a little charge from C1- to C1' increases the total solvation. The gain in structure from the reactants to the IDC can also be seen in the 2-D pictures in Figures 1 and 2. To further support this point, we have calculated the solvent reorganization energy, defined as M,, = ( E,, ) - ( E,, )*, where ( E,, ) is the average water-water potential energy with solute present and ( E, )* is the value for pure water under the same conditions of density and temperature. As AE,, is a partial molar quantity, formal derivations5' show that there is a nonzero contribution from the change in solvent structure, relative to pure solvent, even at infinite dilution. The results for the solvent reorganization energy are shown in Table 3. The last column clearly shows that the solvent disruption energy in SCW increases from the reactant state to the IDC. This greater energetic disruption of solvent is consistent with the greater ordering of solvent about the solute. Notice solvent disruption energies are smaller in SCW than for more ordered AW; however, they are still significant in SCW. Although the decrease in entropy and increase in energy (Figure 11) oppose the formation of the IDC, the intrinsic (gas phase) ion-dipole attraction leads to a minimum in the PMF. In AW, AE is sufficiently unfavorable that the IDC is not present.39 A summary of the key properties in this study is presented in Figure 12. The values of AA, AE and TAS are given for AW41 (bottom) and SCW (top) in the horizontal direction. To explain these changes along the reaction coordinate, we examine the solvation of the reactant and TSC states separately from AW to SCW conditions. This solvation is indicated by the vertical lines. Specifically, the figure lists changes in or ratios of the coordination number about the nucleophile, the number and strength of hydrogen bonds, and the solvation energy for
J. Phys. Chem., Vol. 99,No. 5, 1995 1563
sN2 Reaction in Supercritical Water
1.5
12 10
-
-pr=l.5 p,=2.0
1 .o
8
s li
n h
n
h
W
w
8 om
P t
w
6 4
0.5
2
0
0.0 1
3
7
5
9
r(A) Figure 13. Effect of temperature on the spherical pair distribution function gcro(r) at er = 1.5.
2
I
I
I
I
I
I
Ib
I
--t
pr =0.5
pr= 0.9
SCW (top) compared to AW (bottom). This approach can provide a general understanding about how solvation changes from AW to SCW as a function of charge on the reacting species. As described earlier, the number of water molecules about the C1- nucleophile changes little from AW to SCW in the reactant state, and likewise for the TSC (see Figure 4 and Table 2). Furthermore, the loss of solvation energy is similar for the reactant and TSC states. These two factors cause AE at the TSC to be about the same in AW and SCW. In contrast, the number of hydrogen bonds decreases significantly from AW to SCW, particularly for the TSC, since it is the weaker H bond acceptor. While this does not appear to affect AE, it does influence AS. Since the loss of hydrogen bonding or electrostriction is larger for the TSC than the reactant state, TAS increases from -0.7 in AW42to +3.5 in SCW. The change in free energy for the four steps in the cycle in Figure 12 must sum to zero, or equivalently
+ rASCW
- AAVf]TSC = hASCW
+ LASCW
- AAW]reactants ( 6 )
The values of AA are known in AW and SCW. Thus, the positive change in free energy for the reactant state from AW to SCW is larger than that for the TSC by 3.4 kcdmol. The corresponding changes in E are similar for the reactant and TSC. Thus, the changes in free energy are entropy driven. The change in the function TS from AW to SCW for the TSC is positive and larger than for the reactant state by 4.2 kcal/mol. Physically, this indicates that the loss of electrostriction is more prevalent for the TSC, since it forms weaker H bonds with water. Effect of Temperature. The effect of temperature at supercritical conditions was investigated by comparison of two states: TI = 1 and T, = 1 . 1 at the same reduced density PI = 1.5 for the IDC geometry. The pair correlation function gclo(r) decreases slightly at the higher temperature, reflecting the presence of a smaller number of water molecules surrounding each of the sites of the complex (not shown). As remarked with reference to Table 1, the desolvation caused by increasing temperature is much more pronounced in the case of the less charged atom (Cl', Figure 13) than for the nucleophile, since C1- interacts more strongly with the water. An example of this behavior results from the comparison of the behavior in AW
-pr -pr
1
3
5
7
9
= 1.5
= 2.0
11
13
Figure 14. Effect of solvent density at T, = 1 on (a) the pair distribution function gc&) and (b) the local density eclo(r). vs SCW of gcyo for the reactants (where C1' is in the less polar CH3Cl) and the C1' in the TSC (polar complex), as discussed with reference to Figure 6. This will have important consequences on reactions involving polar and nonpolar species. Other small changes are observed with temperature. For example, the coordination number in the first solvation shell ncl is 7.3 for R" = 1.5aoo at TI = 1.1, compared to 7.5 at TI=l. Also, the numbers of H bonds calculated from the pair energy distribution for the water molecules interacting with the nucleophile are 4.10 and 3.97 at TI = 1.0 and TI = 1.1, respectively, based on the solute-hydrogen rdf's. Thus, the ratio nmh" is 0.56, slightly lower than the value 0.58 at TI = 1.0, at the same density. Effect of Water Density. In order to investigate the effect of water density on the structural properties, we present results from simulations of the IDC at er= 2.0, 1.5, 0.9, and 0.5, for the same reduced temperature T, =l. As shown for gclo(r) in Figure 14, the electrostriction of water molecules in the first and second solvation shells increases significantly as density decreases. Also, as e decreases, the number of water molecules near C1- is nearly constant, as shown in Figure 14b, where the local density eclo(r) is calculated as the product of g&r) and e b , in units of g/cm3. similar features can be seen in the
1564 J. Phys. Chem., Vol. 99, No. 5, 1995
Balbuena et al.
TABLE 4: Coordination Numbers and Number of H Bonds for the IDC at Various Densities nm
geometric
er
ncla
criterion
0.5 0.9 1.5
7.2 7.0 7.5 7.8
4.22 4.11 4.13
2.0
4.21
energeticb criterion 4.40 4.27 4.37 4.42
nHBlnC1
energetic criterion 0.61 0.59 0.58
0.57
The coordination numbers ncl were calculated with a value of R" = 1.5 a00 in eq 5. An upper energetic limit of -7.17 kcaVmol was used in all cases. case of the leaving chloride and also for the methyl group (not shown). The number of nearest neighbors for the chloride nucleophile in the IDC as a function of water density is shown in Table 4. This number decreases only slightly when the density decreases. Clearly, the local density of water about each solute remains constant as shown by Figure 14b, indicating that substantial clustering persists even at the lower densities. Furthermore, the persistence in clustering with decreasing density, is greater for the more ionic species C1 than the less charged Cl'. On the basis of these results, we further elaborate on the three regions of solute-solvent clustering defined by Sun et al., which have been observed in several experimental s t u d i e ~ . ' ~The -~~ description begins with ambient water in the liquidlike region where spectroscopic measures of solvation exhibit linear changes with density.I4 In this liquidlike region, the outer waters are removed from the solute as density decreases, and the spectral shifts may be described with continuum reaction field theories which are relatively linear in density. In the near-supercritical region, from a er of 2 to 0.5, clustering prevents disruption of the inner solvation water molecules, and coordination numbers change little. Here a plateau is observed in solvatochromic shifts for a variety of probes. These clusters present under warm, high-pressure conditions are somewhat more disordered than the usual cold, low-pressure clusters formed in adiabatic expansions. As density is lowered to gaslike values, experimental data14,16-18,20 indicate that water molecules evaporate from the gas phase clusters. The numbers of H bonds ( n m )between the nucleophile and water were calculated from the pair distribution energies and from the solute-hydrogen rdf's for the IDC at each density and are presented also in Table 4, along with the ratio nm/ns for each case. It is seen that this ratio is slightly lower at the higher densities, For the nonpolar solute acetone, it was demonstrated spectroscopically that the degree of acetone-water hydrogen-bonding changes little in this range, in agreement with the present simulations. Conclusions
Simulations in supercritical fluids have evolved rapidly from simple Lennard-Jones mixtures, to polar mixtures, to pure water, to ionic aqueous solutions, and now to chemical reactions in water including ions and hydrogen bonding. The present simulations provide detailed structural information about how the coordination number and hydrogen bonding change with density and temperature for solutes of varying charge. These are key elements for understanding the effect of the solvent on the activation barrier. Under these supercritical conditions the barrier is found to be intermediate between ambient water and gas phase, though much closer to AW. The solute-water correlation functions in SCW show broader peaks and higher valleys than in AW. This higher delocalization in the molecules
is due to the decrease in density and also to the higher temperatures and larger density fluctuations. The coordination number of water about highly polar and ionic solutes in SCW remains almost as high as in ambient water, even for reduced densities of 0.5. The solute-solvent interactions are sufficiently strong such that clusters are prevalent at these low supercritical densities. The persistence in clustering with decreasing density is greater for the more ionic species C1- than the less charged Cl'. Thus, clustering increases as the solute-solvent interaction energy increases relative to the pure solvent energy, as shown clearly in the 2-D pictures. Even though the coordination number remains high in SCW, the number of hydrogen bonds decreases, due to the increase in kT and reduction in density. Thus, hydrogen bond interactions must be converted to nonspecific polar interactions. This loss of hydrogen bonds from AW to SCW decreases as the hydrogen bond strength increases for the series of solutes from Cl' and C1 in the TSC, to C1- in the IDC, to free C1-. The 2-D and spherical radial distribution functions, along with the hydrogen bond energy distribution functions, provide a clear explanation of the origins of AA, AE, and TAS along the reaction coordinate. Because the coordination number remains high in SCW about the C1 nucleophile in the reactant and TSC, AE is similar in AW and SCW. It increases monotonically from the reactants to the TSC, reflecting the weakening of water-solute hydrogen bonds and polar interactions with charge dispersal. The modest decrease in the free energy barrier AA in SCW relative to AW is entropically driven. The hydrogen bonding of water about the reactant state decreases less than that about the TSC from AW to SCW; consequently TAS becomes positive. Another factor that causes the increase in TAS is the decrease in the number and strength of the Cl--water hydrogen bonds from the reactant state to TSC in SCW. In AW there is a decrease in the strength and not in the number. The decrease in entropy from the reactant state to the IDC in SCW is the opposite of the behavior in AW. Upon transfering a little charge from C1- to the C H F 1 dipole, the increase in solvation of the dipole exceeds the loss of solvation about the C1- ion. This trend undergoes a reversal as the reaction proceeds forward, and the entropy changes sign. In future work we will address this reaction at a lower density, er= 0.5, and at higher temperature^.^^ We anticipate that future studies of clusters at gas densities will complement this study of solvation in the other two regions characterized by liquidlike and supercritical densities. Acknowledgement is made to the U.S.Army Research Office for University Research Initiative Grant DAAL 03-92-6-0174 and DAAH 04-93-6-0363 and to the Separations Research Program at The University of Texas. Partial support of this work by a grant from the R. A. Welch Foundation (F-0761) to P.J.R. is also acknowledged. We thank The University of Texas System Center for High Performance Computing and Cray Research Inc. for computational support and Mr. Randy Franklin for help with the graphics. Partial support of this work was also provided by the Army Research Office Contract DAAL03-89-C-0038 with the University of Minnesota Army High Performance Computing Research Center (AHPCRC) and the DoD Shared Resource Center at the AHPCRC. References and Notes (1) Shaw, R. W.; Brill, T.B.; Clifford, A. A.; Eckert, C. A,; Franck, E. U. Chem. Eng. News 1991, 69, 26. (2) Tester, J. W.; Holgate, H. R.; Armellini, F. J.; Webley, P. A.; Killilea, W. R.; Hong, G. T.;Barner, H. E. In ACS Symposium Series;
S N Reaction ~ in Supercritical Water Tedder, D. W., Pohland, F. G., Eds.; American Chemical Society: Washington DC, 1993; Vol. 518, p 35. (3) Gloyna, E. F.; Lee, L. Waste Manage. 1993, 13, 379. (4) Antal, M. J.; Brittain, A.; DeAlmeida, C.; Ramayya, S.; Roy, J. C. In ACS Symposium Series; Squires, T. G., Paulaitis, M. E., Eds.; American Chemical Society: Washington, DC, 1986; Vol. 329, p 77. ( 5 ) Townsend, S. H.; Abraham, M. A.; Huppert, G. L.; Klein, M. T.; Paspek, S. C. Ind. Eng. Chem. Res. 1988, 27, 143. (6) Wood, R. H.; Quint, J. R.; Grolier, J. P.E. J. Phys. Chem. 1981, 85, 3944. (7) Levelt Sengers, J. M. H. In Supercritical Fluid Technology; Ely, J., Bruno, T. J., Eds.; CRC Press: Boca Raton, FL, 1991; p 1. (8) Mesmer, R. E.; Marshall, W. L.; Palmer, D. A.; Simonson, J. M.; Holmes, H. F. J. Solution Chem. 1988, 17, 699. (9) (a) Flarsheim, W. M.; Tsou, Y. M.; Trachtenberg, I.; Johnston, K. P.; Bard, A. J. J. Phys. Chem. 1986, 90, 3857. (b) Flarsheim,W. M.; Bard, A. J.; Johnston, K. P. J. Phys. Chem. 1989, 93, 4234. (10) Xiang, T.; Johnston, K. P. J. Phys. Chem. 1994, 98, 7915. (11) Buback, M.; Crerar, D.; Koplitz, L. M. In Hydrothermal Experimental Techniques; Ulmer, G., Barnes, H. L., Eds.; Wiley: New York, 1987; p 333. (12) Postorino, P.; Tromp, R. H.; Ricci, M.-A.; Soper, A. K.; Neilson, G. W. Nature 1993, 366, 668. (13) Spohn, P. D.; Brill, T. B. J. Phys. Chem. 1989, 93, 6224. (14) Bennett, G.; Johnston, K. P. J. Phys. Chem. 1994, 98, 441. (15) Kim, S.; Johnston, K. P. Ind. Eng. Chem. Res. 1987, 26, 1206. (16) Knutson, B. L.; Tomasko, D. L.; Eckert, C. A.; Debenedetti, P.; Chialvo, A. A. In ACS Symposium Series; Bright, F. V., McNally, M. E. P., Eds.; American Chemical Society: Washington, DC, 1992; Vol. 488; Chapter 5 . (17) Sun, Y.-P.; Fox, M. A.; Johnston, K. P. J. Phys. Chem. 1992,114, 1187. (18) Kajimoto, 0.; Futakami, M.; Kobayashi,T.; Yamasaki, K. J. Am. Chem. SOC.1988, 92, 1347. (19) Johnston, K. P.; Haynes, C. AIChE J. 1987, 33, 2017. (20) Carlier, C.; Randolph, T. W. AIChE J. 1993, 39, 876. (21) Tanger IV,J. C.; Pitzer, K. S. J. Phys. Chem. 1989, 93, 4941. (22) Gupta, R. B.; Johnston, K. P. Ind. Eng. Chem. Res., in press. (23) dePablo, J. J.; Prausnitz, J. M.; Strauch, H. J.; Cummings, P. T. J. Chem. Phys. 1990, 93, 7355. (24) Berendsen, H. J.C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. In Intemleculur Forces; Pullman, B., Ed.; Reidel: Dordrecht, 1981. (25) Walrafen, G. E.; Fisher, M. R.; Hokmabadi, M. S.; Yang, W.-H. J. Chem. Phys. 1986, 85, 6970. (26) Mountain, R. D. J. Chem. Phys. 1989, 90, 1866. (27) Gupta, R. B.; Johnston, K. P. Fluid Phase Equilib. 1994, 99, 135. (28) Cummings, P. T.; Cochran, H. D.; Simonson, J. M.; Mesmer, R. E.; Karaborni, S. J. Chem. Phys. 1991, 94, 5606.
J. Phys. Chem., Vol. 99,No. 5, 1995 1565 (29) Cochran, H. D.; Cummings, P. T.; Karabomi, S. Fluid Phase Equilib. 1992, 71, 1. (30) Cui, S. T.; Hanis, J. G. Chem. Eng. Sci. 1994, 49, 2749. (31) Gao, J. J . Phys. Chem. 1994, 98, 6049. (32) Gao, J. J. Am. Chem. SOC.1993, 115, 6893. (33) Guillot, B.; Guisani, Y. J. Chem. Phys. 1993, 99, 8075. (34) Kalinichev, A. G.; Heinzinger, K. In Advances in Physical Geochemistry; Saxena, S., Ed.; Springer-Verlag: New York, 1992; Chapter 1. (35) Combes, J. R.; Johnston, K. P.; O’Shea, K. E.; Fox, M. A. InACS Symposium Series; Bright, F. V., McNally, M. E. P., Eds.; American Chemical Society: Washington, DC, 1992; Vol. 488; Chapter 3. (36) Kazarian, S. G.; Gupta, R. B.; Clarke, M. J.; Johnston, K. P.; Poliakoff, M. J. Am. Chem. SOC.1993, 115, 11099. (37) Chateauneuf, J. E.; Roberts, C. B.; Brennecke, J. F. In ACS Symposium Series; Bright, F. V., McNally, M. E. P., Eds.; American Chemical Society: Washington, DC, 1992; Vol. 488; Chapter 9. (38) Randolph, T. W.; Clark, D. S.; Blanch, H. W.; Prausnitz, J. M. Proc. Natl. Acad. Sci. U.S.A. 1988, 85, 2979. (39) Balbuena, P. B.; Johnston, K. P.; Rossky, P. J. J. Am. Chem. SOC. 1994, 116, 2689. (40) Archer, D. G.; Wang, P. J. Phys. Chem. Re$ Data 1990, 19, 371. (41) Chandrasekhar, J.; Smith, S. F.; Jorgensen, W. L. J. Am. Chem. SOC.1985, 107, 154. (42) Huston, S. E.; Rossky, P. J.; Zichi, D. A. J. Am. Chem. SOC.1989, 111, 5680. (43) Zwanzig, R. W. J. Chem. Phys. 1954, 22, 1420. (44) Fleischman, S. H., Brooks III, C. L. J. Chem. Phys. 1987,87,3029. (45) Fois, E. S.; Sprik, M.; Paninello, M. Chem. Phys. Lett. 1994, 223, 411. (46) (a) Kollman, P. Chem. Rev. 1993, 93,2395. (b) Pearlman, D. A. J. Phys. Chem. 1994, 98, 1487. (47) Zhu, S.-B.; Wong, C. F.J. Phys. Chem. 1994, 98, 4695. (48) Jorgensen, W. L.; Chandrasekhar, J.; Buckner, J. K.; Madura, J. D. In Annals of the New York Academy of Sciences; Beveridge, D. L., Jorgensen, W. L., Eds.; New York Academy of Sciences: New York, 1986; Vol. 482, p 198. (49) Levelt Sengers, J. M. H. In NATO AS1 Series; Kiran, E., Levelt Sengers, J. M. H., Eds.; Kluwer Academic Publishers: Dordrecht, 1994; Vol. E 273, p 3. (50) Cummings, P. T.; Chialvo, A. A.; Cochran, H. D. Chem. Eng. Sci. 1994,49, 2735. (51) (a) Garisto, F.; Kusalik, P. G.; Patey, G. N. J. Chem. Phys. 1983, 79, 6294. (b) Yu, H.-A.; Karplus, M. J. Chem. Phys. 1988, 89, 2366. (52) Flanagin, L. W.; Balbuena, P. B.; Johnston, K. P.; Rossky, P. J. J. Phys. Chem., submitted.
JP941406Z
