Free Energy, Entropy and Heat Capacity of the ... - ACS Publications

L. Renee Olano and Steven W. Rick. Journal of the American Chemical ... David Rodriguez-Gomez , Eric Darve , Andrew Pohorille. The Journal of Chemical...
11 downloads 0 Views 132KB Size
6884

J. Phys. Chem. B 2000, 104, 6884-6888

Free Energy, Entropy and Heat Capacity of the Hydrophobic Interaction as a Function of Pressure Steven W. Rick† AdVanced Biomedical Computing Center, SAIC-Frederick, NCI-Frederick Cancer Research and DeVelopment Center, Frederick, Maryland, 27102 ReceiVed: March 3, 2000; In Final Form: May 22, 2000

Molecular dynamics simulations of a methane pair in water are used to calculate the thermodynamic properties of the hydrophobic interaction as a function of pressure. Pressure is found to decrease the tendency to form methane aggregates. The entropic contribution to the free energy, which at atmospheric pressure greatly stabilizes aggregation, is highly pressure dependent. As the pressure increases, the entropic stabilization steadily decreases until, at 7 kbar, the entropy of the contact pair is equal to the entropy of the solvent separated pair. The heat capacity change between the contact and solvent separated pair is shown to be large and positive at 1 atm, as is characteristic of hydrophobic processes. At higher pressures, the heat capacity change is zero, indicating that two of the significant properties of the hydrophobic effect, the large entropy decrease and heat capacity increase are lost at high pressures. The free energy, volume and entropy changes are consistent with the corresponding changes for the pressure denaturation of proteins.

I. Introduction

(A‚A)aq h 2(A)aq

The hydrophobic interaction, the tendency of nonpolar species to aggregate in water, is believed to be one of the major forces stabilizing proteins.1,2 Methane pairs in water provide a good computational model for studying the details of the hydrophobic interaction. Under the solvent conditions in which proteins are stable in the native state (T ) 298 K, P ) 1 atm, no salt or other denaturants), computational studies have shown that the methane pairs have a tendency to aggregate, with the pairs in direct contact about 1 kcal/mol more stable than the isolated pair.3-20 Aggregation is found to be entropically favorable and enthalpically unfavorable and is a way to minimize the entropic penalty, since adding nonpolar solutes to water increases the entropy.10,14,18,19 In the thermodynamics of protein folding, the entropic gain from the hydrophobic interaction partially offsets the entropy loss in forming the native protein.21 Under conditions in which proteins denature, either the hydrophobic interaction must change or another force in protein stability must become significant. How the hydrophobic interaction is changed can be examined by calculating the free energy of association under increasingly denaturing conditions, namely the application of pressure. Protein denaturation,

native h denatured

(1)

and the hydrophobic effect, as measured either as the transfer of a nonpolar molecule, A, from a liquid hydrocarbon environment to water,

(A)H h (A)aq

(2)

or the dissociation of two nonpolar molecules in water, † Present address: Department of Chemistry, University of New Orleans, New Orleans, LA 70148, and Chemistry Department, Southern University of New Orleans, New Orleans, LA 70126.

(3)

demonstrate similar trends for the temperature dependence.21 The three processes have been written so that the right side of the equation has more exposure between nonpolar groups and water. The ability of different salts to effect the solubility and aggregation of hydrophobic molecules also correlates with the order of the salts’ denaturing capabilities (the Hofmeister series).2,22 The salt and temperature dependences establish a solid empirical relationship between protein stability and the hydrophobic effect, which the pressure dependence appears to contradict.23 Pressures of about 5 kbar can cause some proteins to denature24-29 and the free energy of the transfer process for methane decreases with pressure, indicating that the hydrophobic effect decreases with pressure.1 Pressure denaturation may be due, at least in part, to a weakened hydrophobic interaction.30,31 Simulations of methane in water find that the tendency to form aggregates decreases with pressure.31-34 However, one simulation of methane pairs in solution finds the hydrophobic interaction is stronger at 5 kbar than at 1 atm.30 In other respects, the pressure dependence of the hydrophobic interaction and protein denaturation are not so similar.23 The volume and heat capacity change for the hydrophobic interaction and protein denaturation have a different pressure dependence. At low pressures, the volume change, ∆V, of protein unfolding is small and positive while at pressures above 1-2 kbar ∆V becomes large and negative.23,35 Because V ) (∂G/∂P)T,N, a positive ∆V at low pressures implies that moderate pressure can stabilize the folded state and this is true for some proteins (see ref 36 and references therein). For the hydrophobic transfer process, the opposite is true. At atmospheric pressure, the transfer of nonpolar molecules from organic solvent to water has a negative ∆V,1,37,38 while at pressures above 2 kbar, ∆V is positive.24,26 For methane pairs, the volume change for the contact pair relative to the separated pair might be expected to follow the same trend as the ∆V for the transfer process but the results are mixed. In the simulation of Payne et al., ∆V is positive at low pressures and negative above 5 kbar.30 The ∆V

