Temperature and Pressure Dependence of ... - ACS Publications

Article pubs.acs.org/JPCB

Temperature and Pressure Dependence of Methane Correlations and Osmotic Second Virial Coefficients in Water Henry S. Ashbaugh,*,† Katie Weiss,‡ Steven M. Williams,† Bin Meng,† and Lalitanand N. Surampudi† †

Department of Chemical and Biomolecular Engineering, Tulane University, New Orleans, Louisiana 70118, United States Alfred University, Alfred, New York 14802, United States

‡

S Supporting Information *

ABSTRACT: We report methane’s osmotic virial coefficient over the temperatures 275 to 370 K and pressures from 1 bar up to 5000 bar evaluated using molecular simulations of a united-atom description of methane in TIP4P/2005 water. In the first half of this work, we describe an approach for calculating the water-mediated contribution to the methane−methane potential-of-mean force over all separations down to complete overlap. The enthalpic, entropic, heat capacity, volumetric, compressibility, and thermal expansivity contributions to the water-mediated interaction free energy are subsequently extracted from these simulations by fitting to a thermodynamic expansion over all the simulated state points. In the second half of this work, methane’s correlation functions are used to evaluate its osmotic second virial coefficient in the temperature− pressure plane. The virial coefficients evaluated from the McMillan−Mayer correlation function integral are shown to be in excellent agreement with those determined from the concentration dependence of methane’s excess chemical potential, providing an independent thermodynamic consistency check on the accuracy of the procedures used here. At atmospheric pressure the osmotic virial coefficient decreases with increasing temperature, indicative of increasing hydrophobic interactions. At low temperature, the virial coefficient decreases with increasing pressure while at high temperature the virial coefficient increases with increasing pressure, reflecting the underlying hyperbolic dependence of the virial coefficient on temperature and pressure. The transition between a decreasing to increasing pressure response of the osmotic virial coefficient is shown to follow the response of the methane−methane contact peak to changes in pressure as a function of temperature, though a universal correlation is not observed.

■

folding,25−33 resulting in cold and pressure induced denaturation.18,34 While the methane−methane PMF in principle can be determined from scattering experiments, the poor solubility of methane in aqueous solution makes measurement of the scattering intensity untenable. Resultantly, direct experimental verification of the forces between methanes mediated by water predicted from molecular simulation is lacking. Alternatively, the sign and magnitude of hydrophobic interactions can be inferred through the osmotic second virial coefficient (B2). Specifically, B2 depends on molecular-level solute pair correlations in aqueous solution through the McMillan−Mayer35 theory expression

INTRODUCTION Nonpolar gases have long stood as model solutes for studying hydrophobic effects since they lack the strong polar and charge interactions with water that complicate interpretation of the dissolution of more complex species.1−8 These gases have been employed to scrutinize the characteristic thermodynamics of hydrophobic hydration on the molecular-scale, from the dominant unfavorable entropies that resist dissolution in water at room temperature, to large positive hydration heat capacity increments, to the observation of minima in their aqueous solubility with increasing temperature.4,9−12 Results for nonpolar gas hydration in turn have informed interpretation of the thermodynamics of self-assembly in aqueous solution, like micellization and protein folding,7,13,14 and the curious temperature dependencies of macromolecular systems and biomolecular assemblies, including the conformational collapse of thermoresponsive polymers and cold denaturation of proteins.15−20 Methane−methane interactions in aqueous solution have similarly served as a starting point for interpreting molecularlevel hydrophobic interactions.21−24 Molecular simulations of the water-mediated interaction between methane pairs in solution, embodied in their potential-of-mean force (PMF), have been used to analyze the temperature and pressure response of hydrophobic contact formation in protein © 2015 American Chemical Society

B2 =

1 2

∫ [1 − gss(r)] dr

(1)

where gss(r) is the solute−solute (s) radial distribution function (RDF) in aqueous solution. When the RDF exhibits prominent solute−solute contacts, B2 is expected to be negative, indicative of attractive interactions in solution. Conversely, unfavorable contacts are expected to make B2 more (or even net) positive, indicative of repulsive interactions in solution. While it is Received: March 2, 2015 Revised: April 24, 2015 Published: May 1, 2015 6280

DOI: 10.1021/acs.jpcb.5b02056 J. Phys. Chem. B 2015, 119, 6280−6294

Article

The Journal of Physical Chemistry B

4.0 simulation package.45 The temperature and pressure were maintained using the Nosé−Hoover thermostat46,47 and the Parrinello−Rahman barostat,48 respectively. Water and methane were modeled using the TIP4P/200549 and OPLS/UA50 potentials, respectively. Water’s internal holonomic constraints on the oxygen, hydrogen, and auxiliary site positions were maintained using the LINCS algorithm. 51 Nonbonded Lennard-Jones (LJ) interactions were truncated beyond 10 Å, while particle mesh Ewald summation52 was used to evaluate long-range electrostatic interactions. Cross interactions between LJ sites were evaluated using Lorentz−Berthelot mixing rules.53 A 2 fs time step was used to integrate the equations of motion in all the simulations discussed below. To evaluate methane−methane pair correlations in water, mixtures of 40 methanes in a box of 4000 water molecules were simulated. Following 5 ns of equilibration, simulations were conducted for 100 ns to evaluate equilibrium averages. A first set of calculations was performed from 275 to 370 K in 5 K increments at 1 bar to evaluate the temperature dependence of these correlations. A second set of calculations were performed at pressures of 500, 1000, 1500, 2000, and 2500 bar and temperatures of 280, 300, 320, 340, and 360 K to evaluate the pressure dependence of these correlations. Additional simulations at 3000, 3500, 4000, 4500, and 5000 bar were performed to extend the pressure range considered at 300 K over an even wider range. As a result, methane−methane pair correlations were determined at 50 different state points to characterize their temperature and pressure dependence. To evaluate the concentration dependence of methane’s excess chemical potential in water, mixtures of 0, 10, 20, 30, and 40 methanes in 4000 water molecules were simulated. Following 5 ns of equilibration, simulations were conducted for 100 ns to evaluate equilibrium averages. Simulations were performed at 1 bar and temperatures of 280 to 370 K in 10 K increments to evaluate the temperature dependence on the concentration dependence of the excess chemical potential. Since determination of the excess chemical potential at elevated pressures is statistically more unreliable, however, we confined these simulations to only examination of the temperature dependence at low pressure. The excess chemical potential of methane in water was evaluated using the test particle insertion formula54

