J. Phys. Chem. 1995,99, 14632-14640
14632
Temperature Effects on the Hydrophobic Hydration of Ethane R. L. Mancera and A. D. Buckingham* Department of Chemistry, University of Cambridge, Lensfield Road, Cambridge CB2 IEW, U.K. Received: January 3, 1995; In Final Form: April 26, 1995@ A series of molecular dynamics simulations have been performed on a system comprising 256 water molecules and 1 ethane molecule at 279, 302, 320, and 352 K. The structural analysis revealed the existence of orientational restrictions for water molecules with respect to the solute in the hydration shell, while a hydrogen bond analysis showed a slight enhancement of the tetrahedral network. This behavior is consistent with the negative entropy and enthalpy functions associated with the hydrophobic hydration of nonpolar solutes near room temperature. The dynamical analysis shows that hydration-shell water molecules have slightly lower diffusion coefficients; on the other hand, calculations of hydrogen bond lifetimes reveal that this region's bonding network breaks up faster than that of the bulk. The overall results suggest that the energetic environment in the hydration shell is similar to that in the bulk and reinforce the idea of a reduction in the number of low-energy paths for molecular orientational rearrangements in the vicinity of nonpolar solutes.
1. Introduction Hydrophobic effects'-6 are believed to play a major role in many fundamental processes in biology and chemistry, like the folding and stability of proteins and nucleic acids, the conformation and association equilibria of small molecules in water, the formation of micelles and biomembranes, and the low solubility of hydrocarbons and inert gases in water. Hydrophobic effects are usually divided into two categories: hydrophobic hydration and hydrophobic interaction. Hydrophobic hydration refers to the thermodynamic, structural, and dynamic properties of very dilute aqueous solutions of nonpolar solutes. Hydrophobic interaction refers to the solvent-induced forces that drive nonpolar solutes together when placed in an aqueous environment. The low solubility of nonpolar molecules in water results from the fact that the transfer of these solute molecules from the gas phase to a dilute aqueous solution involves a large increase in the standard Gibbs free energy (AGO). AGO is actually the combination of a small decrease in the standard enthalpy (AHo) and a large decrease in the standard entropy (ASo) near room along with a negative temperature, as first observed by volume of solution and an increase in the heat capacity of s ~ l u t i o n . ~ -These ' ~ properties are remarkable if we think that the large heat of vaporization and surface tension of water indicate that it has rather a strong cohesive energy so that the creation of a cavity to accommodate a solute molecule would require a considerable input of energy (enthalpy). Although some of this energy should be recovered when the solute is subsequently placed in the cavity, it could be argued that the whole process involves replacing the strong water-water attractive forces by the much weaker water-solute ones. This would imply a net positive enthalpy change, which in tum could explain the observed positive Gibbs free energy. However, for many nonpolar solutes, the solubility decreases when the temperature is raised near 25 "C7with a decrease in the rate of change in the Gibbs free energy so that A W tums out to be negatiue according to the following equation: d(AG"1T)IdT = -AHolP The solubility often passes through a minimum not far above room temperature, - I 4 where AHo vanishes, gradually taking 'Abstract published in Advance ACS Absrrucrs, July 1, 1995.
0022-365419512099-14632$09.00/0
increasing positive values on further heating. In view of this situation, the relation AGO = AH" - TAS" shows that at 25 "C the transfer of the solute into the solution must involve a large negative entropy change. Most of the structural t h e o r i e ~ l - ~proposed .~ to account for the above-mentioned thermodynamic observations rely on the idea that a nonpolar solute increases the degree of hydrogen bonding ("icelikeness") of water molecules in the solvation shell so as to maximize the water-water attraction and avoid the solute's hard core, with a negative contribution to the enthalpy and a substantial decrease in the entropy of the water. This negative entropy change should then come as a consequence of some ordering or restructuring of the water molecules in the first hydration shell of the solute in order for them to maintain their hydrogen bonds. This ordering has been interpreted in terms of restrictions in the orientational freedom of the water molecule^^^ and of slowed-down dynamics;I6 it is considered later. This negative entropy change is believed to give rise to the decrease in solubility with increasing temperature. A statistical mechanical treatment of this model was developed by Nkmethy and Scheraga." It is known that nonpolar gases form crystalline clathrate hydrate^'^,'^ at high pressure and at temperatures, usually below the melting point of ice, in which each guest gas molecule is surrounded by a polyhedral cage in the lattice of host water molecules. It has been argued'-3 that the solvation shell of nonpolar solutes in water should resemble these cages at least in part, since there would be a tendency for the hydrogenbonding groups of the solvent water molecules to avoid being directed toward the solute, with each of these water molecules straddling a portion of the convex surface of the solute to form hydrogen bonds with three other water molecules in the shell. Such an orientation allows the number of strong water-water hydrogen bonds in this region to be close to, if not larger than, that in bulk water. The total binding energy of a water molecule in the solvation layer would therefore be typically at least as large as in the bulk. It is believed that the mobility of water in this hydrophobic hydration shell should be slightly reduced in comparison to bulk water, in part because of the reduced number of hydrogen-bonding rearrangements available. Another picture of hydrophobic hydration comes from the directional, hydrogen-bonding properties of water. As the hydration shell adjusts to accommodate the nonpolar solute, water molecules rearrange around the solute in order to maintain 0 1995 American Chemical Society
Hydrophobic Hydration of Ethane
J. Phys. Chem., Vol. 99, No. 40, 1995 14633
energies than those in the bulk liquid, and the bond energy the maximum possible number of hydrogen bonds. No suggestion is made that hydration water need be different from bulk distribution revealed a slight increase of hydrogen bonding in water in terms of an enhanced interaction but just that its the first hydration-shell water molecules. The translational and orientational freedom is restricted due to the inability of the rotational motions of the hydration-shell water molecules were solute to participate in hydrogen bonds with water, which makes slower (by at least 20%) than those in pure bulk water. Similar certain orientations of the hydration water molecules energetiresults were obtained by Rossky and K a r p l u ~when ~ ~ analyzing cally unfavorable. the water molecules near the nonpolar groups in a simulation of a dipeptide in water and by Rossky and ZickP in a simulation Hydrophobic i n t e r a ~ t i o n ~ Ocan - ~ ~be understood as coming of two Lennard-Jones particles in water situated at a distance from a thermodynamic force in aqueous solution that drives corresponding to a solvent-separated configuration. Rapaport hydrocarbon molecules (and other “hydrophobic” species) and Scheraga22 arrived at the same conclusions in their together to a greater extent than if these solutes were dissolved simulations of four nonpolar particles in water and found that in a more accommodating (nonpolar) solvent, since the surface within the frst solvation shell there was a tendency for the water area they would present to the solvent as a cluster is less than molecules to be oriented so that their 0-H bonds were almost that when they are separate. The entropy change on bringing tangential to an imaginary sphere surrounding the solute, as two solute molecules together should be positive and should opposed to the near-tetrahedral organization of the neighbors thereby lower the free energy of solution. In other words, the of a molecule in pure water, the results being consistent with a association of nonpolar solutes can be regarded as a partial clathrate-type picture. They observed a small degree of reversal of the thermodynamically unfavorable process of slowing down of the molecular motion in the hydration shells solution. -6 as measured by the slow molecular exchange rate between the Computer simulation has emerged as the principal tool to shell and bulk components of the solvent when compared to examine the microscopic properties of aqueous solutions of that in reference water. Guillot et aL4I found in their simulation nonpolar solutes, since their low solubility in water makes it of methane in water a clathrate-like arrangement of the difficult to obtain information from experiment. Much work solvation-shell water molecules. A simulation of methane in a has been done on the computer simulation of hydrophobic polarizable model of water by van Belle and Wodak30revealed h y d r a t i ~ n , ~and ~ - we ~ ~now review the main conclusions. a structuring of the hydration-shell water molecules and slowing An early Monte Carlo (MC) simulation study by Owicki and of their translational and orientational migration. Vaisman et S ~ h e r a g aof~ a~ solution of 1 methane in 100 water molecules analyzed the distance dependence of water structure around revealed that the interaction energy of a pair of water molecules model solutes and, in the case of methane, claimed that hydrogen in the first solvation shell was more negative and more sharply bonds in the first hydration shell are slightly longer and more distributed than that in the bulk liquid, supporting the idea that bent than in bulk water, with just a slight decrease in the number nonpolar solutes increase the degree of hydrogen bonding or of hydrogen bonds per water molecule. The methane coordination structure in their hydration shell. Recent MD test-particle insertion simulations by Guillot and number they computed (-23) is similar to the number of water G u i s ~ a n probed i ~ ~ the temperature dependence of the hydromolecules in the polyhedral cage of the methane phobic hydration of inert gases and methane in SPCE water, Also found is some preference for the first hydration-shell water building solubility saturation curves of these molecules over hydrogens to be closer than the oxygens to the center of the the temperature range 273 to 640 K. Increasing the temperature methane, while the situation is reversed in the second shell, had a strong influence on the solute-water distribution functions consistent with hydrogen bonding between the shells. Similar results were obtained in MC studies by Swaminathan et ~ 2 1 . ~ ~ for all the solutes studied. The amplitude of the first peak of the solute-oxygen radial distribution function (RDF) diminished and Bolis and Clementi,38 who found a slightly distorted and when the temperature increased, although its position did not defective clathrate structure around the nonpolar solute but with shift. At the same time, the coordination number decreased a coordination number (-19-20) closer to that expected for a almost linearly with the water density as long as the temperature pentagonal dodecahedral cage. Alagona and Tani3’ found did not exceed approximately 570 K. They concluded that water clathrate-like structures in the solvation shell of argon. Pangali molecules in the first hydration shell, which are oriented in a MC simulation of a system of two Lennard-Jones et tangentially to the solute at room temperature, began to reorient particles in water, found that hydration-shell water molecules toward the bulk above the boiling point, a clear manifestation engage in stronger hydrogen bonding than bulk molecules and of the “melting” of the icelike structure of water around the pointed out that there is one less water-water interaction at nonpolar solutes. They argued that the different degrees of short range for the shell molecules. In a MC simulation of one solubility of inert gases in water are governed by the soluten-butane molecule in water, J ~ r g e n s e nfound ~ ~ that water. water interaction energy. The bigger the solute, the greater its molecules in all regions of the solution had the same bonding polarizability and the more negative the value of the soluteenergy and hydrogen profiles (numbers and angles) as those of water interaction energy, which tends to increase the solubility. bulk water; however, the coordination numbers (which are On the other hand, the larger the solute, the more negative the somewhat less than 4)were slightly lower in the first hydration entropy of cavity formation, which causes a decrease in the shell. solubility and the minimum solubility exhibited by the solutes. Molecular dynamics (MD) simulations have also shown The energetic term, then, promotes solubility of large solutes, definite orientational preferences reminiscent of structure in while the entropic term disfavors it. clathrates for the water molecules in the first hydration layer There have been some interesting recent experimental findings around a nonpolar solute. This is what Geiger et ~ 2 1 . found ’~ in by means of neutron scattering and NMR regarding the structure a simulation of two Lennard-Jones particles in water. They of aqueous nonpolar solutions. Broadbent and N e i l ~ o n ~ ~ found that the water molecules in the first hydration shell were performed neutron diffraction experiments to study the hydration oriented in such a way that one of the four tetrahedral bond structure around argon in water at elevated pressure and found directions pointed away from the center of the solute. The a relatively strong first hydration shell with only weak longerinteraction energy distributions of these water molecules were range correlations, as revealed by the argon-water radial more sharply distributed and slightly shifted toward lower
Mancera and Buckingham
14634 J. Phys. Chem., Vol. 99, No. 40, 1995 distribution function. Soper et a1.45-47performed neutron scattering experiments to calculate pair correlation functions in tetraalkylammoniumand alcohol-water solutions and found the existence of a disordered cage structure around the methyl group(s) of these molecules, though the local tetrahedral coordination found in pure water remained. Although they found a preference for 0-H bonds to lie approximately tangential to circumscribing spheres around the nonpolar groups, they claim there is no structural or dynamical evidence of an "icelike" structure. Apparently, there is little difference in terms of water structure in going from pure bulk water to water in these systems. They found just a small increase in structure in the tetraalkylammonium solution compared to those of both pure water and the alcohol-water mixtures but little difference between the structures of the alcohol and tetraalkylammonium solutions. These last results support the idea that water molecules in the hydration shell are stabilized with respect to the bulk but suggest that the structural effect of hydrophobic hydration at room temperature is less important than commonly accepted. In an NMR and neutron-scattering study by Brad1 et aL4* of an undercooled solution of tetramethylammonium bromide, the distortion of the hydrogen-bonded water network by the dissolved cations was linked to complex changes in the dynamics of the hydration-shell water molecules. It was found that changes in the mobility of the solvent relative to that of pure water (retardation at high temperatures but facilitation at low temperatures) depend strongly on the state of the hydrationshell hydrogen-bond network; the dynamic picture is not consistent with the picture of a clathrate-like organization of the solvent around a nonpolar solute, and this was corroborated by a highly transient and short-lived hydration cage even in the undercooled phase. Encouraged by these reports and by our studies on the hydrophobic aggregation of ethane,32we performed a series of molecular dynamics simulations of the hydrophobic hydration of ethane to study the structural and dynamic properties of this system at different temperatures. In section 2 we explain briefly the computational details and intermolecular potentials used in our simulations. We then proceed in section 3 to present and discuss our results, and we make some concluding remarks in section 4.
2. Intermolecular Potentials and Computational Method The intermolecular potentials used in all the simulations to be described were the four-site TIP4P model of wateF9 and the OPLS model of ethane.50 The TTp4P potential is successful in predicting the thermodynamic, structural, and dynamic properties of bulk water at both low and high temperatures and p r e s s ~ r e s . ~ ~The . ~ ' -TIP4P ~ ~ and OPLS potentials have been used extensively in simulations of bulk ethanes0 and methanewater giving a good description of the properties of these systems.
Molecular dynamics simulations were conducted in constant ( N , V, E> ensembles using the program Moldy.56 This program
TABLE 1: Summarv of Simulation9 _
_
_
~
T/K
Tdps
TR/ps
ts/fs
PMPa
U/(H mol-')
279.0iz0.5 301.6f0.5 320.2f0.5 351.9iz0.5
90 60 40 30
100 100 100 100
0.50 0.50 0.50 0.25
1 iz5 34iz5 59i5
-9674iz3 -9348*3 -9076f3 -8706iz 3
11Oiz5
T = temperature, TE = equilibration time, TR = run time, ts = timestep, P = pressure, U = average potential energy.
makes use of a modification of the Beeman integration algorithms7 due to Refson and P a ~ l e y . Cubic ~ ~ periodic boundary conditions were applied throughout. In all the simulations a time step of 0.5 fs was used, excect for the one carried out at 352 K in which a time step of 0.25 fs was used. A real space cutoff of 12 A was applied to short-range interactions, with a long-range correction for particles at larger separations. The Ewald sum method was applied to evaluate long-range electrostatic interaction^.^^ The useful data generated during the simulations were stored every 50 time steps. The system consisted of 1 ethane molecule in a bath of 256 water molecules. The initial configuration generated by the program placed the molecules regularly at evenly spaced intervals and assigned them random orientations. An average density of 1.0 g/cm3 was used for all the simulations. Constant volume ensembles were used throughout this work for compatibility with the hydrophobic aggregation studies.32 The required temperatures were obtained by scaling the rotational and translational components of the kinetic energy individually for the 2 molecular species every 20 time steps over a period of 10 ps. We allowed further periods of between 30 and 90 ps for equilibration, with the shorter equilibration periods for the systems at higher temperatures. After these equilibration periods, an additional period of 100 ps in all simulations was allowed for the collection of data. Since a rather long cutoff and short time steps were used, no rescaling of the velocities was needed during the data collection period.
3. Results and Discussion Thermodynamic data for the four simulations can be seen in Table 1. These results were, in general, consistent with separate runs performed on systems of four ethane molecules at approximately the same temperatures reported here.32 Selfdiffusion coefficients for water and ethane were calculated from the slopes of plots of mean-square displacement against time, following the Einstein relation
D = ( ( r - r,J2)/61 where r - ro is the displacement in time t and the average is over all particles of the same species and all choices of time origin over periods of 90 ps, considering that the time ranges of the plots were of 10 ps. The results of are in Table 2. Since only one ethane molecule was present, the statistics in the calculation of many properties were poorer than in the case of water.
TABLE 2: Diffusion Coefficients" T/K
DE$iO-9m2 S-I
279 302 320 352
1.5 rt 1.0 3.3 f 1.0 1.9 f 1.0 4.4 iz 1.0
m2 s-I 2.5 & 0.2 4.4 f 0.2 4.8 iz 0.2 7.4 rt 0.2
DB/iO-9m2 s-I 2.5 f 0.2 4.4 rt 0.2 4.8 iz 0.2 7.5 rt 0.2
m2 s-] 2.1 iz 0.5 4.2 iz 0.5 4.7 & 0.5 7.0 iz 0.5
~ ~ d i m2 o s-I- ~ 2.3 (278 K) 3.2 (301 K) 5.5 (321 K) 6.7 (347 K)
D ~ J ~ m2 O -s-]~ 1.313 (278 K)60 2.3 (298 K)60.61 3.575 (318@ 'K ) 6.0 (348 K)61
T = temperature, DE^ = ethane diffusion coefficient, DW = water diffusion coefficient, DB= bulk water diffusion coefficient, DS = hydrationshell water diffusion coefficient, D w ~= pure water diffusion coefficient, and Dwe = experimental water diffusion coefficient. The errors in the diffusion coefficients in our hydrophobic aggregation study (ref 32) were incorrectly reported and are actually of the same magnitude as the ones reported here.
Hydrophobic Hydration of Ethane For the purpose of this work, we had to make the distiction between hydration-shell and bulk water molecules. Hydrutionshell molecules are defined as those water molecules whose oxygens were located up to a maximum of 5.5 A from either of the methyl groups of the ethane molecule, according to the position of the first minimum in the methyl-oxygen radial distribution function. Bulk molecules are all the other water molecules in the system. It can be seen that at the four temperatures studied, the average diffusion coefficients of the molecules in the hydration shell of the ethane molecule are consistently smaller than those of the molecules in the bulk. This agrees with previous findings. 16.36-40 Furthermore, the higher diffusion coefficients of water with respect to the solute suggest that the hydrationshell molecules are not tightly bound to the ethane. The water molecules become more mobile with increasing temperature; however, it is interesting that at 320 K the diffusion coefficient of ethane actually dropped, and that of water had just a small increase, before rising again at higher temperature. This kind of behavior for the diffusion coefficient of the hydrophobic solute has been found consistently in recent simulation^,^^.^^ where the diffusion of the solutes decreased at those temperatures where a higher hydrophobic aggregation was observed. It should be pointed out that our values for the diffusion coefficient of water at all temperatures are higher than the experimental values, also shown in Table 2. This agrees to some extent with other simulations of pure TIP4P water, where the diffusion coefficients have been found to be higher than the experimental values,60,61and those calculated with the SPCE model, which gives better agreement.53 Our values are higher than those calculated in pure water by Watanabe and Klein53and Brodholt and Wood52 but lower than those calculated by Ferrario and Tani.55 In Table 2 we also report for comparison diffusion coefficients for pure water, obtained from simulations6*under similar conditions over periods of 50 ps after appropriate equilibration periods. It can be seen that the diffusion coefficients for both hydration-shell and bulk water in the mixtures are actually somewhat larger than those in pure water in the temperature range studied. Further work is being done to investigate the nature of the behavior of these diffusion coefficients. Radial distribution functions (RDF) were calculated from the data collected over 100 ps after full equilibration and are shown in Figure 1. It should be pointed out that these functions are consistent with those found in an MC study by Jaros~kiewicz,~~ being the only simulation so far reported of ethane in water. In Figure la, the oxygen-oxygen (0-0) RDF shows that the water in the solutions becomes less structured as the temperature is increased (as seen by the broadening of the peaks), which is what is expected in a liquid when its thermal energy is increased, along with an increase in its mobility, as Table 1 shows. There is no significant difference in the average structure and dynamics of the water in these solutions with respect to those measured in pure bulk The 0-0 RDF agrees well with the experimental structure of pure water,@although recent measurements of the structure of water have revealed the existence of a less high first peak.65 Parts b and c of Figure 1 show the average ethane(methy1)-water structure; the methyl-oxygen (Me-0) RDF agrees quite well with the experimental Me-0 RDF determined in an aqueous solution of methan01.4~ The positions of the peaks in the Me-0 and methyl-hydrogen (Me-H) RDF reflect the existence of a hydration shell in which the dipole vectors are less radial than in the bulk. This structure has been found consistently in the past,20-43-62 and it has been suggested that it could allow for hydrogen bonds to be
J. Phys. Chem., Vol. 99, No. 40,1995 14635
I , 0
1
2
3
4
5
5
7
8
a
I 10
r/A
1.4
,
1
r/A
Figure 1. (a) Oxygen-oxygen radial distribution functions showing how water becomes less structured as temperature is increased. (b) Methyl-oxygen and (c) methyl-hydrogen radial distribution functions showing a nearly temperature-independent hydration structure.
maintained more easily in the hydration shell, similar to what occurs in clathrate hydrates. Since these hydration structures seem to be insensitive to temperature in the range studied, as opposed to the broadening in the 0-0 RDF as temperature is increased, an increasingly high entropic price is being paid to maintain the hydration shell around the solute. It should be expected, however, that at higher temperatures, thermal disruption would begin to dominate and the solvation structure would become increasingly less structured. The simulated ethane-
Mancera and Buckingham
14636 J. Phys. Chem., Vol. 99, No. 40, 1995
TABLE 3: Average Water (H-Bond) Coordination Numbers temperature/K hydration shell bulk 279 3.69 f 0.15 3.72 f 0.05
.. .
-
02
0.1s
-
0.1
-
OIlS
-
302 520 352
O L 1.4
Bond k W h 0.04
I
t
B d saW,dcg
Figure 2. Bulk water relative distributions of hydrogen bond lengths (a) and angles (b), showing the broadening in the distributions as temperature is increased. No significant differences were observed in the case of the hydration-shell water.
water structure is consistent with the hypothesis that hydrophobic hydration is entropically costly near room due to a structuring of the water molecules in the vicinity of the nonpolar solute. However, it is still unclear what the molecular mechanism driving this structuring process is and we present further evidence on this in our water analysis. 3.1. Hydrogen Bond Analysis. We proceeded to analyze the structural properties of the hydrogen-bonded network in the hydration-shell and bulk regions. For this purpose, water molecules were considered to be near-neighbors if their oxygens were 53.5 8, apart. The hydrogen bond between near-neighbor molecules was then selected as that having the minimum H**.O distance among the four possible combinations of H and 0 distances. A hydrogen bond angle was defined as the angle between the 0-H bond vector of one water molecule and the O..*H hydrogen bond vector between this molecule's H and its near-neighbor's 0. One- (not shown) and two-dimensional plots (see below) for the distribution of hydrogen bond lengths and/ or hydrogen bond angles in the system were made to establish an appropriate structural definition of a hydrogen bond. We defined a hydrogen bond as having a maximum length (H..*O) of 2.5 8, and a hydrogen bond angle between 130" and 180". With this definition, a hydrogen bond is considered to be "stronger" the shorter its length and the more linear its angle. Figure 2a shows the distribution of the hydrogen bond lengths (up to 2.5 A, as previously defined) at all temperatures for the bulk water regions. The hydration-shell water analysis showed
3.61 f 0.17 3.51 f 0.18 3.38 & 0.19
3.64 & 0.05 3.57 i0.06 3.47 f 0.06
no significant differences at all temperatures (apart from statistical noise due to a poorer sampling). As temperature is increased, there is a broadening of the distributions, as expected from the additional energy in the system, although the most populated length (1.8-1.9 A) remains the same. Figure 2b shows the distribution of hydrogen bond angles (between 130" and 180°, as previously defined) again at all temperatures for the bulk regions. Again, at all temperatures there were no significant differences with the hydration shell. We again find that as temperature is increased, there is a broadening of the distributions, and in this case the most populated angle moves from about 165" at 279 K to about 162" at 352 K. The main observation is that there are no significant differences in the hydrogen bond properties between the hydration shell and bulk regions. Further evidence comes from two-dimensional plots of the distribution of hydrogen bond lengths and the associated angles, without considering the previous definition of a hydrogen bond so as to screen the whole confrgurational space available. Figure 3 shows the differences in these distributions at 279 and 352 K for the hydration shell and bulk regions. As the temperature goes from 279 to 352 K, for both regions we observe a decrease in the low-energy peak (precisely in the region of short and linear hydrogen bonds: distances up to 2.5 8, and angles between 130" and 180", as in our definition of a hydrogen bond) and an increase in the high-energy peak (for long and nonlinear hydrogen bonds). Comparing the plots at the same temperature between the hydration-shell and bulk regions, we see slightly higher and sharper low-energy peaks for the hydration-shell regions. The evidence from Figure 3 is that molecules in the hydration shell are able to make hydrogen bonds not only with the same structural properties as those in the bulk but also slightly "stronger". We then calculated the average hydrogen-bond coordination numbers (defined here as the maximum number of hydrogen bonds that a water molecule can engage in) for the hydration-shell and bulk water molecules. This was done by summing the total number of hydrogen bonds in which hydration-shell and bulk molecules were engaged at every 25 fs and dividing by the number of water molecules in the hydration shell and in the bulk at that moment. The coordination numbers presented in Table 3 are the time-averages over the whole simulation. The average coordination numbers of the hydration-shell molecules are consistently slightly smaller than that of the bulk molecules, the increasing temperature just having the effect of reducing these averages as the average number of nearest-neighbors decreases (see also the 0-0 RDF in Figure l a to corroborate this, since the area under the first peak decreases slightly with increasing temperature). The distribution of the number of hydrogen bonds (coordination numbers) for both the hydration-shell and bulk molecules can be seen in Figure 4. As the temperature is increased, there is a tendency for the molecules to have smaller (1 -3) coordination numbers and greater departure from tetrahedral coordination, which is what is expected as thermal motion increases. At 279 K the fraction of molecules with tetrahedral coordination is larger for the hydration shell at the expense of higher coordination numbers (5); however, as temperature is increased,
J. Phys. Chem., Vol. 99, No. 40, 1995 14637
Hydrophobic Hydration of Ethane a
b 0.006-
0 008
0.005
0,007
1 ::E
.-
0.006
a
.e
0.001
EP
352K
6
0.004
$
0.003
8
0.002
g
0.001
2
0
‘Fi
o
15
15
C
d
Figure 3. Relative distributions of hydrogen bond lengths and angles for (a) hydration shell at 279 K, (b) hydration shell at 352 K, (c) bulk at 279 K, and (d) bulk at 352 K, showing how hydrogen bcinds become longer and less linear as temperature is increased. Hydrogen bonds in the
hydration shell have slightly higher and sharper low-energy (short and linear) distributions.
-0 0
1
2
3
nc-
4
5
b
0
1
2
3
4
5
6
4
5
6
- w w
Nu-,
n7.
0
1
2
3
ceC.-?c”*um*r
4
5
6
0
1
2 5 c m -
Nu-
Figure 4. Relative distributions of water hydrogen-bond coordination numbers for hydration shell (m) and bulk (0)molecules at (a) 279, (b) 302, (c) 320, and (d) 352 K, showing an increased distribution of lower
coordination numbers with increasing temperature, the: effect being greater in the hydration shell. hydration-shell molecules have an increasingly higher probability than bulk molecules of having a small coordination number (1-3) at the expense of the latter now having a higher probability of tetrahedral coordination. So it seems that the increase in thermal energy is more capable of inhibiting tetrahedral coordination in the hydration-shell region than in the bulk, whereas we see an increase of the low coordination numbers (1-3) in the hydration shell with respect to the bulk. For higher coordination numbers, it is obvious that hydrationshell molecules should have smaller values due to geometric factors, since the number of available near-neighbors is also smaller than in the bulk. As pointed out by Flossky and K a r p l u ~in~ ~a molecular description of the hydration of a dipeptide, since the average coordination numbers in the hydration shell are roughly the same as in the bulk, a more developed hydrogen-bonding network should exist among water
molecules in the hydration shell to allow for bonding with a lower fraction of near-neighbor molecules. These results show that hydration-shell molecules are capable of engaging in normal hydrogen bonding, although the presence of the nonpolar solute does have an effect on the way temperature affects the changes in coordination of these hydration-shell molecules. This reinforces the picture that when a nonpolar molecule is dissolved in water, water molecules are still able to maintain a strong (perhaps even stronger) hydrogenbonded network, thus contributing to the net negative enthalpy change observed. Since the stability of this bonding network exists at the expense of reduced entropy, it is reasonable to believe that hydrogen bonds in the hydration-shell region should be more effectively disrupted than in the bulk region by an increase in temperature. In other words, the entropic contribution to the free energy of solution should be increasingly important with temperature increase, as observed experimentally by the decrease in solubility with increasing temperature near room temperature. However, the evidence presented so far still does not give a clear picture of the origin of the net negative entropy change of this process. We therefore extended our analysis to the structural correlation between water molecules in the hydration-shell and bulk regions. 3.2. Water Analysis. We initially, believed that the origin of the negative entropy change should reside in the orientational restraints imposed on water molecules near a nonpolar solute. Since water molecules, interacting mainly by means of hydrogen bonds, can establish an effective hydrogen-bonded network around the nonpolar solute, we decided to look at the orientational correlations between water molecules in the hydrationshell and bulk regions by analyzing the orientations of the water molecules with respect to the solute and other water molecules, respectively. The vector joining the solute (either of the methyl groups of ethane for hydration-shell molecules or a water
14638 J. Phys. Chem., Vol. 99, No. 40, 1995
1 I
Mancera and Buckingham
0.01
0.008
0.008
0.w
1
0
l 20
40
Bo
Bo
1w
120
140
180
180
1W
120
140
180
180
dw
0.014
0.012
0.01
0.m
0.m
I 0
20
40
Bo
Bo
Wn 0.02
270K 3026
0.018
.._. .....
.-
0.016 0.014
fd
0.012
0’01 0.m 0.008
0.004
0
20
40
Bo
80
1w
I 120
140
180
180
erwg
Figure 5. (a,b) Relative normalized distributions of the orientation angles relative to the solute for hydration-shell and bulk water molecules: (a) a and (b) p (only hydration-shell water distributions are identified), showing the “straddling” conformation of water molecules in the hydration shell. ( c ) Relative normalized distributions of the orientation angles (aand p) for bulk water molecules relative to other water molecules beyond the second hydration shell, showing the local tetrahedral arrangement.
molecule’s oxygen for bulk molecules) to the oxygen subtends an angle a with the water dipole moment vector and an angle ,8 with a vector perpendicular to the, H-0-H molecular plane. Parts a and b of Figure 5 show the normalized distributions of the a and /? angles for the hydration-shell and bulk molecules for all temperatures studied. It is clear that the differences are
significant, although nothing clear can be said about the effect of temperature on the distributions of a or ,!? on the hydrationshell molecules because there is too much noise; however, in the case of the bulk molecules, there is a uniform distribution of orientations for both a and ,!?at all temperatures. Vaisman et al.43did the same analysis for the a angle in the hydration of methane and found similar results. The peaks near 70” for the a angle and the preferred nearly parallelhtiparallel orientations of the ,!? angle correspond to a nearly tangential orientation o f water molecules in the first hydration shell, as already mentioned in our analysis of the ethane-water RDF. This type of a,tructural arrangement of water molecules in the vicinity of a n’onpolarsolute allows them to straddle the surface of the solute and maintain a nearly tetrahedral hydrogen-bond coordination. Figure 5c !shows also the normalized distributions of the a and angles for water molecules in the vicinity of another water molecule, instead of the solute. For this analysis, water molecules more than 8.5 A away from any of the methyl sites of ethane were selected as origins for the angle measurements, and the distribution of the orientations of near-neighbor water molecules (0-0 distances up to 3.5 A) was measured. This allowed for water molecules beyond the second hydration shell to be analyzed with near-neighbor water molecules always outside the first hydration shell. The distribution of a angles reveals the existence of two characteristic peaks: a narrow peak around 127” corresponding to the average orientation of dipole vectors of water molecules which act as hydrogen “donors” and a broader peak around 45” corresponding to the average orientation of dipole vectors of those water molecules which act as hydrogen “acceptors”. Since hydrogen bonds tend to be nearly linear, water molecules acting as hydrogen “donors” have less orientational freedom with respect to hydrogen “acceptor” water molecules. The distribution of the ,!? angles shows the existence of two preferred orientations: a parallel (or antiparallel, depending on the direction of the vector perpendicular to the H-0-€1 molecular plane) orientation (broad peak around 0”/180”) for water molecules acting as hydrogen “acceptors” and a nearly perpendicular orientation (narrow peak around 92”) for water molecules acting as hydrogen “donors”. These angular distributioris result from the local tetrahedral structure of liquid water. It is clear from the previous analyses that hydration-shell molecules might have only a restricted orientational space available for them. To prove this point, we made twodimensional plots of the distributions of a and /? for the hydration-shell water molecules and for the near-neighbor molecules of bulk water (as defined in the last paragraph). We present the results at 302 K in Figure 6; the plots at the other temperatures showed no significant differences. These plots give a clear qualitative picture of the orientational restrictions which are imposed on the water molecules in the vicinity of a nonpolar solute. In the case of hydration-shell water molecules, the orientational space available is consistent with the picture of an average clathrate-like structure, as explained before, that is often assumed in the hydration of nonpolar s o l ~ t e s . ~ - ~ ~ ’ ~ ~ ’ ~ Bulk water molecules show the typical nearly tetrahedral structure reminiscent of ice. The major contribution to the welldefined contours in the region defined by 100” < a < 150” and 70” ,8 > 120” comes from the hydrogen-bond “acceptor” molecules The differences in the “allowed” orientations the water molecules can adopt in the hydration shell and in the bulk
J. Phys. Chem., Vol. 99, No. 40, 1995 14639
Hydrophobic Hydration of Ethane
1
1 ,
0.7
50
100
150
02
ddcp
01
,
,.
0.000506 ----. 0.000404 ..... 0,000303 0.000202 -.-.0.000101 -.-.-
-
-
.............. ................... ...... _..... -. ......................... :;.-.-.......................... , , ._, I
OO
200
4w
Bw
low
800
1m
1100
-
-
1BW
1nM
2ow
Timclfr
P ....
..........
0.3
50
100
150
ddcg
Figure 6. Relative distributions at 302 K of the a and p angles for (a) hydration-shell molecules with respect to the solute and (b) bulk water molecules with respect to other water molecules, showing the differences in the orientational space available for water molecules in the hydration shell and in the bulk.
are striking, and it is possible that the orientations of the hydration-shell molecules are more restrained than those in the bulk. 3.3. Hydrogen Bond Dynamical Analysis. As entropy decrease in the hydration shell has been suggested as the cause of lower mobilities of these molecules, we investigated the behavior of the hydrogen-bonded network in terms of rates of bond breaking for both hydration-shell and bulk molecules. For this purpose, we followed the method outlined by Haughney et aI.,66 where the fraction of broken hydrogen bonds in time is calculated from a histogram containing the number of bonds that break after a given time. These bond-breaking rates depend on the time intervals used to check whether a hydrogen bond has broken or not, although the trends should be the same. We allowed intervals of 25 fs before inspecting the list of all existing hydrogen bonds. Other time intervals showed the same behavior reported here. Parts a and b of Figure 7 show the time decay functions of the fraction of unbroken hydrogen bonds at all temperatures studied for hydration-shell and bulk regions, respectively. We note that for both regions, as temperature is increased, hydrogen bonds break up more rapidly, as expected from the increased thermal energy, and they break up more rapidly in the hydration-shell region than in the bulk region, the difference becoming less notable as temperature is increased. We suppose that the translation and rotation of any water molecule occur by a "concerted" mechanism with other molecules so that as hydrogen bonds between pairs of water molecules are gradually being strained and eventually broken, other hydrogen bonds are simultaneously and gradually being formed with other near-neighbor molecules. The reorientation
0.2
%._
0.1
0 '
......
........ ......... ................. . . ....
200
100
800
EW
...
............ ....... ......
1WO
__
12w
1400
1BW
lB00
2wO
T m i
Figure 7. Hydrogen-bond rupture decay functions for (a) hydration shell and (b) bulk, showing an increased rate of hydrogen-bond breaking in the hydration shell with respect to that in the bulk.
of the molecules involved should then involve a significant degree of correlation. We believe that the differences in mobility between hydration-shell and bulk molecules (as already described by their diffusion coefficients) should be attributed to differences in their capability to carry out such a concerted rearrangement of hydrogen bonding, which in the case of the hydration-shell molecules is limited due to a reduced number of low-energy paths. The solute's main influence resides not only in imposing geometrical constraints on the water molecules, but also in its own diffusion, as it rattles in its solvation cage before breaking it at some point and moving to a "new" solvation cage. On the other hand, the orientational rearrangements of hydration-shell molecules as they diffuse should avoid placing a hydrogen-bonding group directed toward the nonpolar solute. Consequently, the combined diffusion of the solute and the water molecules under their conformational restrictions might induce a faster molecular reorientation mechanism (thus requiring hydrogen bonds to break faster) involving a larger degree of cooperativity. 4. Conclusions
We have performed a series of ( N , V, E> molecular dynamics simulations of ethane in water at different temperatures near room temperature to try to obtain a clearer picture of the molecular mechanisms responsible for the observed thermodynamic properties of the dissolution of nonpolar solutes in water and the effect of temperature on them. Our structural analysis showed the existence of an almost temperature-independent ethane-water structure as revealed by
14640 J. Phys. Chem., Vol. 99, No. 40, I995
radial distribution functions, as opposed to the water-water structure which showed a reduction in structure with increasing temperature. We detected the presence of a different orientational space for water molecules in the hydration shell when compared to those in the bulk. An analysis of the hydrogenbonded network showed a slight structural enhancement of the tetrahedral network, as revealed by the existence of nearly identical hydrogen-bond structural profiles between hydrationshell and bulk water molecules, with a slight shift toward shorter and more linear hydrogen bonds in the case of hydration-shell molecules. This suggests that the energetic environment in the hydration shell is not much different from that in the bulk and reinforces the idea of a reduction in the number of low-energy paths for molecular orientational rearrangements in the vicinity of nonpolar solutes. The dynamic analysis showed that hydration-shell molecules have slightly lower diffusion coefficients than bulk molecules. Furthermore, calculations of hydrogen-bond rupture profiles surprisingly revealed that hydrogen bonds in the hydration-shell region break up faster than those in the bulk region. A dynamic picture combining the solute’s and hydration-shell water’s “concerted” diffusion, along with the water’s orientational restrictions, helps explain this phenomenon in terms of faster molecular reorientations in the diffusion process, thus warranting further studies on the mechanisms of diffusion in aqueous solutions of nonpolar solutes.
Acknowledgment. We are indebted to Dr. N. T. Skipper for his advice and comments, and R.L.M. gratefully acknowledges the National University of Mexico, the British Council, the CVCP of the Universities of the U.K., and the Cambridge Overseas Trust for their support through various grants. We also thank DGSCA (Mexico) and ULCC (U.K.) for generous awards of computer time on the Cray Y-MP4 and the Convex C3800, respectively. References and Notes (1) Tanford, C. The Hydrophobic Effect: Formation of Micelles and Biological Membranes; Wiley: New York, 1973. (2) Franks, F. In Water, a Comprehensiue Treatise; Franks, F., Ed.; Plenum Press: New York, 1975; Vol. 4, Chapter 1. (3) Ben-Naim, A. Hydrophobic Interactions; Plenum Press: New York, 1980. (4) Muller, N. Acc. Chem. Res. 1990, 23, 23. (5) Costas, M.; Kronberg, B.; Silveston, R. J. Chem. SOC.,Faraday Trans. 1994, 90, 1513. (6) Butler, J. A. V.; Reid, W. S. J. Chem. SOC.1936, 1171. (7) Butler, J. A. V. Trans. Faraday SOC. 1937, 229. (8) Miller, K. W.; Hildebrand, J. H. J. Am. Chem. SOC. 1968, 90,3001. (9) Shinoda, K.; Fujihira, M. Bull. Chem. SOC.Jpn. 1968, 41, 2612. (10) Shinoda, K. J . Phys. Chem. 1977, 81, 1300. (11) Dec, S. F.; Gill, S. J. J . Solution Chem. 1985, 14, 417, 827. (12) Privalov, P. L.; Gill, S. J. Adu. Protein Chem. 1988, 39, 191. (13) Privalov, P. L.; Gill, S. J. Pure Appl. Chem. 1989, 6 I , 1097. (14) Gill, S . J.; Dec, S. F.; Olofsson, G.; Wadso, I. J. Phys. Chem. 1985, 89, 3758. (15) Lazaridis, T.;Paulaitis, M. E.J . Phys. Chem. 1992, 96, 3847. (16) Geiger, A.; Rahman, A.; Stillinger, F. H. J. Chem. Phys. 1979, 70, 263.
Mancera and Buckingham (17) Nbmethy, G.; Scheraga, H. A. J. Chem. Phys. 1962, 36, 3401. (18) Tse, J. S.; Klein, M. L.; McDonald, I. R. J . Phys. Chem. 1983, 87, 4198. (19) Rodger, P. M. J. Phys. Chem. 1990, 94, 6080. (201 Paneali. C.: Rao. M.: Beme. B. J. J. Chem. Phvs. 1979. 71. 2982. (21) RavTshanker, G.; Mezei, M.; Beveridge, D. L. Furaday Symp. Chem. SOC.1982, 17, 79. (22) Rapaport, D. C.; Scheraga, H. A. J . Phys. Chem. 1982, 86, 873. (23) Watanabe, K.; Andersen, H. C. J. Phys. Chem. 1986, 90, 795. (24) Jorgensen, W. L.; Buckner, J. K.; Boudon, S.; Tirado-Rives, J. J. Chem. Phys. 1988, 89, 3742. (25) Wallqvist, A. J . Phys. Chem. 1991, 95, 8921. (26) Laaksonen, A.; Stilbs, P. Mol. Phys. 1991, 74, 747. (27) Smith, D. E.; Zhang, L.; Haymet, A. D. J. J. Am. Chem. SOC.1992, 114, 5875. (28) Smith, D. E.; Haymet, A. D. J. J . Chem. Phys. 1993, 98, 6445. (29) Dang, L. X. J . Chem. Phys. 1994, 100, 9032. (30) van Belle, D.; Wodak, S. J. J. Am. Chem. SOC. 1993, 115, 647. (31) Skipper, N. T. Chem. Phys. Lett. 1993, 207, 424. (32) Mancera, R. L.; Buckingham, A. D. Chem. Phys. Lett. 1995,234, 296. (33) Owicki, J. C.; Scheraga, H. A. J. Am. Chem.’Soc. 1977, 99,7413. (34) Swaminathan, S.; Hamson, S. W.; Beveridge, D. L. J. Am. Chem. SOC. 1978, 100, 5705. (35) Pangali, C.; Rao, M.; Beme, B. J. J . Chem. Phys. 1979, 71,2975. (36) Rossky, P. J.; Karplus, M. J. Am. Chem. SOC.1979, 101, 1913. (37) Alagona, G.; Tani, A. J . Chem. Phys. 1980, 72, 580. (38) Bolis, G.; Clementi, E. Chem. Phys. Lett. 1981, 82, 147. (39) Jorgensen, W. L. J . Chem. Phys. 1982, 77, 5757. (40) Rossky, P.J.; Zichi, D. A. Faraday Symp. Chem. SOC. 1982, 17, 69. (41) Guillot, B.; Guissani, Y.; Bratos, S. J. Chem. Phys. 1991,95,3643. (42) Guillot, B.; Guissani, Y. J. Chem. Phys. 1993, 99, 8075. (43) Vaisman, I. I.; Brown, F. K.; Tropsha, A. J . Phys. Chem. 1994, 98, 5559. (44) Broadbent, R. D.; Neilson, G. W. J. Chem. Phys. 1994,100,7543. (45) Tumer, J.; Soper, A. K. J. Chem. Phys. 1994, 101, 6116. (46) Finney, J. L.; Soper, A. K.; Tumer, J. Z . Pure Appl. Chem. 1993, 65, 2521. (47) Soper, A. K.; Finney, J. L. Phys. Rev. Lett. 1993, 71, 4346. (48) Bradl, S.; Lang, E. W.; Tumer, J. Z.; Soper, A. K. J. Phys. Chem. 1994, 98, 8161. (49) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (50) Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. J . Am. Chem. SOC.1984, 106, 6638. (51) Jorgensen, W. L.; Madura, J. D. Mol. Phys. 1985, 56, 1381. (52) Brodholt, J.; Wood, B. Geochim. Cosmochim. Acta 1990,54, 261 1. (53) Watanabe, K.; Klein, M. L. Chem. Phys. 1989, 131, 157. (54) Jorgensen, W. L.; Gao, J.; Ravimohan, C. J. Phys. Chem. 1985, 89, 3470. (55) Ferrario, M.; Tani, A. Chem. Phys. Lett. 1985, 121, 182. (56) Refson, K. MOLecular Dynamics Simulation Code, Moldy; Department of Earth Sciences, Oxford University: Parks Road, Oxford OX1 3PR, 1988, 1992, 1993. (57) Beeman, D. J . Comput. Phys. 1976, 20, 130. (58) Refson, K.; Pawley, G. S. Mol. Phys. 1987, 61, 669. (59) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1987. (60) Mills, R. J. Phys. Chem. 1973, 77, 685. (61) Kryniclu, K.; Green, C. D.; Sawyer, D. W. Faraday Discuss. Chem. SOC.1978, 66, 199. (62) Mancera, R. L.; Buckingham, A. D. Unpublished results. (63) Jaroszkiewicz, G. A. Chem. Phys. Lett. 1984, JOB, 443. (64) Soper, A. K.; Phillips, M. G. Chem. Phys. 1986, 107, 47. (65) Soper, A. K. J . Chem. Phys. 1994, 101, 6888. (66) Haughney, M.; Ferrario, M.; McDonald, I. R. J. Phys. Chem. 1987, 91, 4934. JP950003B