10.1021/jp000841s CCC: $19.00 © 2000 American Chemical Society Published on Web 06/29/2000

Hydrophobic Interaction as a Function of Pressure

J. Phys. Chem. B, Vol. 104, No. 29, 2000 6885

from the study of Hummer et al., has the opposite behavior, negative at low pressures and positive at pressures above 5 kbar.31 The heat capacity change, ∆CP, upon denaturation is insensitive to pressure and is large and positive whereas the ∆CP for the transfer of nonpolar solutes from organic solvents to water is large and positive at 1 atm and zero around 5 kbar.26 To examine the thermodynamics of the hydrophobic interaction in further detail, previous work done at atmospheric pressure using a polarizable fluctuating charge potential for water and methane will be extended to higher pressures.19 The potential of mean force between two methane molecules at high pressures (3, 5 and 7 kbar) and also the entropy, heat capacity and volume changes of association, will be calculated. II. Methods The methane-water system is simulated using the fluctuating charge (FQ) method for treating electronic polarization, which has been applied to many systems.19,39-42 In this model the partial charges on atomic sites respond to changes in their environments and flow between sites on a molecule to minimize the energy. The FQ model for water, which uses a rigid geometry based on the TIP4P water model,43 accurately reproduces many of the equilibrium and dynamical properties of water at 1 atm and 298 K.40 The methane molecules are treated with the FQ model as well, also with rigid bond lengths and angles.19 The methane potential was constructed to give charge values for the isolated molecule in good agreement with electronic structure calculation and to give the correct solvation free energy. The water-methane dimer energy also agrees well with the electronic structure calculations at the fourth-order Møller-Plesset (MP4) level44 and the methane-methane dimer energy is close to the Jorgensen-Madura-Swenson LennardJones energy, which has been shown to be an accurate potential for methane.45 The simulations were done in the isothermal-isobaric (constant T,P,N) ensemble, by coupling to a pressure bath and a Nose´-Hoover temperature bath.46-50 The simulations were done at pressures of 3, 5 and 7 kbar and at a temperature of 298 K, except for additional simulations at 5 kbar which were done at 283 and 313 K in order to calculate the entropy and the heat capacity. All simulations used a 1 fs times step, 256 solvent molecules, the Ewald sum for long-ranged electrostatic interactions and bond constraints were enforced using the SHAKE algorithm.51 The use of Ewald sums has been shown to be particularly important for simulating polarizable systems.52 The potential of mean force calculations were done using umbrella sampling, with a quartic restraining potential, 1/4K(rCC - r0)4 on the methane carbon distance, rCC, with 6 windows centered an angstrom apart from r0 ) 4 to 9 Å with a force constant, K, equal to 40 (kcal/mol)/Å.4 Each window was simulated for 5.0 ns and at the higher temperature, 313 K, each window was simulated for 3.0 ns. An additional set of simulations at 298 K and 1 atm were done to examine the volume change. The potential of mean force was calculated from the biased data using the weighted histogram method.53,54 Error estimates are calculated using a block averaging method in which the data is split into one nanosecond blocks. Uncertainties are reported as (2σ. III. Results The methane pair potential of mean force, w(r), is shown in Figure 1, which compares w(r) at 3, 5 and 7 kbar with the previous result at 1 atm.19 All the curves have a contact pair minimum, around 3.8 Å, a solvent separated pair minimum, around 6.5 to 6.8 Å, and the beginning of another minimum at

Figure 1. Methane pair potential of mean force at P ) 1 atm (solid line), 3 kbar (dashed line), 5 kbar (dotted line) and 7 kbar (dot-dashed line).

TABLE 1: Free Energy Difference between the Contact Pair (CP) and Solvent Separated Pair (SSP), CP and Barrier (B), and CP and the Value at 8 Å pressure

w(rCP) - w(rSSP) (kcal/mol)

w(rCP) - w(rB) (kcal/mol)

w(rCP) - w(8 Å) (kcal/mol)