problematic to construct osmotic membranes permeable to water and impermeable to methane and other nonpolar gases, exact thermodynamic relationships between B2 and the concentration dependence of the solute hydration free energy have been developed and used to extract virial coefficients for short alkanes and inert gases from experimental solubility data.36−38 Estimates of B2 for individual solutes, however, can widely vary in both sign and magnitude.36,37 This may be attributed in part to experimental uncertainties in the solute’s partial molar volumetric properties and their dependence on pressure.39,40 Evaluation of B2 for hydrophobes in water from RDFs determined from simulation has been problematic as well.22,28,36 This difficulty stems in part from the necessity of introducing upper bounds on the integral in eq 1 resulting in poor convergence. As new evaluation techniques have been introduced and computers have advanced that allow longer simulations of larger systems, confidence in the convergence of the B2 integral has grown. Recent simulation studies of methane41 and argon-sized hard sphere solutes42 in water have recently demonstrated that B2 systematically decreases with increasing temperature, indicative of growing, endothermic hydrophobic attractions. Moreover, a recent simulation study of krypton interactions in water at 300 K has demonstrated excellent agreement between B2’s determined using the simulation RDF and from the dependence of krypton’s excess chemical potential on concentration,43 giving confidence that simulations produce thermodynamically consistent results for B2. Motivated by the mounting reliability of B2’s determined from simulation and the general interest in understanding hydrophobic interactions as a function of both temperature and pressure, we set out here to systematically evaluate virial coefficients over a wide range of thermodynamic state points using molecular simulations. In the first half of this work, we present simulation results for the RDF between methanes in aqueous solution over the range of liquid temperatures and pressures from atmospheric conditions to several thousand atmospheres. A new interpolative methodology is introduced to extract the water-mediated contribution to the PMF between methanes in solution that bridges the range of separations passively observed during simulation (i.e., separations beyond methane’s excluded volume) down to complete overlap between methanes. The enthalpic, entropic, heat capacity, volumetric, compressibility, and thermal expansivity contributions to the water-mediated interaction free energy are subsequently extracted from the simulated set of state points, allowing us to contrast thermodynamic contributions to hydrophobic interactions at solute contact against complete overlap. In the second half of this work, we evaluate methane B2’s in aqueous solution using an alternative formulation of the McMillian−Mayer integral presented above useful for obtaining convergence in a closed ensemble simulation. The reliability of the results are determined by comparison with virial coefficients evaluated from the concentration dependence of methane’s excess chemical potential in water. The B2’s determined from simulation are subsequently correlated using a standard thermodynamic expansion over the simulated range of state points to permit their rapid analytical evaluation over the temperature−pressure plane.

μmex = −kT ln

⟨V exp( −Δφ /kT )⟩ ⟨V ⟩

(2)

where kT is the product of Boltzmann’s constant and the absolute temperature, Δφ is the change in the system’s potential energy upon addition of a single methane to the aqueous mixture at a random point in space, V is the instantaneous volume of the simulation, and the angled brackets ⟨···⟩ indicate averages conducted over system mixture configurations generated without the inserted methane included. Test particle calculations were conducted on simulation configurations saved once every picosecond by inserting 5000 test methanes randomly in the simulation box. To assess the impact of the water model used on the partial molar volume and expansivity of methane in solution, a smaller set of simulations were performed of a single OPLS methane in 350 TIP3P and TIP4P waters55 to compare against results for TIP4P/2005. Simulations were conducted from 280 to 370 K in 10 K increments at 1 bar. Following at least 1 ns of equilibration, simulations were conducted for 25 ns to obtain

■

MOLECULAR SIMULATIONS Isothermal−isobaric molecular dynamics44 simulations of water and aqueous mixtures were conducted using the GROMACS 6281


Article


processes has been discussed in the literature.22,27,29−32,56 Here we focus on the impact of temperature and pressure on the solvent-mediated contribution to the methane−methane PMF over a broad range of thermodynamic state points and methane separation. The PMF can determined from the logarithm of the RDF and divided into the direct interactions between methanes (m) in a vacuum and indirect water-mediated interactions in solution as

equilibrium averages. The partial molar volume was determined as the difference in the average volume obtained with and without the added methane, i.e., ⟨VNw+1⟩ − ⟨VNw⟩.

■

RESULTS AND DISCUSSION Water’s Contribution to the Methane−Methane PMF. Methane−methane RDFs harvested from simulations of 40 methanes dissolved in 4000 water molecules as functions of temperature at 1 bar and pressure at 300 K are reported in Figure 1. These pair correlations display many of the previously

w −kT Sngmm(r ) = φmm(r ) + ωmm (r )

(3)

Here r is the separation between methanes, φmm(r) is the direct w methane−methane LJ interaction, and ωmm (r) is water’s indirect contribution to the PMF. The solvent-mediated interaction is frequently referred to as the cavity correlation function in discussions of the RDF, i.e., w ωmm (r ) = −kT Snymm (r ).21 While φmm(r) diverges to infinity due to the excluded volume between two interacting methanes at short separations, no divergence is observed for ωwmm(r), which is finite and well-behaved all the way down to total methane overlap (r → 0). The solvent contribution to the methane−methane PMF at 300 K and 1 bar is displayed in Figure 2. Quantitative

Figure 1. Methane−methane radial distribution function in water evaluated from molecular simulation and fits of the simulation results for the solvent contribution to the PMF (eq 7). (a) Methane-methane radial distribution function at temperatures of 280, 300, 320, 340, and 360 K at a pressure of 1 bar. (b) Methane-methane radial distribution function at pressures of 1, 1000, 2000, 3000, 4000, and 5000 bar at a temperature of 300 K. Symbols are defined in the figure legend. Simulation errors are smaller than the figure symbols.

Figure 2. Water-mediated contribution to the methane−methane PMF at a temperature of 300 K and pressure of 1 bar. Results are shown for the PMF obtained directly from the simulation radial distribution function (eq 3), the free energy at methane overlap, and for the fit of eq 5 bridging the simulation and overlap free energies. Our PMF is compared to that reported by Levy and co-workers in ref 33 obtained independently by free energy perturbation. Symbols are defined in the figure legend. Levy et al.’s results are shifted down from the present results by −2 kJ/mol for clarity.

documented trends for nonpolar gas pair interactions in water. Specific features include a methane−methane contact peak at rcon ≈ 3.9 Å, a solvent-separated peak at rss ≈ 7.0 Å, and a contact barrier at rbar ≈ 5.5 Å. The methane−methane contact peak at 1 bar grows with increasing temperature while the solvent-separated peak mildly decreases (Figure 1a), signifying growing hydrophobic interactions that drive aggregation and phase segregation with rising temperature. The contact and solvent-separated peaks grow in tandem with increasing pressure at 300 K (Figure 1b), following the closer packing of the solvent with increasing pressure. This intuition will be shown below not to hold at elevated temperatures, however, as a result of the expansivity at contact. The equilibrium between the contact and solvent-separated states and its relation to hydrophobic assembly/disassembly

agreement is observed between the solvent contribution obtained here and that reported by Levy and co-workers33 from free energy perturbation calculations performed for methane pairs simulated over a series of fixed separations. Below 3.15 Å, the RDF we observe is zero, making evaluation of the solvent contribution to the PMF over all separations by passive simulation observation untenable due to methane’s excluded volume repulsion. Levy et al.’s perturbation calculations, 33 however, show the solvent contribution smoothly decreases well into the inaccessible excluded volume region down to separations of at least 2.4 Å. In the case of complete methane overlap (r = 0 Å), the solvent contribution to the PMF is expected to drop down to ωwmm(0) = −26 kJ/mol 6282


Article

The Journal of Physical Chemistry B (Figure 2), determined as the free energy difference between two overlapping methanes less twice the hydration free energy of an individual methane (see Appendix for discussion of evaluation of the overlap free energy). We note that this overlap free energy is significantly lower than the free energy of transferring methane from solution to a vacuum (the excess hydration free energy of methane in water at these conditions is ∼9.6 kJ/mol). The reason for this is that when two LJ methanes are brought into coincidence the attractive interactions with the solvent remain so that this process more closely resembles removal of the repulsive Weeks−Chandler− Andersen (WCA) core57 of a single methane from solution (Appendix). This overlap free energy appears reasonable based on visual extrapolation of the solvent-mediated PMF reported in Figure 2 down to r = 0. Based on the continuity of the cavity correlation function and its first and second derivatives for hard solutes into the excluded volume region, Pratt and Chandler21 suggested the solvent contribution to the PMF can be interpolated down to the overlap limit using a cubic function w w ωmm (r )|cubic = ωmm (0) + ar + br 2 + cr 3

(4)

Here ωwmm(0) is the methane overlap free energy, independently evaluated from the simulations discussed in the Appendix, while a, b, and c are parameters that smoothly join the excluded volume and the simulation observed region just outside repulsion. Following ideas put forward by Ashbaugh and Pratt,4,58 we propose that eq 4 can be knitted together with the simulation results for the solvent interaction (eq 3) using a hybrid functional form w w w ωmm (r ) = f (r )ωmm (r )|cubic + [1 − f (r )]ωmm (r )|sim

