J . Phys. Chem. 1988, 92. 5257-5261 systematically vary with the group in the periodic table.) Any proportional correction would maintain the linear relation between AAfHo and V, found here. Pauling's electronegativity scale7 was derived from the heats of formation of diatomic molecules and was based on a qualitative argument concerning ionic character in the bonds. While it was of great utility in organizing the data on heats of formation on diatomic molecules, it could never be extended quantitatively to polyatomic molecules. Yuan did not attempt to use V, for estimation. He simply observed that there was a nearly linear relation between V, and Pauling's values of electronegativity X,.This was not a proper comparison since energies of bonds are proportional to the square of differences in X , for the two elements forming a bond. In fact a comparison of V, with XX2 gives a better linear relation, but it is far from quantitative. However, the biggest discrepancy between X, and V, is in the value assigned to hydrogen. It is very close to the value of carbon in the Pauling scale and almost all other scales but very different from the V, for carbon. Intuitively V, or something proportional to V, seems a quite reasonable measure of the strength of the covalent bond, and so it is perhaps not too surprising that it affords such an excellent correlation. The relation that we have found for AfHo(RX/MeX) permits us to estimate AfHo(RX) for alkyl derivatives when AfHo(MeX) is known. It also allows us to estimate in turn group contributions
5257
for C-X(H)m-3(C)m. This in turn with other groups that are known18 permit us to estimate AfHo(RX) for any single R of almost any complexity. Conversely, if AfHo(RX) are available for any one more complex molecule, eq 5 together with other group values permits us to estimate AfHo(RX) for all other molecules, including CH3X. In subsequent papers we will discuss the extension of this approach to molecules in which X is a group such as Et, i-Pr, t-Bu, C2H3, C3H5, C6H5, C6HsCH2, RCO, HOCO, NC, RO, O2N0, ONO, RCOO, ON, OzN, CN, RS, R S 0 2 , etc. Acknowledgment. This work has been supported by grants from the National Science Foundation (CHE-84-0376 1) and the US. Army Research Office (DAAG29-85-K-0019). Registry No. CH3F, 593-53-3; CH,OH, 67-56-1; C H Q , 74-87-3; CH3NH2,74-89-5; CH3Br, 74-83-9; CH3SH, 74-93-1; CHJ, 74-88-4; CH3CH3, 74-84-0; CH4,74-82-8; CZHSOH, 64-1 7-5; CZHSCI, 75-00-3; C2HSNH2,75-04-7; C2H5Br,74-96-4; CzHSSH,75-08-1; C2H51,75-03-6; CH3CHzCH3, 74-98-6; i-C3H7F, 420-26-8; i-C3H70H, 67-63-0; iC3H7C1,75-29-6; i-C3H7NHz,75-31-0; i-C3H7Br,75-26-3; i-C3H7SH, 75-33-2; i-C3H71, 75-30-9; CH(CH&, 75-28-5; t-C4HgOH, 75-65-0; t-C4H,C1, 507-20-0; t-C4H9NH2, 75-64-9; t-C4HgBr, 507-19-7; tC4HgSH, 75-66-1; t-C.+HJ, 558-17-8; CH,C(CH,)j, 463-82-1. (18) Benson, S. W. Thermochemical Kinetics, 2nd ed.;Wiley: New York, 1976.
Hypernetted Chain Closure Reference Interaction Site Method Theory of Structure and Thermodynamics for Alkanes in Water Toshiko Ichiye and David Chandler* Department of Chemistry, University of California, Berkeley, California 94720 (Received: March 3, 1988)
The solvation structure and thermodynamics of hydrophobic molecules in aqueous solution are studied by the extended reference interaction site method using the hypernetted chain closure (HNC-RISM) and by molecular dynamics simulations. The solvation structure predicted by HNC-RISM agrees only moderately with simulation results. In addition, thermodynamic data from the theory are in poor agreement with experiment, although this may be due at least partially to the potential functions used. These results indicate that HNC-RISM is not a quantitatively accurate theory for nonpolar solutes in aqueous solution. Since problems have already been noted for ionic solutions, this reinforces the need for an improved closure for realistic continuous potential models. A method for obtaining potentials for reduced models is also presented, and the use of HNC-RISM to implement this method is discussed.
I. Introduction The solvation of hydrophobic molecules in water is of great importance in many problems of biological and chemical interest. Consequently, a large number of computer simulations' and more analytic theoretical studies2 have been made on this problem. In this article we study integral equation predictions and molecular dynamics results for the structure and thermodynamics of molecules in aqueous solution. Although a study of hydrophobic solvation using the reference interaction site method has been made by Pratt and Chandler," many new developments since then make a new study of interest. First, an extension of the RISM integral equation using the H N C closure by Rossky and co(1) See, for example: (a) Pangali, C.; Rao, M.; Berne, B. J. J. Chem. Phys. 1979, 71,2975; 1979, 71, 2982. (b) Swope, W. C.; Andersen, H. C. J. Phys. Chem. 1984.88.6548. (c) Joreensen. W. L.: Gao. J.: Ravimohan. C. J. Chem. Phys. 1985,'89,'3470. id) StLatsma, T. P.'; Berendsen, H. J. C.; Postma, J. P. M. J. Chem. Phys. 1986,85, 6720. (2) See, for example: (a) Pratt, L. R.; Chandler, D. J . Chem. Phys. 1977, 67, 3683. (b) Pratt, L. R.; Chandler, D. J . Solution Chem. 1980, 9, 1. (c) Tanaka, H. J . Chem. Phys. 1987,86, 1512. (3) See, for example: (a) Chandler, D.; Andersen, H. C. J . Chem. Phys. 1972, 57, 1930. (b) Lowden, L. J.; Chandler, D. J . Chem. Phys. 1974, 61, c**o
J'LLO.
(4) Pratt, L. R.; Chandler, D. J . Chem. Phys. 1980, 73, 3430, 3434. (5) Hirata, F.; Rossky, P. J. Chem. Phys. Lett. 1981, 83, 329.
0022-3654/88/2092-5257$01.50/0
workers, referred to here as HNC-RISM, appears applicable to physically realistic models of molecules having strong polar interactions. Indeed, it has been used with varied success in the study of liquid water structure6 and ionic solutions.' Second, energy parameters for hydrocarbonss and waterg have been developed that give good structural and thermodynamic properties in computer simulations. These developments have motivated us to make several new applications of the HNC-RISM theory. We have three main objectives: First, we compare experimental thermodynamic data for hydrophobic solvation with our HNCRISM results. This is important in light of a study that has been made attempting to predict the conformational equilibria of butane in aqueous solution by using HNC-RISM.'O Also a limited (6) Pettitt, B. M.; Rossky, P. J. J . Chem. Phys. 1982, 77, 1451.
(7) (a) Pettitt, B. M.; Rossky, P. J. J . Chem. Phys. 1986, 84, 5836. (b) Hirata, F.; Rossky, P. J.; Pettitt, B. M. J. Chem. Phys. 1983, 78, 4133. (c) Chiles, R. A,; Rossky, P. J. J . Am. Chem. SOC.1984, 106,6867. (d) Kuharski, R. A,; Chandler, D. J. Am. Chem. SOC.1987, 91, 2978. (8) Jorgensen, W. L.; Madura, J. D.; Swensen, C. J. J . Am. Chem. SOC. 1984, 106, 6638. (9) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J . Chem. Phys. 1983, 79, 926. (10) Zichi, D. A.; Rossky, P. J. J . Chem. Phys. 1986, 84, 1712. (11) The calculations were performed on the Berkeley CRAY XMP by adapting a code written by R. W. Impey.
0 1988 American Chemical Society
5258
The Journal of Physical Chemistry, Vol. 92, No. 18, 1988
comparison of simulation and experiment for a more developed parameter set than used in the butane calculation gives rather poor agreement for hydrocarbons in aqueous solution.7c This poor agreement can be juxtaposed with the rather good results found from simulations for each of the pure liquid phase^.^,^ It seems, therefore, that these systems are still not well understood, and further study seems warranted. Second, we examine how well the HNC-RISM can predict structural features of this system, independent of the parametrization. To accomplish this, we have also carried out molecular dynamics simulations with the same intermolecular potentials. These calculations, reported herein, illustrate problems with HNC-RISM not entirely different from those already encountered in studies of ionic s o l ~ a t i o n . ~ ~ ~ ~ - ~ As a third objective, we suggest a way of going from a representation in which sites are located at every atom, for which there are known interatomic potentials, to one in which many atoms are represented by a single site. The idea here is to develop simplified reduced models that retain the correct structural features of the solvation of the more complex full-site representation. We have in mind applications for long hydrocarbon chains, such as found in amphiphilic molecules. A brief discussion of the methods used is given in section 11. In section 111, we present structural and thermodynamical results for propane in water. The focus on propane is motivated by our interest in reduced models for hydrocarbon chains in which three successive carbons are replaced by a single site. In section IV, we describe the reduced model and present results for it. Concluding remarks are given in section V. 11. Methods
The HNC-RISM formulation has been described el~ewhere,~ so only the resulting equations are given here. For a solute at infinite dilution, the Ornstein-Zernike-like equation of Chandler and Andersen3 for the radial distribution is
Ichiye and Chandler TABLE I: Model Parameters 104~2,
methane CHa propane CH3 propane CH, water 0 water H
A6 kcal/mol
885.9 880.1 593.5 58.0 8.7 X lod0
-3 167.1 -2482.1 -1673.7 -525.0
A
bond length, propane . . water a
B2,
AI2 kcal/mol
CH,-CH1, 1.53 0-H, 0.9372
0 0 -0.80 +0.40
bond angle, deg CH1-CH,-CH,. 112 H-0-H,i04.5-
another derivation based on the assumption of Gaussian fluctuations for the solvent gives13 a similar result: -@fisolvG=
p E j d r c,,(r) ass
+
S,S’
The superscripts H N C and G indicate the hypernetted chain approximation and Gaussian model, respectively. We will use both equations for comparison with experiment. The molecular dynamics results are for a system of 249 water molecules plus one propane with periodic boundary conditions.’ The simulation was run for 100 ps with a time step of 0.001 ps after an equilibration period of 10 ps. The average temperature was 25 OC. The molecular geometry and integration potentials used are Jorgensen’s TIPS modeli4 for liquid water as modified by Pettitt and Rossky6 and Jorgensen’s OPLS model for propane.8 Zichi and Rossky’O used similar models for their study of hydrocarbons in water. The interaction potential is of the TIPS format:
+ -BaBs +r6-
where u,(r) is the pair potential between a solute and solvent sites a and s. The procedure used for solving these equations is the same as that employed by Kuharski and Chandler.7d The pure solvent correlation functions xss,(r)(see eq 1) can be calculated from simulation, experimental data, or HNC-RISM for the pure liquid. Although the HNC-RISM and simulation gsst(r)’sare qualitatively like experimental results, there are quantitative differences between all of the results. However, we have found that the effects of these differences are often relatively small in calculations of thermodynamic properties of the solute in solvent. Results given here are for Xsst(r)from HNC-RISM, unless otherwise noted. The excess chemical potential of solvation, hLs, can be computed from a closed form expression in which standard coupling parameter integration is applied to HNC-RISM theory.I2 This expression is
0
Used in water-water interactions.
AaAs ua, = rI2 where the Greek and Roman subscripts refer to interaction sites on the solute and solvent, respectively, p is the bulk number density , is the matrix of the solvent at temperature T = ( @ / k B ) - ’w(r) of intramolecular site-site correlation functions, and X(r) is the matrix of site-site density pair correlation function of the pure solvent. The HNC-like closure is
q, electrons
zmzs
r
(4)
The parameters are summarized in Table I. 111. Results for Solvation of Propane in Water
Thermodynamic of Solvation. The thermodynamic properties of propane in water were calculated from the HNC-RISM results. The excess chemical potential of solvation was calculated via eq 3a and 3b at several temperatures (Figure la). The entropy (Figure lb) and enthalpy (Figure IC) of solvation were calculated via finite differences. Finally, a measure of the hydrophobic effect,15dAHI= AbLpr- 3ApMo where Pr and Me refer to propane and methane, respectively, is shown in Figure Id. In each case, the expression based on a Gaussian bath (eq 3b) gives results that are closer to experiment than the expression that is consistent with the H N C closure (eq 3a). The poor agreement of the calculations with experiment is perhaps surprising, especially considering the success of the earlier work by Pratt and Chandler.4 However, there are several important differences in the two calculations. First, the PrattChandler values of u and t for the solute-water interactions are significantly different from those used here. Second, in the Pratt-Chandler theory, the xoo(r) was obtained from experiment and not theory. Third, the closure adopted by Pratt and Chandler is of the Percus-Yevick (PY) from rather than HNC. In effect, the difference in closure means that the interactions between the solute and the proton sites were ignored on the Pratt-Chandler calculations. This is true because in the PY closure, cOH(r)is zero when umH(r) = 0. If we use these same conditions in our calcu(13) Chandler, D.; Singh, Y.; Richardson, D. M. J . Chem. Phys. 1984,81.
s,s‘
where the circumflexes denote Fourier transforms. Alternatively, (12) Singer, S.; Chandler, D. Mol. Phys. 1985, 55, 621
1975.
(14) Jorgensen, W. L. J. Am. Chem. SOC.1981, 103, 335, 341, 345. (15) Ben-Naim, A. J . Chem. Phys. 1971, 54, 1387; 54 3696; 1972, 57, 5257; 1972, 57, 5266. Ben-Naim, A,; Wilf, J.; Yaacobi, M . J . Phys. Chem. 1973, 77, 95. Ben-Naim, A,; Yaacobi, M. J . Phys. Chem. 1974, 78, 170. Ben-Naim, A. Water and Aqueous Solutions; Plenum: New York, 1974.
Alkanes in Water: The HNC-RISM Theory
I
50.
I
;
0.
h
M
$
Y
1
I
“4
10.1
-10.
I
1
40.
0.
I
The Journal of Physical Chemistry, Vol. 92, No. 18, 1988 5259
1 I
0
A
A
I
I
I
I
1
I
I
I
I
I
r.
u
-20.
A
1
-30.
VI
a -40.
‘
-50.
0.
I b {
0 A
i
0
0 1
10.
I
1
20.
30.
T
I
I
50.
60.
I
40.
(“C)
C
I
-
0 0 o
-10.
-8>
r
‘
I -5. -
B
v
-10.
-
O
-
1
1
I
I
d
A
A
0
\
3
I
0
-
0
0
A
0
0
0
0
I
i
Lo
-15.
‘
0.
I
10.
I
20
I
I
30.
40.
1 -
50.
60.
T (“C) Figure 1. Thermodynamic data for propane in aqueous solution as a function of temperature: (a) chemical potential of solvation; (b) entropy of solvation; (c) enthalpy of solvation; (d) hydrophobic interaction parameter, 6AHI. The 0 are from the HNC-RISM approximation for Ar (eq 3a), the 0 are from the Gaussian model for A@ (eq 3b), and the A are from experiment (Wilhelm, E.; Battino, R.; Wilcock, R. J. Chem. Rev. 1977, 77, 219).
lations (i.e., the Pratt-Chandler e and u, and neglect H sites), the results are in close agreement with experiment, particularly for the H N C chemical potential (Figure 1). Structure of Solvation. Since the discrepancies between the HNC-RISM and experimental results for thermodynamic properties could lie either in the potential or in the HNC-RISM formalism itself, we compared the g,(r)’s from HNC-RISM with a molecular dynamics simulation of propane in water using the same potential parameters. The site representation for the solute is shown in Figure 2. The two equivalent end methyl groups are referred to as CH3 and the central methylene as CH2. Radial distribution functions for propane in water at 25 O C are shown in Figure 3. The HNCRISM results are shown as dashed lines. The distribution of the
Figure 2. Models of propane. Solid lines indicate (schematically) van der Waals interaction spheres. Dashed line denotes the interaction sphere for the reduced molecule. The center of that sphere is the position of the virtual site labeled X .
water 0 around the CH3 groups is very similar to that around a single-site solute, while the distribution around the CHI group has a split in the first peak due to solvent approaching from a side in which the CH, groups block closer approach. There is an unphysical nonzero density at small r for the HNC-RISM CHI-H radial distribution. The molecular dynamics results are shown as solid lines in Figure 3. Although many of the overall features are the same, there are some discrepancies. Most important is that HNC-RISM predicts radial distributions are out of phase with those of molecular dynamics. Also, the spiitting of the first peak predicted by HNC-RISM for the CHI-0 distribution is not apparent in the simulation. Since there are discrepancies between X,.(r) from HNC-RISM and molecular dynamics, we also used xss,(r)from molecular dynamics in eq 1. The results did not change significantly. Use of a PY type closure in eq 1 gave gao(r) similar to but somewhat better than the HNC closure with peak positions half-way between those of H N C and those of simulation. However, the PY gmH(r) is unphysically large (-0.4) at 0 < r < uaH.This behavior occurs because caH(r)= 0 with the PY closure, as indicated above, unless auxilliary interactions are added to exclude unphysical overlaps. [Through geometric considerations, one may augment the PY closure with auxilliary sites and distances of closest approach. That procedure would give nonzero c , ~ ( ~ ) ’ s . ] Discussion of Results for Propane in Water. Although we have chosen to focus on structural predictions and not parameter optimization, there are several indications that the latter may be necessary. The parameters for the solute-solvent interactions are obtained by combining rules for pure solute and pure solvent parameters [we use the rule given in eq 4; other combining rules, cas = I / 2 ( u a + us) and e,, = ( C , C ~ ) ~give / ~ , similar results]. However, we note that Jorgensen’s Monte Carlo study of the thermodynamics of solvation of hydrocarbons in TIPS4P water gives rather poor agreement with experiment.Ic It is possible that the combining rules are inappropriate for a nonpolar solute in polar, hydrogen-bonded water. For instance, the van der Waals radius of water has been optimized to give a good structure for pure water and is in fact much larger than the value from gasphase studies,16Le., the value used in the Pratt-Chandler studies. The gas-phase value may be more appropriate for a nonpolar solute-water interaction. Our concern herein, however, is on how well HNC-RISM treats the statistical mechanical manifestations of nonpolar solute-water interaction. The molecular dynamics results show that HNCRISM underestimates the distance of closest approach and gives some incorrect features for g&). However, HNC-RISM does give a nontrivial prediction for caH(r)even though there are no a-H interactions for this system. A nonzero caH(r)may be necessary to describe some of the ordering effects attributed to the hydrophobic effect. The large differences in the thermodynamic properties from use of the two different expressions for the chemical potential (eq (16) Hirshfelder, J. 0.; Curtiss, C . F.; Byrd, R. B. Molecular Theory of Gases and Liquids; Wiley: New York, 1954).
Ichiye and Chandler a site-site model with known geometry and parameters such as Jorgensen’s OPLS8 is replaced with one where each site now represents several sites of the original model. The new sites of the reduced model may be chosen at the geometric center or center of mass of the sites to be replace. It may correspond to the location of one of the original sites or it may be a virtual site, Le., one with no potential energy interactions with the solvent. (Such virtual sites are a special case of “auxiliary sites”3b,17,18 with the distance of closest approach set equal to zero.) The first step is to obtain g,(r), the radial distribution functions of the solvent sites s about the site v of the reduced model in the presence of all of the sites. Whether v is a real site or a virtual site, obtaining g, may be accomplished in several ways. One way is to carry out a simulation. This is the most accurate way, though perhaps the most inconvenient computationally. Alternatively, one can add a virtual site to a RISM type calculation. One condition that must hold for a rigorous theory is that for all real sites a
a
I
00
I
I
1
I
1.5
HNC-RISM,
simulation
&(r) with virtual site = has(r)without virtual site
+
0.0 0.
/
2.
I
I
I
4.
6.
8.
10
r (Angstroms)
2.0
t 0.5
I
I
2.0
LL
M
I
I
,+simulation
I
I/ I
l
I
l
1
I
I
I
I
I
1
CH,-H HNC-RISM
dl
simulation
1.0
t-
00 0
real
sites
HNC-RISM\,r
1 ‘, ~
-
I
CH,-0
I 0.0
For a PY closure with continuous potentials, c,(r) is zero for all r, so the Ornstein-Zernike-like equations (eq 1) for real sites are unchanged, as shown for the case of auxiliary sites with zero hard-sphere ~ - a d i i . ’ ~The % ’ ~radial distribution function for u is then simply a function of cms(r)[or hms(r)]:
Propane in Water
r
1
(5)
/
2
1
I
- - - +
4
I
8
I
0
10
r (Angstroms)
Figure 3. Radial distribution functions between propane sites and water sites at T = 298 K: (a) CH,-O; (b) CH,-H; (c) CH,-O; (d) CH2-H. In each figure, the dashed line corresponds to results from the HNCRISM equation, and the solid corresponds to results from molecular
dynamics. 3a and 3b) is disturbing. The H N C expression is the result of thermodynamic integration,12 while the Gaussian equation can be explained in terms of a physical modelI3 If in fact the physical model is good, the difference between the two equations may reflect the inaccuracies of the H N C closure. It is not possible to settle the issue unambiguously by using our results alone. A simulation study of the thermodynamics could provide the necessary information. Overall, these results suggest that predictions of the HNCRISM theory should be viewed with caution. IV. Reduced Models Another issue we are concerned with is the development of reduced solvent-solute potential models that correctly retain aspects of the solvation pair structure. In such a reduced model,
where OUo1is the aP element of &-I. However, for an HNC-like closure, c,(r) is not necessarily zero for all r. In a rigorous theory, however, all the c,(r)’s must behave in such a way that eq 5 holds. But eq 5 does not hold exactly for the H N C closure. The nonzero contribution to c J r ) in the H N C theory arises from so-called “unallowed” graphs.”~’* With these remarks in mind, let us now consider a ”reduced molecule” or “reduced set of sites” where the collection of real interaction sites (three for propane) are replaced by a single interaction site. The location of that site is at the point of the virtual site in the full molecule or collection of sites (see Figure 2 ) . Let H,(r) denote the pair correlation function between the center of the reduced model for the molecule (or collection of sites) and the solvent interaction sites. Our criterion for forming the solvent-solute pair potential, U,(r), between the reduced model and the solvent is to pick that potential such that (7)
This condition together with the Chandler-Andersen equation, ( l ) , and the H N C closure, (2), yields -PUYs(r)= -ffus(r) + ~ f f v s 4 x - 1 1 s ~ + s ( rIn) [ 1
+ ff,(r)l
(8)
S‘
We have studied the utility of the HNC-RISM theory in implementing these reduced models with several calculations devoted to the aqueous solvation of propane. First, we examined whether c,(r) = 0 is a useful approximation. This was done by computing the has(r)functions for propane in water via molecular dynamics as described earlier and then inserting these correlation functions into eq 6. The g J r ) that results is compared with that of simulation in Figure 4. It is seen that the approximation c,(r) = 0 produces a virtual site correlation function with main peaks that are too broad and too low. One might hope that HNC-RISM could improve these results since c,(r) # 0 in that approximation. To study this possibility, we performed two sets of HNC-RISM calculations. In one, X,(r)’s obtained from simulation were used as input; in the other, the HNC-RISM correlation functions for (17) Ladanyi, B. M.; Chandler, D. J . Chem. Phys. 1975, 62, 4308. (18) Cummings, P. T.; Gray, C. G.; Sullivan, D. E.; J . Phys. A : Math. Gen. 1981, 14, 1483. (19) Chandler, D.; McCoy, J. D.; Singer, S. J. J . Chem. Phys. 1986, 85, 5911. 5911.
The Journal of Physical Chemistry, Vol. 92, No. 18, 1988 5261
Alkanes in Water: The HNC-RISM Theory
1.0
Propane in Water 2.0 I
0.5
I
c
I
I
1
I
I
I
I
I
I
I f ! I
0.0
I
- )
1.5
-0.5
‘
0.
I
I
2.
4.
I
6.
1
8.
I 10.
r (Angstroms)
,/
1
i - _
0.0 0.
4.
2.
I
I
8.
8.
Figure 5. Potential energy function between the reduced model site u of propane and the 0 of water at T = 298 K. The solid line is from from HNC-RISM, and the dashed line is the average Lennard-Jones potential (eq 9). 10.
r (Angstroms)
Figure 4. Radial distribution function between the reduced model site v of propane and the 0 of water at 298 K. The results from HNC-RISM are shown as a dashed line, from eq 6 as a dotted line, and from molecular dynamics as a solid line.
the pure solvent pair correlations were used. The results for both calculations were nearly identical, and the latter are illustrated in Figure 4. The widths and heights of peaks in g,(r) seem accurately predicted by the HNC-RISM theory, but the HNCRISM gus(r)’sare significantly out of phase with those of simulation. These errors due to the H N C closure are perhaps disappointing, but they are not unlike those found in earlier applications of the HNC-RISM theory.’ Finally, we have used eq 8 to compute the effective pair potential in the reduced model. We use H,(r) obtained from the HNCRISM theory. Our result is shown in Figure 5. It is interesting to note that the potential so obtained is quite similar to one found by simply averaging the orientations of the propane molecule:
In this equation, the virtual site v is at the origin, the a t h site of the propane molecule is at &.)(a), the oxygen of a water molecule is at r, and J d 3 denotes the normalized integration over orientations of the propane molecule. In summary, the basic procedure for calculating effective potentials involves (1) calculating h,(r) from u,(r), ( 2 ) identifying H,(r) with h,(r) (eq 7), and (3) using a closure to obtain Uus(r) from H,(r). In light of the inaccuracies in HNC-RISM theory noted in section IV, however, a few remarks are in order: First, the H,(r) obtained from simulation are difficult to Fourier transform since there are barely two peaks. Various methods for tapering or estimating the long-range behavior could be used but would introduce uncertainties. Increasing the box size would lead to an extremely inefficient procedure. Second, the procedure is internally consistent since the HNC-RISM theory is used to obtain h,(r) from u,(r) and the U J r ) from H,(r). Finally, the primary
idea of this procedure is the identification of H,(r) with h,(r) (eq 7 ) . The procedure could be implemented with other better closures if available. Moreover, the results using the H N C closure for both steps 1 and 3 do appear physically reasonable since the minimum of Uys(r)agrees with that of Uus(r)(eq 9) and the repulsive wall of Uus(r)is less steep than that of u , ( r ) .
V. Conclusions The results presented here show that the predictions of HNCRISM theory are in only rough agreement for the distribution functions with molecular dynamics and for thermodynamic measurables with experiment. The discrepancies with experimental values for thermodynamic data may be due to errors in the distribution functions and/or the interaction parameters between the solute and water. Errors in HNC-RISM have been previously noted, particularly in the number densities and thermodynamics of ionic solution^.'^ Since the problems seemed to be exacerbated by increasing charge of the solute, the hope was that HNC-RISM would at least be able to treat nonpolar solutes in aqueous solution. However, the results presented here show that HNC-RISM is not quantitatively accurate even for this case. Other uses of the RISM theory but without the H N C closure have been more accurate. For example, studies of the structures of nonassociated and also the phenomenon of hydrophobic s o l v a t i ~ nhave ~ ~ ~been treated relatively accurately with the RISM theory. Thus, the general framework of the theory seems capable of better results than those we have found herein. It may be that improvements can be obtained with the new closures recently derived with density functional theory.Is This possibility awaits further research. Despite our uncertainties in the accuracy of the HNC-RISM theory, the application of it to develop a reduced model for alkane chains yields a sensible result that is readily interpretable in terms of simple orientational averaging. Therefore, we expect that this application yields a reliable result that can be used to simplify models employed in simulations of long chains. Acknowledgment. This research has been supported by a grant from the NIH. Registry No. H 2 0 , 7732-18-5; propane, 74-98-6.