1 atm 3 kbar 5 kbar 7 kbar

-0.68 -0.51 ( 0.02 -0.44 ( 0.03 -0.46 ( 0.04

-0.81 -0.87 ( 0.04 -0.97 ( 0.05 -0.93 ( 0.02

-1.04 -0.96 ( 0.02 -1.07 ( 0.06 -1.06 ( 0.01

9 Å. All curves are shown relative to their values at 8 Å. The w(r) at 1 atm was only calculated out to 8.5 Å. The high-pressure calculations used an additional umbrella sampling window centered at 9 Å. Structure in the w(r) beyond 8 Å has been seen in other studies at 1 atm.17,18 Calculating w(r) beyond 10 Å would require a bigger simulation box with more molecules, since this would be greater than half the current box length. The value of the contact pair (CP) minimum is relatively insensitive to pressure, while the solvent separated pair (SSP) becomes more stable with pressure (Table 1). The barrier height between the CP and the SSP increases and moves to smaller distances as the pressure is increased. The SSP minimum decrease with pressure, reaching a minimum at 5 kbar, was reported in previous studies.30,31 However, the study of Payne et al., found that the contact pair minimum, relative to w(r) at 8 Å, also decreases with pressure, and w(rCP) - w(rSSP) appears to be not strongly pressure dependent.30 The present results for w(rCP) - w(rSSP) are more consistent with the results of Hummer et al.31 The entropy of association can be calculated directly as the difference between the free energy and enthalpy

-T S(r) ) w(r) - 〈U(r)〉 - P〈V(r)〉

(4)

where 〈‚‚‚〉 denotes the ensemble average, U(r) is the potential energy, P is the pressure and V(r) is the volume.10 The entropy can also be calculated through a finite difference temperature derivative of the free energy

-S(r) )

w(r,T+∆T) - w(r,T-∆T) 2∆T

(5)

which involves calculating w(r) at two additional temperatures but can be more efficient numerically than the direct method of eq 4.10 Figure 2A shows the entropy change from ref 19, calculated from eq 4, and the high-pressure results at 3 and 7 kbar, calculated from eq 5. Figure 2B compares the results at 5 kbar calculated from eqs 4 and 5. The entropies are given

6886 J. Phys. Chem. B, Vol. 104, No. 29, 2000

Rick

Figure 3. Heat capacity as a function of methane separation at 1 atm (solid line) and 5 kbar (dotted and dashed line). The solid and dotted lines were calculated from the numerical second derivative of the free energy (eq 6) and the dashed line was calculated from the numerical first derivative of the enthalpy (eq 7).

to find CP(r) at 1 atm (Figure 3 and Table 2). The heat capacity can also be found from

CP(r) )

Figure 2. Entropic contribution, -TS(r), to the potential of mean force, at (A) 1 atm (solid line),19 3 kbar (dotted line) and 7 kbar (dashed line) and (B) 5 kbar, calculated by eq 5 (solid line) and eq 4 (dotted line).

TABLE 2: Entropy and Heat Capacity Differencesa CP(rCP) CP(rCP) -T(S(rCP) - -T(S(rCP) S(rSSP)) S(8 Å)) CP(rSSP) CP(8 Å) (kcal/mol) [(kcal/mol)/K] [(kcal/mol)/K] pressure (kcal/mol) 1 atm 3 kbar 5 kbar 7 kbar

-4.0a -3.2 ( 0.9b -1.4 ( 0.7a -2.2 ( 2.3b 0.1 ( 1.0 b

-3.2a -0.52c -0.68c -2.1 ( 0.8b -0.2 ( 0.4a 0.37 ( 0.12c 0.29 ( 0.17c -1.2 ( 0.9b -0.07 ( 0.16d -0.04 ( 0.07d 0.3 ( 0.8b

a a calculated by eq 5; b calculated by eq 4; c calculated by eq 6; and d calculated by eq 7.

relative to their values at 8 Å. The entropies from the different methods are in good agreement and provide an independent estimate of the entropy at this pressure. The differences between the two methods are comparable to the statistical uncertainties (see Table 2). As can be seen in Figure 2 and Table 2, the entropy change shows a strong dependence on pressure. At 1 atm, entropy greatly stabilizes the contact pair. At 7 kbar, the entropy change is roughly zero. In all cases the entropy decreases substantially at the barrier between the contact pair and the solvent separated pair, which moves from around 6 Å at 1 atm to 5 Å at 7 kbar. The heat capacity change for the process can be calculated using a finite difference second derivative of the free energy,

( )

CP(r) ) -T

∂2w(r) ∂T2



P,N

w(r,T+∆T) - 2w(r,T) + w(r,T-∆T)

-T

(∆T)2

(6)

The w(r) at three different temperatures from ref 19 were used

( ) ∂H(r) ∂T

P,N



H(r,T+∆T) - H(r,T-∆T) 2∆T

(7)

where H(r) is the enthalpy (〈U(r)〉 - P〈V(r)〉). A third expression for CP involving the fluctuations in the enthalpy (CP(r) ) 〈∆H(r)2〉kT2) was not used because of the large numerical uncertainties in 〈∆H(r)2〉.51 At 5 kbar, CP was calculated using both eqs 6 and 7. At 1 atm, only eq 6 was used since the enthalpy at different temperatures was not available. The heat capacity change is large and negative at 1 atm, in agreement with other hydrophobic processes.2 At 5 kbar, the heat capacity change is smaller and zero or slightly positive, depending on what method is used. A heat capacity change of about zero agrees with the heat capacity change for the transfer of nonpolar solutes from organic solvents to water at 5 kbar.26 The value of the heat capacity change at 1 atm is highly approximate, due to a large degree of uncertainty from the numerical second derivative. The results at 5 kbar were calculated from simulations more than twice as long than the previous simulations at 1 atm,19 so the free energies should be more accurate. In addition, the level of agreement between the heat capacities calculated the two different ways lends validity to the qualitative conclusions regarding the pressure dependence of the heat capacity. The average volume, 〈V(r)〉, as a function of solute separation, relative to the value at 8 Å, is shown in Figure 4A,B. At 1 atm, the volume of the contact pair is negative, in agreement with the conclusion of Payne et al.30 At higher pressures, the volume of the solvent separated pair is about equal to the contact pair volume and at 5 kbar may be slightly smaller. The volume changes are small, as is indicated from the pressure dependence of w(r). The volume difference, ∆V, between the CP and the SSP can be estimated from ∆w ()w(rCP) - w(rSSP)) at different pressures. At 5 kbar, ∆V can be estimated from (∆w(P)7 kbar) - ∆w(P)3 atm))/(4 kbar), which, from the data in Table 1, is equal to 0.5 mL/mol, less than the numerical uncertainty in 〈V(r)〉. The volume changes are small enough so that the contribution of the volume term, P〈V(r)〉, to the enthalpy is small at all pressures. IV. Conclusion The free energy of methane association in water has been studied as a function of pressure. Pressure stabilizes the solvent

Hydrophobic Interaction as a Function of Pressure

Figure 4. Average volume as a function of solute separation at A) 1 atm (solid line) and 3 kbar (dotted line) and B) 5 kbar (solid line) and 7 kbar (dotted line).

separated pair relative to the contact pair, in support of a mechanism of pressure protein denaturation involving incorporation of water molecules into the protein.31 Entropy, which plays a large role in the hydrophobic interaction and greatly stabilizes the contact pair at 1 atm,10,14,18,19 shows a strong pressure dependence (Figure 2A,B). At 7 kbar, the entropy change is zero and the hydrophobic interaction changes from entropically dominated to entropically indifferent. The heat capacity change was shown to be large and negative at 1 atm, which is characteristic of many hydrophobic processes.2 At 5 kbar, the heat capacity change is no longer negative (Figure 3). Thus, two important characteristics of the hydrophobic effect, the large entropy decrease and the large heat capacity increase, are lost at pressures around 5 kbar. The entropy and heat capacity changes can be related to the structure of pure water, which also show a strong pressure dependence.55-57 The structure of pure water shows a large temperature dependence at 1 atm, but the effects of temperature become smaller at pressures of 1 to 6 kbar.57 Therefore, many of the properties of water should be less temperature dependent at higher pressures and the heat capacity and entropy should be larger at 1 atm than at higher pressures. The tetrahedral structure, which gives water its unique hydrogen bonding arrangement with a high entropy, is lost at high pressure.55-57 In conventional views of solvation, the addition of solutes disrupts this hydrogen bond arrangement and leads to a decrease in the entropy.1,2,58 At higher pressures, the hydrogen bond arrangement of the pure solvent is already disrupted and the addition of solutes does not lead to a large entropy change. The results presented here help to clarify the relationship between the pressure dependence of protein denaturation (eq 1), the transfer of nonpolar solutes from a liquid hydrocarbon

J. Phys. Chem. B, Vol. 104, No. 29, 2000 6887 to an aqueous environment (eq 2) and hydrophobic dissociation (eq 3). For all processes (protein denaturation,24-26,28,29 the transfer process,1 and the present and previous results for aqueous methane31-34) the application of pressure around 5 kbar drives the equilibrium of eqs 1-3 toward the right. A volume change, ∆V, at 1 atm for methane pair dissociation is found to be positive. The volume change for protein denaturation is also positive23,35 but the ∆V for the transfer process is negative.1,37,38 At higher pressure, ∆V is negative for protein denaturation23,35 and positive for the transfer process.24,26 At 5 kbar, a zero or slightly negative ∆V is found (for the contact pair to solvent separated pair transition) and so, both at low and high pressures, the ∆V results are more consistent with the volume changes for protein denaturation than the hydrophobic transfer process. The positive ∆V at 1 atm implies that the hydrophobic effect becomes stronger at pressures just above 1 atm1,30 and some proteins are stabilized by moderate pressure increases.36 The results presented here for the entropy change, which becomes zero at 7 kbar, can be compared to the entropy of denaturation. The slope of the pressure-temperature protein stability curve becomes zero around 300 K and 4 kbar, depending on the protein and the pH.26,28 From the Clapeyron equation, (dP/dT)∆G)0 ) ∆S/∆V, the zero slope implies a zero entropy change. The ∆S for protein denaturation, at 1 atm, is due not just to the hydrophobic contribution but to a number of factors, including the large gain in conformational entropy, which is why ∆S is positive. Unlike temperature denaturation, proteins denatured by pressure retain some native structure27,59 so the conformational entropy changes will not be as large. The biggest structural change for pressure denaturation may be the incorporation of water molecules in protein interiors. On the basis of the results presented here, the ∆S for this should be near zero, which is consistent with the small ∆S for protein denaturation at high pressures. Acknowledgment. This project has been funded in whole or in part with Federal funds from the National Cancer Institute, National Institutes of Health, under Contract No. NO1-CO56000. References and Notes (1) Ben-Naim, A. Hydrophobic Interactions. Plenum: New York, 1980. (2) Dill, K. A. Biochemistry 1990, 29, 7133. (3) Pratt, L. R.; Chandler, D. J. Chem. Phys. 1977, 67, 3683. (4) Pangali, C.; Rao, M.; Berne, B. J. J. Chem. Phys. 1979, 71, 2975. (5) Swaminathan, S.; Beveridge, D. L. J. Am. Chem. Soc. 1979, 101, 5832. (6) Ravishanker, G.; Mezei, M.; Beveridge, D. L. Faraday Symp. Chem Soc. 1982, 17, 79. (7) Rapaport, D. C.; Scheraga, H. A. J. Phys. Chem. 1982, 86, 873. (8) Watanabe, K.; Andersen, H. C. J. Phys. Chem. 1986, 90, 795. (9) Jorgensen, W. L.; Buckner, J. K.; Boudon, S.; Tirado-Rives, J. J. Chem. Phys. 1988, 89, 3742. (10) Smith, D. E.; Zhang, L.; Haymet, A. D. J. J. Am. Chem. Soc. 1992, 114, 5875. (11) Smith, D. E.; Haymet, A. D. J. J. Chem. Phys. 1993, 98, 6445. (12) van Belle, D.; Wodak, S. J. J. Am. Chem. Soc. 1993, 115, 647. (13) Dang, L. X. J. Chem. Phys. 1994, 100, 9032. (14) Jungwirth, P.; Zahradnik, R. Chem. Phys. Lett. 1994, 217, 319. (15) Head-Gordon, T. Chem. Phys. Lett. 1994, 227, 215. (16) New, M. H.; Berne, B. J. J. Am. Chem. Soc. 1995, 117, 7172. (17) Lu¨demann, S.; Schreiber, H.; Abseher, R.; Steinhauser, O. J. Chem. Phys. 1996, 104, 286. (18) Lu¨demann, S.; Abseher, R.; Schreiber, H.; Steinhauser, O. J. Am. Chem. Soc. 1997, 119, 4206. (19) Rick, S. W.; Berne, B. J. J. Phys. Chem. B 1997, 101, 10488. (20) Young, W. S.; Brooks, C. L., III. J. Chem. Phys. 1997, 106, 9265. (21) Baldwin, R. L. Proc. Natl. Acad. Sci. U.S.A. 1986, 83, 8069. (22) Sacco, A.; De Cillis, F. M.; Holz, M. J. Chem. Soc., Faraday Trans. 1998, 94, 2089.

6888 J. Phys. Chem. B, Vol. 104, No. 29, 2000 (23) Kauzmann, W. Nature 1987, 325, 763. (24) Brandts, J. F.; Oliveria, R. J.; Westort, C. Biochemistry 1970, 9, 1038. (25) Hawley, S. A. Biochemistry 1971, 10, 2436. (26) Zipp, A.; Kauzmann, W. Biochemistry 1973, 12, 4217. (27) Weber, G.; Drickamer, H. G. Q. ReV. Biophys. 1983, 16, 89. (28) Gross, M.; Jaenicke, R. Eur. J. Biochem. 1994, 221, 617. (29) Jonas, J.; Jonas, A. Annu. ReV. Biomol. Struct. 1994, 23, 287. (30) Payne, V. A.; Matubayashi, N.; Murphy, L. R.; Levy, R. M. J. Phys. Chem. B 1997, 101, 2054. (31) Hummer, G.; Garde, S.; Garcı´a, A. E.; Paulaitis, M. E.; Pratt, L. R. Proc. Natl. Acad. Sci. 1998, 95, 1552. (32) Wallqvist, A. Chem. Phys. Lett. 1991, 182, 237. (33) Wallqvist, A. J. Phys. Chem. 1991, 95, 8921. (34) Wallqvist, A. J. Chem. Phys. 1992, 96, 1655. (35) Harpaz, Y.; Gerstein, M.; Chothia, C. Structure 1994, 2, 641. (36) Sun, M. M. C.; Tolliday, N.; Vetriani, C.; Robb, F. T.; Clark, D. S. Protein Sci. 1999, 8, 1056. (37) Gross, D. S.; Schnier, P. D.; Rodriguez-Cruz, S. E.; Fagerquist, C. K.; Williams, E. R. Proc. Natl. Acad. Sci. U.S.A. 1996, 267, 1483. (38) Scheraga, H. A. ACS Symp. Ser. 1994, 568, 641. (39) Rick, S. W.; Berne, B. J. J. Am. Chem. Soc. 1996, 118, 672. (40) Rick, S. W.; Stuart, S. J.; Berne, B. J. J. Chem. Phys. 1994, 101, 6141. (41) Liu, Y.-P.; Kim, K.; Berne, B. J.; Friesner, R. A.; Rick, S. W. J. Chem. Phys. 1998, 108, 4739. (42) Banks, J. L.; Kaminski, G. A.; Zhou, R.; Mainz, D. T.; Berne, B. J.; Friesner, R. A. J. Chem. Phys. 1999, 110, 741.

Rick (43) Jorgensen, W.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (44) Szcze¸ s´niak, M. M.; Chałasin´ski, G.; Cybulski, S. M.; Cieplak, P. J. Chem. Phys. 1993, 98, 3078. (45) Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. J. Am. Chem. Soc. 1984, 106, 6638. (46) Andersen, H. C. J. Chem. Phys. 1980, 72, 2384. (47) Ciccotti, G.; Ryckaert, J. P. Comput. Phys. Rep. 1986, 4, 345. (48) Martyna, G. J.; Tobias, D. J.; Klein, M. L. J. Chem. Phys. 1994, 101, 4177. (49) Nose´, S. Mol. Phys. 1984, 52, 255. (50) Hoover, W. G. Phys. ReV. A 1985, 31, 1695. (51) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, U.K., 1987. (52) Rick, S. W. In Simulation and Theory of Electrostatic Interactions in Solution; Pratt, L. R.; Hummer, G., Ed.; American Institute of Physics: Melville, NY, 1999; pp 114-126. (53) Ferrenberg, A. M.; Swendsen, R. H. Phys. ReV. Lett. 1989, 63, 1195. (54) Kumar, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; Rosenberg, J. M. J. Comput. Chem. 1992, 13, 1011. (55) Wu, A. Y.; Whalley, E.; Dolling, G. Mol. Phys. 1982, 47, 603. (56) Okhulkov, A. V.; Demianets, Y. N.; Gorbaty, Y. E. J. Chem. Phys. 1994, 100, 1578. (57) Bellisent-Funel, M.; Bosio, L. J. Chem. Phys. 1995, 102, 3727. (58) Stillinger, F. H. Science 1980, 209, 451. (59) Zhang, J.; Peng, X.; Jonas, A.; Jonas, J. Biochemistry 1995, 34, 8631.