Figure 3. Water-mediated contribution to the methane−methane PMF obtained by interpolation between the simulation results and methane-pair overlap free energies following eq 5. (a) Solvent contribution of the PMF at temperatures of 280, 300, 320, 340, and 360 K at 1 bar pressure. (b) Solvent contribution of the PMF at pressures of 1, 1000, 2000, 3000, 4000, and 5000 bar at a temperature of 300 K. Symbols are defined in the figure legends. The inset figures show PMF details for radii greater than 4 Å.

(5)

In this expression, f(r) is a switching function that interpolates between eq 4 and the simulation results over the lower and upper limits of R1 and Ru ⎧ 1, r ≤ Rl , ⎪ ⎪ (R − r )2 (R u − r )3 − 2 , Rl < r < R u , f (r ) = ⎨ 3 u 2 (R u − R l)3 ⎪ (R u − R l ) ⎪ 0, r ≥ Ru ⎩

to the PMF at overlap deepens to an even greater extent (Figure 3b). With increasing pressure the solvent PMF near contact at ∼4 Å is effectively fixed, while the barrier between the contact and solvent separated minima is drawn to lower separations (Figure 3b inset), in difference to the results at fixed pressure. As with the RDFs reported in Figure 1, the present results are in general agreement with previously reported trends observed for the temperature and pressure dependence of the water-mediated contribution to the methane−methane PMF.25,28,30,31 In difference to those reports, however, Figure 3 demonstrates that with high quality information regarding the PMF outside the inaccessible excluded volume and the limiting value at complete methane overlap, eq 5 is able to construct a reasonable interpolation between these two limits that varies in a smooth, systematic manner. Given the range of thermodynamic states simulated, it should be possible to break the solvent contribution to the PMF into its constituent thermodynamic components that fully describe methane−methane interactions over a wide range of temperatures and pressures. Under the assumption that second derivative properties of the interaction free energy (i.e., the heat capacity, compressibility, and expansivity) are independent of temperature and pressure, water’s contribution to the methane−methane PMF can be described using the following thermodynamic expression:18

(6)

The parameters a, b, and c in eq 4 are obtained here by minimizing the mean square difference between eq 5 and ωwmm(r)|sim over the range R1 < r < Ru. Figure 2 demonstrates that eq 5 smoothly interpolates between the complete overlap limit and simulation regime assuming R1 = 3.5 Å and Ru = 5.0 Å. Encouragingly, the interpolated function agrees quantitatively with Levy et al.’s simulation results down to 2.4 Å, supporting the accuracy of our interpolation procedure. The solvent contribution to the PMF from r = 0 Å to 12 Å as a function of temperature at 1 bar and pressure at 300 K are reported in Figure 3, corresponding to the same state points in Figure 1. As can be seen, this hybrid functional form provides an excellent description down to complete overlap, showing systematic changes in the solvent contribution to the interaction free energy well into the excluded volume region. At 1 bar (Figure 3a), the solvent PMF from overlap to contact deepens with increasing temperature, while the barrier between the contact and solvent separated minima (or RDF peaks in Figure 1) shifts to higher separations (Figure 3a inset). Over the simulated pressure range at 300 K, the solvent contribution 6283


Article


Figure 4. Thermodynamic components of the water-mediated contribution to the methane−methane PMF determined from fitting of eq 7 over the 50 simulated state points. The reference temperature and pressure used for the fitting is T0 = 300 K and P0 = 1 bar. (a) Solvent contribution to the methane−methane interaction enthalpy, hwmm,0, at the reference temperature and pressure. (b) Solvent contribution to the methane−methane interaction entropy, swmm,0, at the reference temperature and pressure. (c) Solvent contribution to the methane−methane interaction heat capacity, cwmm. (d) Solvent contribution to the methane−methane interaction volume, υwmm,0, at the reference temperature and pressure. (e) Solvent contribution to the methane−methane interaction isothermal compressibility, κwmm = −∂υwmm/∂P|T. (f) Solvent contribution to the methane−methane interaction thermal expansivity, αwmm = −∂υwmm/∂T|P. Each of these functions from 0 to22 Å are reported in tabular form in the Supporting Information.

can be observed from the excellent reproduction of the methane−methane RDFs reported in Figure 1. While methane−methane interaction enthalpies, entropies, heat capacities, volumes, and compressibilities have been previously reported, we believe this is the first time they have all been collectively determined over such a broad range of temperatures and pressures. This also permits evaluation of the interaction expansivity as well. Moreover, the present results allow direct comparison between solvent contribution to methane−methane interactions at contact (rcon ∼ 3.9 Å) and complete overlap, which effectively corresponds to removal of the repulsive WCA volume of a single methane from solution as noted above. For the remainder of this section, we consider the impact of temperature and pressure on the enthalpic and entropic components of the water-mediated interaction free energy. These thermodynamic functions are obtained through appropriate derivatives of eq 7 as

w w w ωmm (r , T , P) = hmm,0 (r ) − Tsmm,0 (r ) w w + cmm (r )[(T − T0) − T Sn(T /T0)] + υmm,0 (r )(P − P0) 1 w w − κ mm(r )(P − P0)2 + αmm (r )(T − T0)(P − P0) (7) 2

Here T0 and P0 correspond to the reference temperature and pressure, respectively, chosen to be 300 K and 1 bar. The six functional components of the solvent-mediated PMF are hwmm,0(r), the solvent contribution to the interaction enthalpy at T0 and P0; swmm,0, the solvent contribution to the interaction entropy at T0 and P0; cwmm(r), the solvent contribution to the interaction heat capacity; υwmm,0(r), the solvent contribution to the interaction volume at T0 and P0; κwmm(r), the partial isothermal compressibility of interaction; and αwmm(r), the partial thermal expansivity of interaction. In difference to standard thermodynamic definitions of the compressibility and expansivity, which are determined by the pressure and temperature derivatives of the logarithm of the volume, the partial compressibility and expansivity derivatives are determined here by the absolute derivative (i.e., κwmm(r) = −∂υwmm(r)/ ∂P|T and αwmm(r) = −∂υwmm(r)/∂T|P) since the methane− methane interaction volume can change sign as a function of distance. The thermodynamic components of the water-mediated interaction free energy, obtained by least-squares fitting of eq 7 to the solvent-mediated interaction free energies over all 50 simulated state points, are reported in Figure 4. These functions are also reported in tabular form from r = 0 Å to 22 Å in the Supporting Information. The quality of the fit is excellent, as

w (r , T , P ) = − T 2 hmm

w (r , T , P)/T ∂ωmm ∂T

= =

w (r ) hmm,0

−

1 w w T0(P − P0) κ mm(r )(P − P0)2 − αmm 2

T , P) + +

w (r , Tsmm

P

w (r , ωmm

w (r )(T cmm

T , P)

w (r )(P − P0) − T0) + υmm,0

(8a)

and 6284


Article

The Journal of Physical Chemistry B w smm (r , T , P ) = −

=

w smm,0 (r )

+

w (r , T , P ) ∂ωmm ∂T

w cmm Sn(T /T0)

−

interacting pair. The product of the temperature and the entropy at contact is nearly independent of temperature near contact; however, the entropy at overlap has a very strong temperature dependence. In particular, while the overlap entropy strongly favors the interaction at low temperature, the overlap entropy becomes increasingly repulsive with increasing temperature appearing to change sign at ∼362 K at 1 bar. As a result of the distinct temperature dependencies of the entropy at contact and overlap, the net entropic contribution to the interaction free energy appears to display a favorable minimum value in the vicinity of the contact separation for temperatures above 320 K. Taken together, the enthalpic and entropic components of the solvent-mediated interaction free energy as a function of temperature indicate that the molecular-scale interactions that drive two methanes toward contact (r = rcon) can be independent of the sign of their values upon desolvation of a methane core (r = 0 Å), as illustrated in the case of the entropy (Figure 5b), which displays distinct temperature dependencies between these two characteristic separations. These observations challenge the validity of interpretations of hydrophobic interactions in terms of hydrophobic hydration thermodynamics, which in turn are used to justify surface area based correlations parametrized to hydration free energies to model molecular scale hydrophobic interactions.59 The enthalpic and entropic components of the solvent interaction free energy as a function of pressure at 300 K are reported in Figure 6. The repulsive enthalpic barrier between contact and the solvent separated minimum is pushed down with increasing pressure, falling to a value of nearly zero at the

P

w (r )(P αmm

− P0)

(8b)

The enthalpic and entropic components of the solvent interaction free energy as a function of temperature at 1 bar are reported in Figure 5. The interaction enthalpy (Figure 5a) is

Figure 5. Enthalpic (a) and entropic (b) components of the watermediated contribution to the methane−methane PMF at temperatures of 280, 300, 320, 340, and 360 K at a pressure of 1 bar. Symbols are defined in the figure legend.

comparable to that reported by Garde et al.30 and Dias and Chan,31 with a broad repulsive barrier in the vicinity of the contact minimum (rcon = 3.9 Å). This barrier has been attributed to lost dispersion interactions between methane and water as two methanes are brought into contact. As the methanes are squeezed into the excluded volume region, however, attractions are regained and the enthalpy becomes attractive. Since the impact of methane−water attractions largely cancels in the evaluation of the overlap free energy (Appendix), the net enthalpic attraction is not a result of methane−water interactions but increased favorable interactions of water with itself. Indeed, the hydration enthalpy of repulsive methane-sized solutes (hard spheres and WCA particles) has been observed to be positive at ambient temperatures and above.4,12 It is not surprising then that desolvation of a single methane’s repulsive core upon complete overlap is enthalpically favorable. The repulsive enthalpy at contact is overcome by a more favorable, positive entropy at contact, leading to a net attraction (Figure 5b). This has been interpreted as a signature of the hydrophobic interaction resulting from the release of configurationally restricted water molecules in the first shell of the

Figure 6. Enthalpic (a) and entropic (b) components of the watermediated contribution to the methane−methane PMF at pressures of 1, 1000, 2000, 3000, 4000, and 5000 bar at a temperature of 300 K. Symbols are defined in the figure legend. 6285


Article


improved convergence. Sample smeared methane−methane RDFs at 300 K and 1 and 5000 bar are reported in Figure 7

largest pressures simulated (Figure 6a). While both Garde et al.30 and Dias and Chan31 observe a slight shift in the repulsive enthalpic barrier to shorter separations with increasing pressure, they did not report a suppression of the barrier magnitude with increasing pressure. We attribute these differences most likely to differences in the water models used in each study, TIP3P for Garde et al.30 and TIP4P for Dias and Chan,31 compared to TIP4P/2005 presently. Following the increased favorability of enthalpic contact interactions with increasing pressure, the overlap enthalpy becomes increasingly more attractive with increasing pressure as well. While the interaction entropy is largely favorable at 1 bar, the entropy at contact becomes less attractive with a growing barrier between the solvent separated and contact separations at higher pressures (Figure 6b). Dias and Chan31 reported comparable changes in the entropic interaction with increasing pressure. In difference to contact, the overlap entropy becomes increasingly favorable with increasing pressure, further highlighting the distinctions between interactions at contact and core desolvation noted above. While a number of thermodynamic second derivative properties have been assumed to be constant in fitting the PMF to eq 7 over a broad range of state points, this approximation only introduces slight differences in the properties reported in Figures 5 and 6 when finite difference derivatives of the simulation results are used instead with the trends described above preserved. More advantageously, eq 7 permits a ready interpolation between state points that cannot be easily obtained from finite difference derivatives. Osmotic Second Virial Coefficient of Methane in Water. The derivation of the McMillan−Mayer expression for the osmotic second virial coefficient (eq 1) presumes methane correlations are evaluated in an open ensemble,35,38 while presently we employed the closed isothermal−isobaric ensemble. Although RDF normalization differences between ensembles can be overcome, the evaluation of the B2 integral is still challenging due to long-range correlations that make convergence of eq 1 difficult to ascertain over finite integration domains. Recognizing that the choice of the molecular center used to measure correlations between species in eq 1 and associated Kirkwood−Buff integrals60 does not change the ultimate result, Lockwood and Rossky61 advocated using a continuous distribution of sites rather than more traditional nuclear centers to smear out RDFs. Smearing washes out packing correlations and improves convergence of the correlation function integrals. Here we smear methane−methane RDFs using a normalized isotropic step function centered on the carbon nucleus of each methane ⎧ 3/(4πλ 3), r ≤ λ , p(r , λ ) = ⎨ ⎪ r>λ 0, ⎩

Figure 7. Impact of correlation smearing applied to the methane− methane radial distribution function determined from simulation at 300 K and pressure of 1 and 5000 bar. Smearing was carried out via eq 10 using the isotropic, spherical smearing function (eq 9) with a radius of λ = 2.5 Å. A smearing radius of λ = 0 Å indicates the unsmeared distribution function Symbols are defined in the figure legend. The results at 5000 bar are shift up by +1 for clarity.

using λ = 2.5 Å. Despite the fact that the unsmeared (λ = 0 Å) RDF peaks grow with increasing pressure, the smeared RDFs are relatively featureless with a single repulsive dip at complete overlap followed by a mollified attractive correlation peak centered near the contact separation of the unsmeared RDFs. While eq 1 is independent of λ, we adopt a smearing radius of λ = 2.5 Å in the following discussion since the radius-of-gyration of the doubly smeared distribution is 2.74 Å (i.e., ⟨|r′ + r″|2⟩1/2 = 2.74 Å), comparable to the diameter of the water molecules that pack between the methanes and give rise to the packing correlations. We note that choosing significantly larger values of λ can lead to problems in the normalization of RDFs evaluated from finite sized simulation boxes with periodic images. To overcome normalization differences between RDFs evaluated in open versus closed ensembles, Schnell, Vlugt, Kjelstrup and co-workers64−68 considered how occupancy fluctuations in an open observation volume embedded within an even larger closed ensemble scale with the subvolume’s size. It was reasoned that fluctuations within the subvolume asymptotically approach the open system thermodynamic limit for an infinite-sized observation volume. For a spherical observation volume of radius R within a larger closed ensemble, the McMillan-Mayer osmotic second virial coefficient integral can subsequently be reformulated as64

⎪

(9)

where λ is the smearing radius. The smeared RDF is evaluated by the double convolution integral

B2 (R , λ) = 2π

∫|r | 0), while the thick black contour indicates the theta condition (B2 = 0). The thick dashed lines indicate the conditions ∂B2/∂P|T = 0 (purple dashed line) and ∂B2/∂T|P = 0 (green dashed line). These two lines intersect at a predicted saddle point at T = 338 K and P = 2690 bar, just outside the range of simulated state points.

Table 1. Least Squares Fitted Parameters of Eq 14 to the Osmotic Second Virial Coefficients Obtained from the Simulation Radial Distribution Functions and Those Fitted to the Simulation Radial Distribution Functions Using Eq 7 β1 β2 β3 β4 β5 β6

3

(cm K/mol) (cm3/mol) (cm3/mol) (cm3 K/(mol bar)) (cm3 K/(mol bar2)) (cm3/(mol bar))

simulated RDF’s

fitted RDF’s

64533 −215.71 1185.3 −4.2046 −1.4790 × 10−5 0.12255

67939. −226.69 1446.9 −4.0791 −1.2179 × 10−4 0.12376

transition with increasing temperature from B2s that become more attractive with increasing pressure to B2s that become more repulsive with increasing pressure can be linked to the hyperbolic nature of the underlying functional description, tied to the point at which ∂B2/∂P|T = 0. This condition is satisfied for points along the curve

coefficients obtained from eq 14 and those from the simulated or fitted RDFs are 10.6 cm3/mol and 0.5 cm3/mol, respectively. The larger error in the fit of the B2’s obtained from the simulated RDFs reflects the greater noise in their virial coefficients (Figure 9), especially with increasing pressure. Nevertheless, the agreement between the βi’s obtain from fitting to either data set is good (Table 1) with a root-meansquare difference between osmotic virial coefficients obtained from either parameter set of 1 cm3/mol over all the temperatures and pressures simulated. Subsequently, below we confine our discussion below to predictions of eq 14 using parameters describing virial coefficients from the fitted RDFs. The osmotic second virial coefficient of methane in water as a continuous function of temperature and pressure from eq 14 is shown as a contour plot in Figure 11. While B2 is predominantly attractive over most temperatures and pressures, an island of repulsive values is observed at low temperatures (T < 300 K) and pressures (P < 600 bar). This plot reveals an underlying hyperbolic form to the temperature and pressure dependence of B2, similar to descriptions of the transfer free energy of organic species between their pure liquid and water71 and distinct from elliptical free energy descriptions typically used to describe the stability of proteins.18,34 The observed

P = P0 −

β4 2β5

−

β6 2β5T

(15)

Equation 15 marks a nearly vertical boundary between portions of the curve where B2 decreases/increases with increasing pressure (Figure 11). At ambient pressure, the transition temperature is 333 K, which lies between the simulated temperatures of 320 and 340 K, where this transition was anticipated to occur based on our description of Figure 9. Extrapolating just outside the range of simulated temperatures and pressures, the curve ∂B2/∂P|T = 0 intersects the curve ∂B2/ ∂T|P = 0 at the state point T = 338 K and P = 2690 bar (Figure 12), where the determinant of the Hessian matrix is negative, confirming this is a hyperbolic saddle point. It is worth noting that B2 passes through a maximum at fixed pressure when ∂B2/ ∂T|P = 0 where attractions between methanes at elevated pressure are weakest along an isobar. We do not explore this behavior in greater detail here, however. A question that arises is “How is the boundary between B2s that decrease/increase with increasing pressure tied to the methane−methane PMF?” Interestingly the transition temperature at 333 K appears to coincide with the crossing temperature between the osmotic virial coefficient at 1 bar 6289


Article


The partial molar volume of methane obtained from simulation as a function of temperature at 1 bar is compared to an experimental correlation for the volume developed by Shock and co-workers76,77 in Figure 10b. The volume of methane in TIP4P/2005 water compares very favorably with the experimental correlation over a wide range of temperatures, differing from experiment by at most 1 cm3/mol at the highest temperatures simulated. While TIP4P water accurately captures the volume of methane near 300 K, the volume changes much more dramatically with temperature than experiment with a thermal expansivity of methane nearly twice as large. TIP3P water displays the worst agreement with experiment, consistently overpredicting methane’s volume62 over the entire temperature range simulated with an expansivity greater than even that in TIP4P. While not definitive, we may anticipate that TIP3P and TIP4P may exhibit shifts in the B2 transition down to lower temperatures due to their enhanced methane expansivity compared to both experiment and TIP4P/2005. This is consistent with the observation that many of the thermodynamic signatures of methane hydration in other water models are shifted to lower temperatures due to most commonly used water models exhibiting their temperatures of maximum density well below that of real water. The TIP3P and TIP4P water models exhibit density maxima at −74.2 °C and −19.2 °C, respectively, while TIP4P/2005 water (Tmd = 2.6 °C) compares quite favorably with experiment (Tmd = 4 °C).72

Figure 12. Correlation between the osmotic second virial coefficient and the maximum peak value of the methane−methane RDF at contact. The data plotted is obtained from the RDFs fitted to eq 7. Symbols are defined in the figure legend. The red arrow indicates the decrease in B2 with increasing pressure at 320 K. The purple arrow indicates the increase in B2 with increasing pressure at 360 K.

and the gas phase virial coefficient noted in our discussion of Figure 9a. We believe this observation is likely a coincidence, however, since the gas phase virial coefficient is independent of solvent interactions. Solvent effects rather are manifest in packing correlations in the hydrated RDFs. We have plotted B2 versus the maximum peak height in the methane−methane RDF (gmax mm ) (Figure 12), following Koga’s suggestion that these two quantities are correlated.41 While we do observe that B2 tends to decrease with increasing values of the primary contact peak height, reflecting the influence of increasing attractions, we do not observe a universal correlation which all the data fall upon. Nevertheless, as indicted by the arrows (Figure 12), the contact peak height tends to increase with increasing pressure at 320 K and decrease with increasing pressure at 360 K. At 340 K the peak barely increases with pressure, while B2 weakly increases. While the correlation between gmax mm and B2 is not universal, we do conclude the variation in B2 with pressure largely follows gmax mm . max To obtain an estimate of the point at which gmm is independent of pressure, we consider the methane−methane interaction volume at rcon = 3.9 Å. From Figure 4d,f, at methane−methane contact the interaction volume is υwmm,0(rcon = 3.9 Å) = −1.18 cm3/mol and expansivity is αwmm(rcon) = 0.0294 cm3/(mol K). When the interaction volume is negative/ positive we expect the contact peak to grow/shrink with increasing pressure. At the temperature T = T0 − υwmm,0(rcon)/ αwmm(rcon) = 340 K the interaction volume at contact is 0 and the peak magnitude is expected to be indifferent to increasing pressure. This temperature is in reasonable agreement with the transition temperature calculated from eq 15 above (333 K) where the pressure response changes sign at ambient pressure. It has been previously documented that the response of methane hydration thermodynamics is quite sensitive to the water model used.69,72−75 It may be asked, how sensitive is the observed transition from methane becoming more to less attractive with increasing pressure dependent on the water model? While a full investigation of this question would be computationally intensive, the discussion above suggests the response depends on the accuracy of water model at capturing the thermal expansivity of methane and methane pairs in water.

■

SUMMARY An empirical scheme for bridging water’s contribution to hydrophobic pair interactions over all separations from overlap to complete separation was proposed and applied to examine the methane−methane PMF over a broad range of temperatures and pressures. This bridge into the excluded volume region smoothly interpolates between water’s contribution to the methane−methane PMF evaluated by passive molecular simulation and the methane overlap free energy evaluated from particle intersection techniques using a cubic polynomial. The thermodynamic components of the interaction free energies were subsequently extracted over the simulated temperature and pressure range by fitting to a thermodynamic expansion (eq 7) to obtain the solvent-mediated enthalpy, entropy, heat capacity, volume, compressibility, and expansivity of interaction. Comparing methane hydrophobic interactions at contact (r = rcon) to the desolvation of the repulsive core of a single methane at overlap (r = 0), we observe distinct differences between the enthalpic and entopic contributions to the interaction free energy and their associated temperature and pressure dependencies. For example, at 300 K and 1 bar hydrophobic contact is enthalpically unfavorable while removal of the solute’s repulsive core is favorable. In addition, while hydrophobic contact is overwhelming entropically favorable at 1 bar pressure, the net contribution to the hydrophobic contact free energy is largely independent of temperature (i.e., −Tswmm(rcon) ≈ constant) at 1 bar. This stands in contrast to the desolvation entropy, which is large and favorable at low temperature and becomes increasingly less favorable with increasing temperature. These observations strongly question the validity of using simple area based correlations relating individual hydrophobe hydration to model molecular-scale hydrophobic interactions. The quality of the methane−methane RDFs obtained in the first half of this work subsequently permitted evaluation of methane’s osmotic second virial coefficient over the simulated 6290


Article

The Journal of Physical Chemistry B range of temperatures and pressures. We found that the ensemble independent extrapolation procedure proposed by Schnell, Vlugt, Kjelstrup, and co-workers64−68 yielded osmotic virial coefficients from isothermal−isobaric ensemble simulations that varied in a systematic manner with changes in the thermodynamic state and agreed with an independent thermodynamic consistency check based on the concentration dependence of methane’s solvation free energy. Moreover, the osmotic virial coefficients obtained from RDFs fitted to eq 7 were in excellent quantitative agreement those obtained directly from the simulated RDFs, but were smoother, allowing for a more accurate examination of the temperature and pressure dependence of B2. At ambient pressure, methane’s B2 was found to be endothermic in agreement with recent simulation studies, decreasing from positive to negative values with increasing temperatures, crossing the theta point at 300 K. With increasing pressure, however, B2 was found to either become more attractive (decrease) or repulsive (increase) depending on the temperature. The transition between these two regimes follows the dependence of the first peak of the methane−methane RDF on pressure, which in turn was linked to the volume and expansivity of interaction at contact. The underlying functional dependence of B2 on temperature and pressure was found to be hyperbolic, similar to the transfer of organic solutes between water and their pure liquids. Considering the thermal response of the volumetric properties of methane in different water models suggests that the susceptibility of B2 to pressure changes may be shifted in temperature. The accuracy of an individual water model at capturing water’s experimental equation-of-state properties is therefore of potential concern in accurately predicting the interaction properties between methane in solution, in agreement with previous studies of the hydration properties of nonpolar gases.

■

Figure A1. Excess chemical potential of methane and methane-like solutes in water as a function of temperature and pressure used to obtain the hydration free energy of bringing two methanes from infinitely far away to direct overlap. (a) Excess chemical potentials as a function of temperature at 1 bar. (b) Excess chemical potentials as a function of pressure at 300 K. The symbols are defined in the figure, while the solid black lines indicate fits to eq A2. Fit parameters are reported in Table A1. Error bars are smaller than the figure symbols.

APPENDIX

Methane Overlap Free Energy

For a united-atom description of methane, the solvent contribution to the PMF in the overlap limit (r = 0) corresponds to the excess free energy of hydrating a LJ site with a methane−oxygen interaction well-depth twice that of a single methane (i.e., 2εmo) less two times that of an individual methane with a normal interactions, w ωmm (0) = μmex |2εmo − 2μmex |εmo

The methane pair overlap free energy (Figure A1) is even more favorable than that associated with removing a single methane from solution (i.e., ωwmm(0) < −μex m |εmo). Assuming that methane−water attractive interactions only play a perturbative role in determining the hydration free energy, we may expect methane−water attractions in the overlap limit to be equal to twice that of two methanes distantly separated from one another. As a first approximation, then the overlap free energy should be comparable to that of removing a repulsive WCA particle with it’s repulsive core dictated by the single methane− WCA oxygen well-depth from solution (i.e., ωwmm(0) ≈ −μex m |εmo ). As can be seen in Figure A1 this approximation predicts a slightly more negative (favorable) overlap free energy than calculated from eq A1, suggesting methane−water attractions play little role in determination of the solvent-mediated interactions in for strongly overlapping methanes. Near quantitative agreement is observed for two overlapping WCA methanes (i.e., ωwmm(0) ≈ WCA ex WCA μex m |2εmo − 2μm |εmo , results reported in Table A1), confirming the minimal role of attractions in the final stages of methane overlap. The excess hydration free energies of methane, methane with twice the methane-oxygen interaction well-depth, the overlap

(A1)

These free energy components are readily obtained by analyzing our configurations generated during simulations of 40 methanes in a bath of 4000 water molecules using Widom’s test particle insertion (eq 2). The hydration free energy of the single methane is unfavorable for dissolution and closely follows the free energy of a methane at infinite dilution (Figure A1), though with slight differences resulting from methane correlations encapsulated in the osmotic second virial coefficient (Figure 10a). Nevertheless, the methane free energy at a ∼1% mole fraction exhibits the same enthalpic, entropic, and heat capacity signatures associated with methane solvation at infinite dilution. Interestingly, the hydration free energy of a methane with twice the normal methane−oxygen well-depth is negative and favorable for dissolution (Figure A1). Nevertheless, while being soluble, this “attractive” methane still displays many of the characteristic temperature dependencies of hydrophobic hydration. 6291


Article


Table A1. Least Square Fit Parameters of Eq A2 Fitted to the Excess Chemical Potentials of a Methane (εmo), a Methane with Twice the Well-Depth (2εmo), and the Overlap Energya OPLS-UA Methane μex m |εmo μex m |2εmo ωwmm(0) WCA Methane WCA μex m |εmo WCA μex m |2εmo WCA ωex mm(0)|

a

hex m,0 (kJ/mol)

sex m,0 (J/(mol K))

cex m (J/(mol K))

3 υex m,0 (cm /mol)

3 κex m (cm /(mol bar))

3 αex m (cm /(mol K))

−8.75

−59.8

188

36.8

−1.40 × 10−3

3.07 × 10−2

−4

1.86 × 10−2

−32.1

−81.4

176

32.8

−14.6

38.2

−201

−40.8

16.1

−42.6

220.

47.7

18.3 −13.9

−44.4 40.8

216

54.1

−224

−41.4

−3.02 × 10

−3

2.50 × 10

−3.79 × 10−2

−2.97 × 10−3

5.03 × 10−2

−3

5.71 × 10−2

−3.27 × 10

−3

2.67 × 10

−4.34 × 10−2

Results are reported for the OPLS-UA methane and the repulsive WCA core of the OPLS-UA methane.

■

free energy, and the corresponding repulsive WCA cores are well described by the free energy expression18 ex ex − Tsm,0 + cmex[(T − T0) − T Sn(T /T0)] μmex (T , P) = hm,0 1 ex + υm,0 (P − P0) − κmex(P − P0)2 + αmex(T − T0)(P − P0) 2 (A2)

As for eq 7, T0 and P0 correspond to the reference temperature and pressure, chosen here to be 300 K and 1 bar. In this expression, hex m,0 is the excess hydration enthalpy at T0 and P0, ex sex m,0 is the excess hydration entropy at T0 and P0, cm is the ex hydration heat capacity, υm,0 is the excess partial molar volume at T0 and P0, κex m is the solute’s absolute partial compressibility, and αex is the solute’s absolute partial expansivity. Least squares m fits of eq A2 are compared to the simulation free energies as a function of temperature and pressure in Figure A1. The fitted thermodynamic parameters in eq A2 are reported in Table A1. The overlap free energies determined here are used to establish the zero methane separation (r = 0) limit for the cubic interpolative PMF bridge (eq 4).

■

ASSOCIATED CONTENT

S Supporting Information *

The thermodynamic components of water’s contribution to the methane−methane PMF obtained from fits of eq 7 and displayed in Figure 4 are reported in tabular form in the Supporting Information document. The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.5b02056.

■

REFERENCES

(1) Murphy, K. P.; Privalov, P. L.; Gill, S. J. Common Features of Protein Unfolding and Dissolution of Hydrophobic Compounds. Science 1990, 247, 559−561. (2) Hummer, G.; Garde, S.; Garcia, A. E.; Paulaitis, M. E.; Pratt, L. R. Hydrophobic Effects on a Molecular Scale. J. Phys. Chem. B 1998, 102, 10469−10482. (3) Blokzijl, W.; Engberts, J. Hydrophobic Effects - Opinions and Facts. Angew. Chem., Int. Ed. 1993, 32, 1545−1579. (4) Ashbaugh, H. S.; Pratt, L. R. Colloquium: Scaled Particle Theory and the Length Scales of Hydrophobicity. Rev. Mod. Phys. 2006, 78, 159−178. (5) Pierotti, R. A. Scaled Particle Theory of Aqueous and Nonaqueous Solutions. Chem. Rev. 1976, 76, 717−726. (6) Bennaim, A.; Marcus, Y. Solvation Thermodynamics of Nonionic Solutes. J. Chem. Phys. 1984, 81, 2016−2027. (7) Tanford, C. The Hydrophobic Effect: Formation of Micelles and Biological Membranes, 2nd ed.; John Wiley and Sons: New York, 1980. (8) Pollack, G. L. Why Gases Dissolve in Liquids. Science 1991, 251, 1323−1330. (9) Lazaridis, T.; Paulaitis, M. E. Entropy of Hydrophobic Hydration - A New Statistical Mechanical Formulation. J. Phys. Chem. 1992, 96, 3847−3855. (10) Ashbaugh, H. S.; Paulaitis, M. E. Entropy of Hydrophobic Hydration: Extension to Hydrophobic Chains. J. Phys. Chem. 1996, 100, 1900−1913. (11) Ashbaugh, H. S. Entropy Crossover from Molecular to Macroscopic Cavity Hydration. Chem. Phys. Lett. 2009, 477, 109−111. (12) Garde, S.; Garcia, A. E.; Pratt, L. R.; Hummer, G. Temperature Dependence of the Solubility of Non-polar Gases in Water. Biophys. Chem. 1999, 78, 21−32. (13) Kauzmann, W. Some Factors in the Interpretation of Protein Denaturation. Adv. Protein Chem. 1959, 14, 1−63. (14) Dill, K. A. Dominant Forces in Protein Folding. Biochemistry 1990, 29, 7133−7155. (15) Kunugi, S.; Takano, K.; Tanaka, N.; Suwa, K.; Akashi, M. Effects of Pressure on the Behavior of the Thermoresponsive Polymer Poly(N-Vinylisobutyramide) (PNVIBA). Macromolecules 1997, 30, 4499−4501. (16) Graziano, G. On the Temperature-Induced Coil to Globule Transition of Poly-N-Isopropylacrylamide in Dilute Aqueous Solutions. Int. J. Biol. Macromol. 2000, 27, 89−97. (17) Paschek, D.; Nonn, S.; Geiger, A. Low-Temperature and HighPressure Induced Swelling of a Hydrophobic Polymer-Chain in Aqueous Solution. Phys. Chem. Chem. Phys. 2005, 7, 2780−2786. (18) Hawley, S. A. Reversible Pressure-Temperature Denaturation of Chymotrysinogen. Biochemistry 1971, 10, 2436−2442. (19) Dias, C. L.; Ala-Nissila, T.; Wong-ekkabut, J.; Vattulainen, I.; Grant, M.; Karttunen, M. The Hydrophobic Effect and Its Role in Cold Denaturation. Cryobiology 2010, 60, 91−99. (20) Athawale, M. V.; Goel, G.; Ghosh, T.; Truskett, T. M.; Garde, S. Effects of Lengthscales and Attractions on the Collapse of Hydro-

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS H.S.A. gratefully acknowledges extensive fruitful conversations with Prof. Lawrence Pratt on hydrophobic interactions and the evaluation of osmotic virial coefficients. This material is based upon work supported by the Gulf of Mexico Research Initiative and the National Science Foundation under the NSF EPSCoR Cooperative Agreement No. EPS-1003897 with additional support from the Louisiana Board of Regents. We also gratefully acknowledge computational support from the Louisiana Optical Network Initiative. 6292


Article

The Journal of Physical Chemistry B phobic Polymers in Water. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 733−738. (21) Pratt, L. R.; Chandler, D. Theory of the Hydrophobic Effect. J. Chem. Phys. 1977, 67, 3683−3704. (22) Smith, D. E.; Haymet, A. D. J. Free-Energy, Entropy, and Internal Energy of Hydrophobic Interactions - Computer-Simulations. J. Chem. Phys. 1993, 98, 6445−6454. (23) Hummer, G.; Garde, S.; Garcia, A. E.; Pohorille, A.; Pratt, L. R. An Information Theory Model of Hydrophobic Interactions. Proc. Natl. Acad. Sci. U.S.A. 1996, 93, 8951−8955. (24) Pangali, C.; Rao, M.; Berne, B. J. Monte-Carlo Simulation of the Hydrophobic Interaction. J. Chem. Phys. 1979, 71, 2975−2981. (25) Skipper, N. T. Computer Simulation of MethaneWater Solutions. Evidence for a Temperature-Dependent Hydrophobic Attraction. Chem. Phys. Lett. 1993, 207, 424−429. (26) Mancera, R. L.; Buckingham, A. D.; Skipper, N. T. The Aggregation of Methane in Aqueous Solution. J. Chem. Soc., Faraday Trans. 1997, 93, 2263−2267. (27) Rick, S. W. Heat Capacity Change of the Hydrophobic Interaction. J. Phys. Chem. B 2003, 107, 9853−9857. (28) Ludemann, S.; Abseher, R.; Schreiber, H.; Steinhauser, O. The Temperature-Dependence of Hydrophobic Association in Water. Pair versus Bulk Hydrophobic Interactions. J. Am. Chem. Soc. 1997, 119, 4206−4213. (29) Hummer, G.; Garde, S.; Garcia, A. E.; Paulaitis, M. E.; Pratt, L. R. The Pressure Dependence of Hydrophobic Interactions is Consistent with the Observed Pressure Denaturation of Proteins. Proc. Natl. Acad. Sci. U.S.A. 1998, 95, 1552−1555. (30) Ghosh, T.; Garcia, A. E.; Garde, S. Molecular Dynamics Simulations of Pressure Effects on Hydrophobic Interactions. J. Am. Chem. Soc. 2001, 123, 10997−11003. (31) Dias, C. L.; Chan, H. S. Pressure-Dependent Properties of Elementary Hydrophobic Interactions: Ramifications for Activation Properties of Protein Folding. J. Phys. Chem. B 2014, 118, 7488−7509. (32) Zhu, S.; Elcock, A. H. A Complete Thermodynamic Characterization of Electrostatic and Hydrophobic Associations in the Temperature Range 0 to 100 degrees C from Explicit-Solvent Molecular Dynamics Simulations. J. Chem. Theory Comput. 2010, 6, 1293−1306. (33) Payne, V. A.; Matubayasi, N.; Murphy, L. R.; Levy, R. M. Monte Carlo Study of the Effect of Pressure on Hydrophobic Association. J. Phys. Chem. B 1997, 101, 2054−2060. (34) Zipp, A.; Kauzmann, W. Pressure Denaturation of Metmyoglobin. Biochemistry 1973, 12, 4217−4228. (35) McMillan, W. G., Jr.; Mayer, J. E. The Statistical Thermodynamics of Multicomponent Systems. J. Chem. Phys. 1945, 13, 276−305. (36) Watanabe, K.; Andersen, H. C. Molecular-Dynamics Study of the Hydrophobic Interaction in an Aqueous-Solution of Krypton. J. Phys. Chem. 1986, 90, 795−802. (37) Kennan, R. P.; Pollack, G. L. Pressure-Dependence of the Solubility of Nitrogen, Argon, Krypton, and Xenon in Water. J. Chem. Phys. 1990, 93, 2724−2735. (38) Vafaei, S.; Tomberli, B.; Gray, C. G. McMillan-Mayer Theory of Solutions Revisited: Simplifications and Extensions. J. Chem. Phys. 2014, 141, 154501. (39) Alvarez, J. L.; Prini, R. F. Pressure-Dependence of the Solubility of Nitrogen, Argon, Krypton, and Xenon in Water - Comment. J. Chem. Phys. 1992, 96, 3357−3358. (40) Kennan, R. P.; Pollack, G. L. Pressure-Dependence of the Solubility of Nitrogen, Argon, Krypton, and Xenon in Water - Reply. J. Chem. Phys. 1992, 96, 3359−3360. (41) Koga, K. Osmotic Second Virial Coefficient of Methane in Water. J. Phys. Chem. B 2013, 117, 12619−12624. (42) Chaudhari, M. I.; Holleran, S. A.; Ashbaugh, H. S.; Pratt, L. R. Molecular-Scale Hydrophobic Interactions between Hard-Sphere Reference Solutes are Attractive and Endothermic. Proc. Natl. Acad. Sci. U.S.A. 2013, 110, 20557−20562.

(43) Chaudhari, M. I.; Sabo, D.; Pratt, L. R.; Rempe, S. B. Hydration of Kr(aq) in Dilute and Concentrated Solutions. J. Phys. Chem. B 2015, DOI: 10.1021/jp508866h. (44) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications, 2nd ed.; Academic Press: San Diego, CA, 2001. (45) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435−447. (46) Nosé, S. A Unified Formulation of the Constant Temperature Molecular-Dynamics Methods. J. Chem. Phys. 1984, 81, 511−519. (47) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695−1697. (48) Parrinello, M.; Rahman, A. Polymophic Transitions in SingleCrystals - A New Molecular-Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (49) Abascal, J. L. F.; Vega, C. A General Purpose Model for the Condensed Phases of Water: TIP4P/2005. J. Chem. Phys. 2005, 123, 234505. (50) Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. Optimized Intermolecular Potential Functions for Liquid Hydrocarbons. J. Am. Chem. Soc. 1984, 106, 6638−6646. (51) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. E. G. M. LINCS: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463−1472. (52) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N· log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (53) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, 1987. (54) Shing, K. S.; Chung, S. T. Computer-Simulation Methods for the Calculation of Solubility in Supercritical Extraction Systems. J. Phys. Chem. 1987, 91, 1674−1681. (55) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926−935. (56) Rick, S. W. Free Energy, Entropy and Heat Capacity of the Hydrophobic Interaction as a Function of Pressure. J. Phys. Chem. B 2000, 104, 6884−6888. (57) Weeks, J. D.; Chandler, D.; Andersen, H. C. Role of Repulsive Forces in Determining the Equilibrium Structure of Simple Liquids. J. Chem. Phys. 1971, 54, 5237−5247. (58) Ashbaugh, H. S.; Pratt, L. R. Contrasting Nonaqueous Against Aqueous Solvation on the Basis of Scaled-Particle Theory. J. Phys. Chem. B 2007, 111, 9330−9336. (59) Czaplewski, C.; Ripoll, D. R.; Liwo, A.; Rodziewicz-Motowidlo, S.; Wawak, R. J.; Scheraga, H. A. Can Cooperativity in Hydrophobic Association be Reproduced Correctly by Implicit Solvation Models? Int. J. Quantum Chem. 2002, 88, 41−55. (60) Kirkwood, J. G.; Buff, F. B. The Statisitcal Mechnical Theory of Solutions. I. J. Chem. Phys. 1951, 19, 774−777. (61) Lockwood, D. M.; Rossky, P. J. Evaluation of Functional Group Contributions to Excess Volumetric Properties of Solvated Molecules. J. Phys. Chem. B 1999, 103, 1982−1990. (62) Sangwai, A. V.; Ashbaugh, H. S. Aqueous Partial Molar Volumes from Simulation, and Individual Group Contributions. Ind. Eng. Chem. Res. 2008, 47, 5169−5174. (63) Priya, M. H.; Ashbaugh, H. S.; Paulaitis, M. E. Cosolvent Preferential Molecular Interactions in Aqueous Solutions. J. Phys. Chem. B 2011, 115, 13633−13642. (64) Schnell, S. K.; Englebienne, P.; Simon, J. M.; Kruger, P.; Balaji, S. P.; Kjelstrup, S.; Bedeaux, D.; Bardow, A.; Vlugt, T. J. H. How to Apply the Kirkwood-Buff Theory to Individual Species in Salt Solutions. Chem. Phys. Lett. 2013, 582, 154−157. (65) Schnell, S. K.; Liu, X.; Simon, J. M.; Bardow, A.; Bedeaux, D.; Vlugt, T. J. H.; Kjelstrup, S. Calculating Thermodynamic Properties from Fluctuations at Small Scales. J. Phys. Chem. B 2011, 115, 10911− 10918. 6293


Article

The Journal of Physical Chemistry B (66) Schnell, S. K.; Vlugt, T. J. H.; Simon, J. M.; Bedeaux, D.; Kjelstrup, S. Thermodynamics of a Small System in a Mu T Reservoir. Chem. Phys. Lett. 2011, 504, 199−201. (67) Schnell, S. K.; Vlugt, T. J. H.; Simon, J. M.; Bedeaux, D.; Kjelstrup, S. Thermodynamics of Small Systems Embedded in a Reservoir: A Detailed Analysis of Finite Size Effects. Mol. Phys. 2012, 110, 1069−1079. (68) Kruger, P.; Schnell, S. K.; Bedeaux, D.; Kjelstrup, S.; Vlugt, T. J. H.; Simon, J. M. Kirkwood−Buff Integrals for Finite Volumes. J. Phys. Chem. Lett. 2013, 4, 235−238. (69) Docherty, H.; Galindo, A.; Vega, C.; Sanz, E. A Potential Model for Methane in Water Describing Correctly the Solubility of the Gas and the Properties of the Methane Hydrate. J. Chem. Phys. 2006, 125, 074510. (70) Chandler, D. Introduction to Modern Statistical Mechanics; Oxford University Press: London, 1987. (71) Sawamura, S.; Nagaoka, K.; Machikawa, T. Effects of Pressure and Temperature on the Solubility of Alkylbenzenes in Water: Volumetric Property of Hydrophobic Hydration. J. Phys. Chem. B 2001, 105, 2429−2436. (72) Ashbaugh, H. S.; Collett, N. J.; Hatch, H. W.; Staton, J. A. Assessing the Thermodynamic Signatures of Hydrophobic Hydration for Several Common Water Models. J. Chem. Phys. 2010, 132, 124504. (73) Paschek, D. Temperature Dependence of the Hydrophobic Hydration and Interaction of Simple Solutes: An Examination of Five Popular Water Models. J. Chem. Phys. 2004, 120, 6674−6690. (74) Dyer, P. J.; Docherty, H.; Cummings, P. T. The Importance of Polarizability in the Modeling of Solubility: Quantifying the Effect of Solute Polarizability on the Solubility of Small Nonpolar Solutes in Popular Models of Water. J. Chem. Phys. 2008, 129, 7. (75) Krouskop, P. E.; Madura, J. D.; Paschek, D.; Krukau, A. Solubility of Simple, Nonpolar Compounds in TIP4P-Ew. J. Chem. Phys. 2006, 124, 016102. (76) Plyasunova, N. V.; Plyasunov, A. V.; Shock, E. L. Database of Thermodynamic Properties for Aqueous Organic Compounds. Int. J. Thermophys. 2004, 25, 351−360. (77) ORganic Compounds HYDration Properties Database. http:// webdocs.asu.edu/ (accessed April 30th, 2015).

6294


Temperature and Pressure Dependence of ... - ACS Publications

Recommend Documents