Hydrophobicity Varying with Temperature, Pressure, and Salt

Jan 22, 2018 - Temperature-, pressure-, and salt-concentration-induced variations in the solubility of small nonpolar solutes in aqueous solution and ...
1 downloads 8 Views 1MB Size
Article pubs.acs.org/JPCB

Cite This: J. Phys. Chem. B XXXX, XXX, XXX−XXX

Hydrophobicity Varying with Temperature, Pressure, and Salt Concentration K. Koga*,†,‡ and N. Yamamoto‡,§ †

Research Institute for Interdisciplinary Science, Okayama University, Okayama 700-8530, Japan Department of Chemistry, Faculty of Science, Okayama University, Okayama 700-8530, Japan



ABSTRACT: Temperature-, pressure-, and salt-concentration-induced variations in the solubility of small nonpolar solutes in aqueous solution and the corresponding variations in the solvent-induced pair attraction between such solute molecules are investigated. The variations in the solvation free energy of a solute and those in the solvent-induced pair attraction are well reproduced by a mean-field approximation in which the repulsive cores of solute molecules are treated as hard spheres and the mean-field energy of a solute molecule is taken to be the average potential energy that the solute molecule feels in solution. The mechanisms of variation in the solvation free energy and those of variation in the solvent-induced pair potential, with increasing temperature, pressure, and salt concentration, are clarified. Correlations between the solvation free energy and the solventinduced pair potential at a contact distance in temperature, pressure, and salt concentration variations are near linear in any mode of variation, but the slope of the linear relation is dependent on the mode of variation and is determined by a ratio of the solvation thermodynamic quantities characteristic of each mode of variation.



INTRODUCTION The low solubility of nonpolar species in water and the effective attractive forces between such molecules in water together are the defining element of the hydrophobic effect. The decrease in the solubility of nonpolar substances with increasing temperature is a characteristic feature of the hydrophobicity. It is commonly assumed that solute molecules in aqueous solution attract each other more strongly when a solute becomes less soluble in water, i.e., when it becomes more hydrophobic, in response to some change in the thermodynamic state such as a temperature rise, but this correlation has not been critically examined by experiment. The solubilities of typical hydrophobic solutes in aqueous solution have been measured with great accuracy,1−3 while the hydrophobic interactions between such molecules have been scarcely quantified, either as the pair correlation function h(r) or as the osmotic second virial coefficient B. Earlier studies reported that B for small hydrophobic molecules and molecules with hydrophobic groups are close to zero, indicating that the hydrophobic interaction is not as attractive as commonly assumed,4,5 and a recent study using Raman vibrational spectroscopy concluded that the hydrophobic interactions between small hydrophobic groups of alcohols in water are repulsive.6 The idea that the strength of hydrophobic attraction is correlated with the © XXXX American Chemical Society

solubility of that hydrophobic solute has yet to be tested experimentally. So far, much of our knowledge on hydrophobic interactions among typical nonpolar molecules comes from theoretical studies, which can give exact or accurate results for given models. Studies of exactly soluble lattice models show that there exist definite correlations between the solubility and the hydrophobic interaction.7−10 Computer simulation studies of realistic models and other theoretical approaches have given the solute−solute correlation function h(r), the potential of mean force at some distance, or the osmotic virial coefficient B.11−16 Most of the studies confirmed that the lower the solubility the stronger the effective attraction between solute molecules in solution, and a few found the opposite correlation when the temperature or pressure is high enough9 or when the concentration of a cosolvent such as alcohols increases.16 The solubility of a given species in aqueous solution may be changed not only by varying the temperature but also by varying the pressure, the chemical composition in solution, or other thermodynamic variables. The mechanism of the Special Issue: Benjamin Widom Festschrift Received: December 11, 2017 Revised: January 22, 2018 Published: January 22, 2018 A

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

examine the validity of the mean-field approximation in calculating the hydrophobic interaction between methane molecules. The underlying idea in the mean-field theory of liquids, sometimes referred to as the van der Waals picture, is transparent in the sense that it focuses only on the universal key features of simple liquids, simplifying the rest as much as possible. That is, the structure of liquids is primarily determined by the strong intermolecular repulsion at short distances and each molecule is in a deep but uniform negative background potential.19,20 Since it is built upon the common features of simple liquids, the mean-field approximation may well be applied to liquid mixtures or, more specifically, solvation and association of solutes in a solvent. However, when the main component of a solvent is water, the solvent structure is not simply determined by the short-ranged repulsive forces among molecules, either in its bulk region or in the vicinity of a solute molecule. Still, nonpolar solute molecules may be modeled as attracting hard spheres with respect to the solute−solvent intermolecular interactions. What we will examine here and in the following section is validity of the mean-field approximation to the solvation and association of solute molecules in aqueous solution. Another important subject on the validity of the mean-field theory concerns the small and large length scales of solute molecules in a solvent,25,26 but the present work deals with small solute molecules and the issue of length scales is not examined here. The first approximation in the mean-field model is that realistic pair potentials ϕi(r) between a solute molecule (particle) A and other molecules i are now taken to be those of attracting hard spheres, which is the simplest approximation to the short-ranged repulsive part of ϕi(r). Then, the potential energy Ψ(r) that molecule A feels at r in the solvent is +∞ if A and any other molecule i overlap; otherwise, Ψ(r) = ∑i ϕi(|r − ri|), which is mostly negative due to the weak, long-ranged attractive part of ϕi(r). Then, from the potential distribution theorem applied to the solute species,27

solubility change could be different depending on the mode of change of the thermodynamic state. Likewise, the effective strength of the solute−solute attraction may be changed in many ways, and the mechanism of variation could be different for different modes of changing the thermodynamic state. In this study, we consider the temperature, pressure, and salt concentration as varying thermodynamic variables. First, we examine the effect of each on the hydrophobic hydration and the hydrophobic interaction using computer simulation and applying a simple mean-field approximation. It turns out that the mean-field approximation describes very well the solvation and the association of small nonpolar molecules in aqueous solution. Second, we try to elucidate the nature and robustness of the correlation between the solubility and the hydrophobic attraction. To do so, we perform molecular dynamic simulations of aqueous solutions of methane and calculate the solvation free energy of the solute in solutions and the solute−solute pair correlation function, both as functions of the temperature, pressure, and NaCl concentration. We will see that in each mode of change of the thermodynamic state there exists a near-linear relation of the strength of hydrophobic attraction (as measured by the solvent-mediated potential W*(rc) of mean force at a solute−solute contact distance) to the solubility (as measured by the solvation free energy μA* of the solute). It will be shown that the rate of change of W*(rc) with μ*A is determined by a ratio of solvation thermodynamic quantities characteristic of the mode of variation in the thermodynamic state.



HYDROPHOBIC HYDRATION The solubility of a hydrophobic solute A in an aqueous solution β in equilibrium with another phase α (a gas, liquid, or solid) may be measured in terms of the equilibrium ratio (ρβA/ραA)eq of the number densities of A in the two coexisting phases. The ratio (ρβA/ραA)eq is called the Ostwald absorption coefficient when α is a gas phase, and when the gas is dilute, it is practically equal to Σ ≡ ρβA/zA, with the activity zA of A defined so as to be asymptotically equal to ραA as ραA → 0. The solvation free energy μ*A of A in a solution (now phase β) is the excess chemical potential μA − μidA , i.e., the difference between the chemical potential μA of A in the solution and the chemical potential μidA of the hypothetical ideal gas of A with the same density ρβA and temperature T as those in the solution. With that definition, μ*A is a logarithmic measure of the solubility Σ: μA* = −kT ln Σ

Σ = ω⟨e−Ψ / kT ⟩c

(2)

where ω is the probability that the solute molecule A (a test particle) will fit into the solvent at an arbitrary point in equilibrium configurations unaffected by the test particle and ⟨···⟩c is the conditional average taken over those configurations in which the test particle successfully fits. The factorization, eq 2, is exact for solute molecules with a hard-sphere core. The second approximation here is

(1)

⟨e−Ψ / kT ⟩c ≃ e−u / kT

For typical small hydrophobes, argon, methane, and oxygen, in water at 298 K and 1 bar, experimental values of Σ fall in a range from 0.031 to 0.03417 (meaning approximately 4 × 104 water molecules per solute in each solution). As the temperature rises, Σ decreases and, equivalently, μA*/kT increases until either the temperature of the solubility minimum or the boiling point is reached. The mean-field theory of liquids and the potential distribution theorem,18−20 initiated and formulated by Ben Widom in the 1960s, are the very cornerstone of liquid theory, based on which subsequent theories of liquids and inhomogeneous fluids of fundamental importance have been developed.21−24 In this section, we illustrate the validity of a simple version of the mean-field approximation in reproducing the solvation free energy of methane varying with the temperature, pressure, and salt concentration; in the next section, we will

(3)

where u is a mean potential energy that the solute molecule would feel in the solution. There are many ways of evaluating u with varying degrees of simplification. Here we identify u as u = ⟨Ψ(r)⟩ =

∑ ϕi(|r − ri|) i

(4)

which is the mean potential energy that the solute molecule A at r feels when the solvent molecules move subject to their mutual interaction and to their interaction with the solute molecule via the original pair potential ϕi(r). From eqs 1−3, the solvation free energy of solute A is written μA* = −kT ln ω + u B

(5) DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 1. Temperature dependences of the solvation thermodynamic quantities for the model methane in water: (a) the solvation free energy μA*, the solvation energy u, the solvation free energy −kT ln ω of the hard sphere, and the mean-field prediction μ*A (mean-field) = −kT ln ω + u; (b) μA*/kT, u/kT, −ln ω, and μA*/kT(mean-field) = −ln ω + u/kT. Standard deviations are smaller than the size of the symbols.

the value of μ*A at 298 K and 1 bar. Then, the temperature dependence of μA* is accurately given by the mean-field approximation. It has been remarked that the mean-field approximation of this kind should work well for small hydrophobic molecules.25,26 It is seen in Figure 1a that the temperature dependence of u is weaker than that of −kT ln ω: at 298 K, (∂u/∂T)p = 0.0219 kJ/(mol·K) and [∂(−kT ln ω)/∂T]p = 0.0581 kJ/(mol·K). Thus, the increase of μ*A with temperature is largely due to that of the solvation free energy of the hard sphere. It is also noticed that u has a near linear dependence on T: the line along the points of u in Figure 1a is the linear fit line

(In a homogeneous solution, μA*, ω, and u have no r dependence.) We will see below the validity of this simple mean-field approximation, eq 5, when applied to the solvation of a small hydrophobic molecule in pure water and salt solution. The model systems we study are realistic models for aqueous solutions of methane: the “methane” molecules are LennardJones (LJ) particles,28 the “water” molecules are those of the TIP4P/2005 model,29 and the Na+ and Cl− ions are LJ particles with charges.30 First, the solvation free energies μ*A of methane in aqueous solution at varying temperatures T, pressures p, and NaCl concentrations (molality) m are calculated by the testparticle insertion method. Equilibrium configurations of solvent molecules for which the test-particle insertion is carried out are generated by molecular dynamics (MD) simulations. The resulting μ*A are considered “exact” for the model systems. Second, we apply the mean-field approximation to evaluate μA*(mean-field) at the same thermodynamic states as in the exact calculations. The mean potential energy u of the solute molecule A in eq 4 is calculated by MD simulation of the system containing solute A in the solvent, while the probability ω that the test solute particle with a hard-sphere core fits in the solvent is given by the test-particle insertion method with a suitable choice of the hard-sphere diameters for each species. The details of the numerical methods and the model systems are described in the Computational Details section. Figure 1 shows how the solvation free energy μA* and other solvation thermodynamic quantities vary with temperature in the range 238−373 K at a pressure of 1 bar. Note that the plotted data of μ*A and μ*A /kT (black circles in Figure 1) are the results obtained without any approximation. They are very close to the experimental data for methane in water, because the model potential parameters are optimized to reproduce the solubility data.31 The monotonic increase of μ*A with temperature in this range and the increase of μ*A /kT up to its maximum at around 350 K are both well-known for methane in water. Plotted in the same figure are −kT ln ω the solvation free energy of the hard-sphere particle, u the mean potential energy of the LJ methane in solution, and −kT ln ω + u the mean-field solvation free energy μ*A (mean-field) in eq 5. The hard sphere diameter d = 3.17 Å of methane particles was chosen so that the mean-field prediction with that diameter fits

u = k(u0 + u1T )

(6)

with u0 = −2649 K and u1 = 2.633. Thus, the concavity of μ*A (T) comes solely from that of −kT ln ω. In Figure 1b, μ*A /kT, u/kT, and −ln ω are plotted against the temperature. Note that μA*/kT = −ln Σ is practically a logarithmic measure of the solubility (ρβA/ραA)eq and so the maximum of μ*A /kT seen at about 358 K corresponds to the solubility minimum of methane in water. It is to be noted that −ln ω is the largest at about 273 K, i.e., the probability ω that a hard sphere with diameter d chosen for methane particles fits in water is minimal at the low temperature. It is found that the values of −T ln ω are accurately fit by a quadratic function −T ln ω = c0 + c1T + c 2T 2

(7)

with c0 = −2976 K, c1 = 32.042, and c2 = −4.204 × 10−2 K−1. The curve along the data points of −kT ln ω in Figure 1a is the fit curve. Then, the temperature of minimum ω i s Tω = c0/c 2 , w h i c h i s 2 6 6 K , a n d −ln ω = c1 − 2 c0c 2 = 9.67 at T = Tω. At low temperatures below 298 K, −ln ω is nearly independent of T and is approximately equal to c1 − 2 c0c 2 while u/kT varies significantly and, since u is almost perfectly linear in T, u/kT ≃ u0/kT + u1/k with the coefficients given above. Then, at low temperatures around Tm, the temperature of maximum density μA* kT

≃−

( − u 0) + u1 + c1 − 2 c0c 2 T

(low temperatures) (8)

C

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 2. Pressure (p) and NaCl concentration (m) dependences of μA*/kT, u/kT, −ln ω, and μA*/kT(mean-field) = −ln ω + u/kT. Standard deviations are smaller than or as large as the size of the symbols: (a) pressure dependences at T = 298 K and m = 0; (b) NaCl concentration dependences at 298 K and 1 bar.

⎛ ∂ ln ω ⎞ ⎜ ⎟ ≃ 0 ⎝ ∂T ⎠ ρ

The only temperature dependence of μA*/kT is 1/T in the first term. As the temperature increases further, −ln ω decreases, because the probability ω that the test particle fits in the solvent increases as the density of water decreases, and μ*A /kT becomes maximal at some temperature below the boiling point of water. This can be seen from the expression of μA*/kT that follows from eqs 6 and 7 μA* kT



u0 + c0 + u1 + c1 + c 2T T

over the whole range of temperatures when the density of water is fixed to that at 298 K and 1 bar. It was also shown in our earlier study32 that −ln ω varies only by 9% as the temperature rises from 273 to 640 K (d = 3.44 Å and the fixed density is 1.002 g/cm3). Thus, the above relation (eq 11) for the probability ω of finding a small cavity in liquid water is expected to hold for a wide range of temperatures. The second term γVVhs * /kT in the right-hand side of eq 10 is negative when T < Tm, zero at T = Tm, and positive when T > Tm, and the magnitude is much greater than the first term except at temperatures around Tm, where both are null. From eqs 10 and 11, for a wide range of temperatures,

(9)

which is accurate over the whole range of temperatures. This approximate expression gives the temperature Ts of solubility minimum and the maximum of μA*/kT as Ts = μA* kT

u0 + c0 = 366 K, c2

* γVV hs ⎛ ∂ ln ω ⎞ ⎜ ⎟ ≃ ⎝ ∂T ⎠ p kT

= c1 + u1 − 2 c 2(c0 + u0) = 3.92

(12)

Figure 2 shows the pressure and NaCl concentration dependences of μ*A /kT, −ln ω, u/kT, and the mean-field μA*/kT at a temperature of 298 K. The pressure ranges from −1000 bar (a metastable state) to 6000 bar. The mean-field approximation, −ln ω + u/kT, agrees very well with μ*A /kT over the entire range of pressures. As seen in Figure 2a, the pressure dependence of u/kT is much smaller than that of −ln ω: linear fits to the data give [∂(u/kT)/∂p]T = 2.57 × 10−5 bar−1 and [∂(−ln ω)/∂p]T = 1.33 × 10−3 bar−1, respectively, and so

2

which are very close to the experimental results. When the solute−solvent pair potentials are of hard spheres, μ*A /kT = −ln ω, and from the thermodynamic identity (∂f/∂T)p = (∂f/∂T)ρ − (∂p/∂T)ρ(∂f/∂p)T with now f = μA*/kT, we have * γVV hs ⎛ ∂ ln ω ⎞ ⎛ ∂ ln ω ⎞ ⎜ ⎟ = ⎜ ⎟ + ⎝ ∂T ⎠ p ⎝ ∂T ⎠ ρ kT

(11)

(10)

where γV = (∂p/∂T)ρ is the thermal pressure coefficient of the solvent and

* V * = (∂μA*/∂p)T ≃ [∂( −kT ln ω)/∂p]T = V hs

Note that V* is the volume change of the system when a single methane molecule A is added and fixed to any given position at fixed pressure and temperature. From the linear fits, we find V* = 33.6 cm3/mol and Vhs * = 33.0 cm3/mol. The partial molecular volume vp is V* + kTχT with kTχT = 1.1 cm3/mol at 298 K, and therefore, we have vp = 34.7 cm3/mol for the present model of methane in water. The NaCl concentration dependences of μA*/kT, −ln ω, u/kT, and the mean-field prediction of μA*/kT are shown in Figure 2b. The temperature is again fixed to 298 K and the pressure to 1 bar. The salt concentrations (molalities) m are 0, 1, 2, and 3 mol/kg. The mean-field prediction −ln ω + u/kT vs m is again in good agreement with the μA*/kT vs m from the

⎡ ⎤ * = ⎢ ∂( −kT ln ω) ⎥ V hs ∂p ⎣ ⎦T

is the volume change of the system upon placing the hard sphere particle at a fixed position in the solvent, which is related to the partial molecular volume vhs of the hard sphere by vhs = Vhs * + kTχT. The first term (∂ ln ω/∂T)ρ in the right-hand side of eq 10 is identically null if the solvent is a hard-sphere fluid; it would be close to 0 if the solvent is a simple liquid, for the structure of solvent is determined by packing. Now the solvent is water with its characteristic hydrogen-bonding network structure. Nevertheless, we find that D

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 3. Potentials of mean force W(r) between methane molecules in aqueous solutions at 298 K: (a) response to pressure (p = −1000, 1, 1000, 2000, 3000, and 6000 bar and the solvent is pure water); (b) response to NaCl concentration (m = 0, 1, 2, and 3 mol/kg and p = 1 bar).

free energy calculation. Both −ln ω and u/kT are near linear functions of m, and the linear fits give the slopes 0.488 and 0.0262 (mol/kg)−1, respectively. Thus, the variation in the probability ω for cavity formation m is the main factor that determines the rate of change of μA* with m in the salt solution.

methane.34,35 Figure 3b shows the response of W(r) to increasing salt concentration m. As the NaCl concentration increases, the whole potential curve W(r) shifts to lower energies. The rate of variation of W(r) with m at fixed r is most notable at around the contact distance and gradually becomes less pronounced as the distance increases. The overall trend is in qualitative agreement with an earlier simulation result by Ghosh et al.36 With rising temperature, pressure, or NaCl concentration, the first minimum of the methane−methane potential of mean force W(r) becomes deeper; i.e., the hydrophobic interaction between methane molecules in water at a contact distance becomes more attractive, as illustrated in Figure 3 for the pressure and NaCl concentration dependences and as shown in Figure 3b in ref 13 for the temperature dependence. We shall now apply the mean-field approximation to the hydrophobic interactions or, more specifically, to the solvent-induced part W*(r) of the potential of mean force between methane molecules. Suppose the mean-field approximation, eq 5, holds for the solvation of a pair of A particles that are a fixed distance r. Then, the solvation free energy μD*(r) of the dimer is



HYDROPHOBIC INTERACTION The effective interaction between a pair of solute molecules in the solvent is described by the potential of mean force W(r), where r is the distance between the two molecules. The solvent-induced part W*(r) of the potential of mean force is the part remaining after the direct solute−solute interaction ψ(r) is subtracted. It is W*(r) that changes in response to variations in the thermodynamic state of the solution. One may regard W*(r) as the excess of the solvation free energy μ*D(r) of a pair of solute particles that are a fixed distance r apart over twice the solvation free energy μA* of a single solute particle, i.e., W *(r ) = μD*(r ) − 2μA*

(13)

The potential of mean force W(r) is related to the radial distribution function g(r) between solute particles as W(r) = −kT ln g(r) and the solvent-induced part as W*(r) = −kT ln[g(r)eψ(r)/kT]. The latter, W*(r), is independent of ψ(r), and so is g(r)eψ(r)/kT, the fact we made use of to evaluate g(r) of the LJ methane particles in model aqueous solutions from g(r) of the particles interacting with each other via the WCA repulsive part of ψ(r). Figure 3 shows how W(r) changes in response to increasing pressure and NaCl concentration at 298 K. The temperature dependence of W(r) for the same system was reported elsewhere.13 It is seen in Figure 3a that, as the pressure increases, the first minimum of W(r) becomes more negative, i.e., the effective interaction between methane molecules at a contact distance becomes stronger, the second minimum becomes deeper with the peak position shifting to shorter distances, and the height of the energy barrier between the first and second minima slightly increases. These results are consistent with the pressure dependence of g(r) reported by Ashbaugh et al.33 and also are qualitatively in accord with the pressure dependence of W(r) at a contact distance evaluated by Graziano;14 some earlier studies reported different pressure dependences of the first minimum of W(r) for models of

μD*(r ) = −kT ln ω(r ) + u(r )

(14)

where ω(r) is a probability that the dimer of hard spheres fits in the solvent; each sphere is identical to that chosen for the test methane particle and the center-to-center distance is fixed to r. u(r) is the average potential energy that the dimer feels when the solvent molecules move subject to their mutual interaction and to their interaction with each of the A particles of the dimer via the solute−solvent LJ pair potential; note that u(r) does not include the direct interaction ψ(r) between A particles. From eqs 5, 13, and 14, we have the mean-field approximation to W*(r) *(r ) + uattr(r ) W *(r ) = W hs

(15)

with W *(r ) = −kT ln

ω(r ) , ω2

uattr(r ) = u(r ) − 2u

(16)

We focus on W*(rc), where rc = 3.9 Å is the distance at which W(r) is minimal, or equivalently g(r) is maximal, at 298 E

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 4. Temperature dependences of W*(rc), Whs * (rc), uattr(rc), and W*(rc)(mean-field) at 1 bar. The error bars are 0.6 kJ/mol or less for Whs * (rc), 0.3 kJ/mol or less for uattr(rc), and smaller than the size of the symbols for W*(rc).

Figure 5. Pressure and NaCl concentration dependences of W*(rc)/kT, Whs * (rc)/kT, uattr(rc)/kT, and W*(rc)/kT(mean-field) at T = 298 K. The * (rc)/kT, 0.1 or less for uattr(rc)/kT, and smaller than the size of the symbols for W*(rc)/kT: (a) pressure error bars are 0.3 or less for Whs dependences at T = 298 K and m = 0; (b) NaCl concentration dependences at 298 K and 1 bar.

the potential energy difference uattr(rc) = u(rc) − 2u is nearly temperature-independent and merely shifts Whs * (rc) uniformly to the less negative W*(rc). Note that both the mean potential energies u of a single LJ methane particle and u(r) of the pair increase linearly with temperature in such a way that u(rc) − 2u is nearly constant. Shown in Figure 4b are those quantities plotted in Figure 4a divided by kT and so dimensionless. One can see that both W*hs(rc)/kT = −ln[ω(r)/ω2] and uattr(rc)/kT decrease with temperature, but the temperature dependence of uattr(rc)/kT is merely 1/kT. The variations of W*(rc)/kT with pressure and NaCl concentration and those of the other quantities are shown in Figure 5. The temperature is fixed at 298 K, and so kT’s are now constant unlike those in Figure 4b. For both the pressure and salt concentration dependences, the result of the mean-field approximation compares very well with W*(rc)/kT directly obtained from MD simulations. With increasing pressure, uattr(rc) decreases while W*hs(rc) increases, and as a result, the sum, W*(rc), slightly decreases (Figure 5a). Therefore, the mechanism of W*(rc) decreasing with pressure is different from that of W*(rc) decreasing with temperature. The pressure derivatives of uattr(rc) and W*hs(rc) are opposite in sign, and their magnitudes are of the same order. That may be the case for other hydrophobic solutes, too, and then, the pressure derivative of W*(rc) can be either positive or negative

K and 1 bar when the solvent is pure water. The contact distance rc varies only slightly with temperature, pressure, and salt concentration. Figure 4 displays the variations with temperature of W*(rc), the mean-field prediction W*hs(rc) + uattr(rc), and each of its terms. It is shown that W*(rc) decreases monotonically with increasing temperature and the mean-field approximation reproduces the temperature dependence reasonably well. In fact, the agreement is remarkable given that each of −kT ln ω(rc), −kT ln ω2, u(rc), and 2u in eq 16 is much larger in magnitude than W*(rc) and greatly changes with temperature. Also, the hard-sphere diameter is not optimized for the pair of solute particles, and thus the virtually perfect agreement at 298 K is notable. The mean-field approximation gives lower values at temperatures higher than 298 K, and the deviation increases with temperature. The deviation at higher temperatures is possibly due to the fact that the hard-sphere diameter for a single methane particle is chosen at 298 K; the mean-field prediction for W*(rc) would be even better at any temperature had one employed a temperature-dependent hardsphere diameter for a single methane particle. Further investigation is required for verifying this conjecture. The plots of W*hs(rc) and uattr(rc) in Figure 4a indicate that the decrease of W*(rc) with increasing temperature, i.e., the hydrophobic attraction being stronger with temperature, is essentially due to the temperature dependence of Whs * (rc), while F

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 6. Osmotic second virial coefficients B of methane in aqueous solutions: (a) B vs pressure for methane in pure water at 298 K; (b) B vs NaCl concentration (molality) at 298 K and 1 bar. The error bars indicate the maximum and minimum values of B as a function of the cutoff distance (see ref 13).

depending on a subtle balance of those of uattr(rc) and W*hs(rc). Figure 5b shows the variation of W*(rc)/kT with the NaCl concentration. The two components Whs * (rc)/kT and uattr(rc)/ kT both decrease with increasing salt concentration, but the concentration derivative of W*hs(rc) is greater than that of u attr (r c ). That is, the concentration dependence of −ln[ω(r)/ω2] is the main cause for the hydrophobic attraction strengthened with increasing salt concentration. The potential of mean force W(r), or the correlation function g(r), contains the full information on the effective interaction between solute molecules in a solvent, but those for nonpolar solutes in water are hardly accessible by experiment. The osmotic second virial coefficient B is a measure of the effective interaction as measured by the deviation from the van’t Hoff osmotic pressure law and is given by the integrals at the infinitedilution limit 1 lim 2 ρAβ → 0 1 = − lim 2 ρAβ → 0

B=−

reduces the attractive interaction between the hydrophobic solute molecules. The magnitudes and the pressure dependence of B are in good agreement with those at 300 K reported earlier.33 It is found that the pressure of 6000 bar has roughly the same effect on B as that of the NaCl concentration of 1 mol/kg.



CORRELATION BETWEEN THE SOLUBILITY AND THE HYDROPHOBIC ATTRACTION Now we examine correlations between the solubility and the hydrophobic attraction in the variations of the temperature, pressure, and NaCl concentration. Figure 7 shows the

∫ [e−W (r)/kT − 1] dτ ∫ [g(r) − 1] dτ

(17)

where the integration is over the whole space and dτ is the volume element. Experimental data of B for nonpolar solutes are largely unavailable but can be derived from solubility data or equations of state.37−39 Recently, the osmotic second virial coefficient B for methane in water was evaluated from the correlation-function integral, eq 17, with g(r) at short distances obtained by simulation and at large distances supplemented by an exact asymptotic form.13 It was found that B is positive at low temperatures, changes its sign at around 280 K, and is negative at higher temperatures; only at 320 K, it turns to be more negative than the gas-phase B. Ashbaugh et al. also obtained basically the same magnitudes and temperature dependence of B for methane in water using a different method to evaluate the correlation-function integral.33 Shown in Figure 6 are the pressure and NaCl concentration dependences of B. As it does with increasing temperature, B decreases with increasing pressure or increasing concentration of NaCl. These reflect increasing attractive components in W(r) with those changes of the thermodynamic state. In the metastable state of −1000 bar, B is found to be close to 0 or slightly positive, indicating that the expanded water significantly

Figure 7. Correlation between W*(rc) and 2μA*.

correlations between μ*A , a measure of the solubility, and W*(rc), a measure of the strength of the hydrophobic interactions. (The abscissa in the figure is 2μA*, instead of μA*, for the purpose of analyzing the correlations later.) In any case of the temperature, pressure, and NaCl concentration variations, W*(rc) decreases with increasing μA*; i.e., the effective attractions between the nonpolar solute molecules increase as the solubility of that solute in aqueous solution decreases. Also, W*(rc) is near linear in μ*A in the ranges of temperatures, pressures, and salt concentrations. Barbosa and Widom showed that, for the Bethe lattice model of a solution of a low-solubility solute, W*(rc) is exactly linear in μA* 10 G

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B 2 W *(rc) = − (μA* − pv0) Z

ΛT (rc) − 1 = −0.178,

(18)

Λ m(rc) − 1 = −0.202

where Z is the coordination number of the lattice and v0 is the cell volume. [W*(rc) and μ*A are W(1) and ΔG*p , respectively, in their notation.] It is plausible that W*(rc) being linear in μ*A is commonly observed for solutions of a low solubility solute. In the present study, the rate of decrease of W*(rc) with 2μA* is dependent on the mode of variation. The rate for the salt effect is the greatest, the rate for the temperature effect is nearly the same but slightly smaller, and the rate for the pressure effect is significantly smaller than the others. We can see what determines the rate dW*(rc)/d(2μA*) for each mode of variation as follows. Define Λ(r) by dW *(r ) = Λ (r ) − 1 d(2μA*)

and these reflect the ratios of the solvation thermodynamic quantities in eq 21 S*(rc) = 0.82, 2S* (μ M w )*(rc) salt

2(μsalt M w )*

= 0.80

(19)

⎛ ∂μ* ⎞ V * = ⎜⎜ A ⎟⎟ , ⎝ ∂p ⎠T

⎛ ∂μ* ⎞ (μsalt M w )* = ⎜⎜ A ⎟⎟ ⎝ ∂m ⎠T , p

(20) Figure 8. Correlation between B and μ*A .

Likewise, the corresponding solvation quantities S*(r), V*(r), and (μsaltMw)*(r) for a dimer with the A−A distance r are given as the partial derivatives of the solvation free energy μ*D(r). Then, the ratios of the solvation quantity for the dimer to twice that for a single molecule A are S*(r)/2S* = dμD*(r)/d(2μA*) for the temperature variation and so on. However, from eqs 13 and 19, Λ (r ) =

V *(rc) = 0.99, 2V *

The reason why the rate dW*(rc)/d(2μ*A ) for the pressure variation is significantly smaller than the other rates is that the volume change V*(rc) for the solvation of the dimer is smaller only by 1% than twice the volume change V* for the solvation of a single solute particle A. The ratio S*(rc)/2S* is nearly equal to (μsaltMw)*(rc)/2(μsaltMw)*, and thus, the slope dW*(rc)/d(2μA*) for the temperature variation coincides with that for the salt concentration variation. The osmotic second virial coefficient B, too, correlates with the solvation free energy μA* in the three modes of thermodynamic variation, as shown in Figure 8. In any mode

When Λ(rc) < 1, we find the commonly expected correlation between the solubility and the hydrophobic interaction, i.e., the lower the solubility the stronger the attraction; when Λ(rc) > 1, we will see the opposite correlation. Let S* and V* denote the entropy change and the volume change in the solvation process adding a solute molecule at a fixed position in the solvent under a fixed pressure.40 Note that V* is not the partial molecular volume vp but is related to it by V* = vp − kTχT. Likewise, let (μsaltMw)* be the change of the chemical potential μsalt of the salt multiplied by the mass of water Mw in the same solvation process. The solvation quantities S*, V*, and (μsaltMw)* for a single solute molecule A are identical to the partial derivatives of the solvation free energy μA* of A: ⎛ ∂μ* ⎞ S* = −⎜⎜ A ⎟⎟ , ⎝ ∂T ⎠ p

Λ p(rc) − 1 = − 0.012,

dμD*(r ) d(2μ*) A

of variation, B decreases with increasing μA*; i.e., the effective pair attractions between the nonpolar solute molecules become stronger as the solubility of that solute in aqueous solution decreases. This is consistent with the changes of W*(rc) with μA*, as seen in Figure 7. Also, we notice that B varies with μA* almost linearly in the ranges of temperatures, pressures, and salt concentrations. The rates of decrease of B with μ*A are, however, noticeably different for the three different modes of variation. The rate dB/dμA* for the salt concentration variation is by far the greatest, followed by that for the temperature variation, and that for the pressure variation is the smallest. The linear fits to B vs μA* give the slopes

(21)

for each mode of variation. Thus, Λ(r) for the temperature, pressure, and salt concentration variations reflect the ratios of the solvation quantities V *(r ) S*(r ) , , Λ p(r ) = 2V * 2S* (μ M w )*(r ) Λ m(r ) = salt 2(μsalt M w )*

Λ T (r ) =

dB (T ) = −16, dμA*

dB (p) = −4.2, dμA*

dB (m) = −93 dμA*

all in units of cm3/kJ; T, p, and m in the parentheses indicate the thermodynamic variable to be changed. Note that dB/dμA* for the salt concentration variation is much greater than that for the temperature variation, while dW*(rc)/dμ*A for the salt concentration variation is nearly the same as that for the temperature variation (see Figure 7). This is due to the fact that

(22)

where the subscripts T, p, and m denote the temperature, pressure, and concentration variations. The linear fits to the three sets of data in Figure 7 give the rates dW*(rc)/d(2μA*) = Λ(rc) − 1: H

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B B is the integral of e−W(r)/kT − 1 over the entire space; the whole part of the potential curve W(r) lowers with increasing NaCl concentration, as shown in Figure 3b, whereas W(r) decreases mainly around the solute−solute contact distance rc with increasing temperature (see Figure 3b in ref 13).

We calculate the solute−solute radial distribution functions g(r) with sufficient accuracy that the osmotic second virial coefficient B can be obtained from the integral of g(r) − 1 in eq 16. To do this, the LJ potential ϕ(r) for the methane−methane interaction was replaced by the repulsive part ϕWCA(r) of the Weeks−Chandler−Andersen (WCA) potentials23 in the MD simulation so that the methane particles are dispersed in the solvent as they are at the limit of infinite dilution. With this modification in the methane−methane pair potential, the resulting gWCA(r) is found to be independent of the solute concentration and thus we can assume that gWCA(r) is very close to the one at infinite dilution. Then, g(r), which is what we need, is given by g(r) = gWCA(r)e−[ϕ(r)−ϕWCA(r)]/kT. For the integral ∫ [g(r) − 1] dτ to converge, it is necessary that g(r) converges to 1 as r → ∞, but for a closed system of N solute particles, it converges to 1 − c/N, with c being an intensive quantity as derived by Lebowitz and Percus.43 Thus, g(r) is corrected by the factor 1 − c/N. The method for calculating B is the same as what we proposed earlier.13 We calculate the short-range contribution to B from the finite range 0 < r < Rc using g(r) as obtained from simulation and the long-range contribution to B from the range Rc < r < ∞ using an exact asymptotic expression for g(r) between particles with shortrange interactions derived by Evans et al.44 The simulation lengths for production runs at varying temperatures are 400 ns at 238 K and 100 ns at higher temperatures; those at varying pressures are 300 ns or longer; and those at varying salt concentrations are 300 ns. The solvation free energies μ*A of the LJ methane particles and −kT ln ω of hard spheres in water and in NaCl aqueous solutions at infinite dilution were obtained by the test particle insertion method. The test-particle insertion was performed for 105 or 2 × 105 equilibrium configurations, for each of which the test particle was inserted 105 times. The diameter d of the hardsphere methane particle is 3.17 Å, while the diameter di of a solvent atom i is taken to be σ of the LJ potential. Then, the diameter of the excluded sphere of the test hard-sphere particle is d + di. The solvation energies u for a single LJ methane particle and u(rc) for a pair of the particles at fixed distance rc = 3.9 Å were obtained from the MD simulations of the system containing one methane particle or one pair of the particles in the solvent. The simulation length for the production runs is 12 or 20 ns for u and 10 ns for u(rc). The solvation free energies −kT ln ω(rc) for the dimer of hard spheres were obtained by the thermodynamic integration method. The pair potentials ψλ(r) between a single solute particle and a solvent atom i are taken to be of the form λ5ψWCA(r), where ψWCA(r) is the repulsive part of the WCA potentials with σ = (d + di)/2. Then, ψλ(r) will be the hard sphere potential as λ → ∞, and therefore, −kT ln ω(rc) = ∫ ∞ 0 ⟨dΨλ/dλ⟩ dλ, with Ψλ being the sum of the pair potentials ψλ(r) over the two solute particles and solvent atoms. In practice, x = e−λ was used as a charging parameter ranging from 0 to 1; then, ψλ(r) becomes the hard sphere potential as x → 0 and it is identically 0 at x = 1. We evaluated ⟨dΨλ/dλ⟩ at 11 points: x = 0.05, 0.1, 0.2, ..., 0.8, 0.9, 0.95, then fitted a function x(1 − x)(∑8i=0 cixi) to the data, and performed the integration. Error bars (standard deviations) of the mean solvation thermodynamic quantities were obtained from the variance of the five block averages.



COMPUTATIONAL DETAILS Isobaric−isothermal MD simulations were performed for the model systems of pure water, NaCl aqueous solutions, and those containing methane, each in a cubic simulation cell under periodic boundary conditions. The potential function for the intermolecular interactions of water molecules is of the TIP4P/ 200529 and that of methane molecules is of the united-atom Lennard-Jones (LJ) with the following energy and size parameters:28,41 ϵ2 = 1.2301 kJ/mol and σ2 = 3.73 Å. The methane−water pair interaction is modeled by the LJ potential function between the center of a methane molecule and the oxygen site of a water molecule with the following LJ parameters:31 ϵ12 = 1.043 kJ/mol and σ12 = 3.4445 Å. The pair potential between Na+−Na+, Cl−−Cl−, and Na+−Cl− ion pairs is the sum of the LJ and Coulomb potentials with the parameter sets of Joung and Cheatham.30 Each of the water− ion pair potentials is the sum of the LJ potential between oxygen and ion centers and the Coulomb potentials between the partial change sites and the ion center. Each of the methane−ion pair potentials is the LJ potential between methane and ion centers. All of the LJ parameter sets are determined by the Lorentz−Berthelot combining rules except for the water−methane interactions described above. For all of the intermolecular pair potentials, the LJ part is truncated at 9 Å and the long-range Coulomb interaction is calculated by the Ewald sum with a real-space cutoff of 9 Å. The set of the model potentials employed here reproduces with great accuracy the density of pure water,29 the second virial coefficient of a methane gas, and the solubility of methane in water31 over a wide range of temperatures at 1 bar. All of the MD simulations were performed by using the program package GROMACS.42 The time step is 1 fs, the pressure−temperature control method is Nosé−Hoover, and equilibrium configurations are sampled every 50 steps during production runs after an equilibration run of 2 × 105 steps or more. The thermodynamic states (T, p, and m of NaCl) for which the MD simulations and the free energy calculations were performed are as follows: (i) T = 238, 258, 278, 298, 318, 338, 358, and 373 K at 1 bar and m = 0; (ii) p = −1000, 1, 1000, 2000, 3000, and 6000 bar at 298 K and m = 0; and (iii) m = 1, 2, and 3 mol/kg at T = 298 K and p = 1 bar. MD simulations were performed for the systems having both the solvent (water or water + NaCl) and solute (methane) species and for the systems of the solvent species alone. The former systems contain 48 methane molecules in 4000 water molecules (water + methane) with 0, 72, 144, and 216 NaCl ion pairs (m = 0, 1, 2, and 3 mol/kg). They are the systems for which g(r) of the nonpolar solute particles in aqueous solution is obtained. The latter systems consist of either 4000 water molecules (m = 0) or 2000 water molecules with 0, 36, 72, and 108 NaCl ion pairs (m = 0, 1, 2, and 3 mol/kg). Equilibrium configurations of these systems will be used for obtaining the solvation free energies μA* of the LJ methane particles and −ln ω of the hard spheres by the test-particle insertion method. The system sizes are sufficiently large and the solute concentration is sufficiently low that the resulting μA*, −ln ω, and g(r) may be considered to be those at finite dilution.



CONCLUSIONS A simple approximation to the solvation free energy and the effective pair potential between solute molecules in solution is I

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

present mean-field model, too, would fail to reproduce μ*A and W*(rc) in the thermodynamic states close to the critical point. Another subject in the present study is the correlation between the solubility of small hydrophobic molecules in aqueous solution and the strength of hydrophobic attraction. It was found that W*(rc) is a near linear function of μ*A in any mode of variation (Figure 7). The slope dW*(rc)/d(2μA*) = Λ(rc) − 1 is, however, dependent on the mode of variation and is related to the ratio of solvation thermodynamic quantities characteristic of the mode of variation as given in eq 22. It was remarked that the slopes for the temperature and NaCl concentration variations are very close to each other because the ratio S*(rc)/2S* is nearly equal to (μsaltMw)*(rc)/ 2(μsaltMw)*. The rate of change in B with μ*A is also dependent on the mode of variation. The magnitude of dB/dμA* for the NaCl concentration variation is the greatest, and that for the pressure variation is the smallest. The near linear correlation between W*(rc) and μ*A is expected to hold for hydrophobic solutes in aqueous solutions upon temperature, pressure, and salt concentration variations when the thermodynamic states are not far from ambient conditions. In general, however, the linear correlation will not continue to hold beyond a certain range of the thermodynamic states; e.g., μA* of methane turns to decrease with increasing temperature above 373 K along the vapor-pressure curve of water,45 and thus, at some temperature S* = −(∂μ*A /∂T)p = 0 and then, from eqs 19 and 22, the slope dW*(rc)/d(2μA*) turns to be infinity. In fact, such deviations from the near linearity can be seen in earlier theoretical and simulation studies.9,12

applied to understand the hydrophobic hydration and hydrophobic interactions of methane in aqueous solution with varying temperature, pressure, and salt concentration. The theoretical background is the mean-field theory of liquids combined with the potential distribution theorem initiated by Widom.18−20 In the approximation, the exponential measure Σ = e−μ*A /kT of the solvation free energy μ*A is assumed to be the product of the probability ω that the hard-sphere test particle fits into the solvent and exp(−u/kT) the Boltzmann factor of the mean-field potential energy u that the solute molecule would feel if its hard core fits among solvent molecules. The potential energy u is taken to be the average potential energy of the solute molecule A interacting with solvent molecules i via the original pair potentials ϕi(r). The mean-field μ*A is given by eq 5. Likewise, it is assumed that the solvation free energy μ*D(r) for a pair of solute molecules that are a fixed distance r apart is given as the sum of the hard-sphere contribution −kT ln ω(r) and the mean potential energy u(r). Then, the solventinduced part W*(r) of the potential of mean force is given by eqs 15 with 16. With a suitable choice of diameter d for the hard core of a solute molecule, the temperature, pressure, and salt concentration dependences of the solvation free energy μA* are accurately reproduced by the mean-field approximation, demonstrating that the underlying idea in the mean-field theory of liquids is valid for the solvation of small nonpolar solutes in aqueous solution. For each mode of change of the thermodynamic state, the mechanism of variation in μ*A or μ*A / kT is clarified. For example, the temperature dependence of μA*/kT below and around room temperature is mainly due to 1/ kT in u/kT, since u is only weakly dependent on T over the whole temperature range and −ln ω is nearly constant around Tm, the temperature of maximum density of water. However, the maximum of μ*A /kT at about 358 K, which corresponds to the solubility minimum of methane in water, is due to the fact that the rate of decrease in −ln ω is equal to the rate of increase in u/kT at that temperature. The pressure dependence of μA*/ kT and the NaCl concentration dependence of μ*A /kT at 298 K are determined essentially by the behaviors of −ln ω, while u/ kT remains nearly constant with varying pressure and NaCl concentration. Further, the temperature, pressure, and NaCl concentration dependences of the solvent-induced part W*(rc) of the effective pair potential at the contact distance are reasonably well reproduced by the mean-field approximation without any adjustment of the hard-sphere diameter. It demonstrates that the mean-field approximation holds very well for the solvation of a pair of solute molecules, too, for the evaluation of W*(rc) involves differences between two large numbers [−kT ln ω(rc) + 2kT ln ω and u(r) − 2u]. Again, the mean-field approach makes clear the mechanisms of variation in W*(rc) with increasing temperature, pressure, and salt concentration; e.g., the temperature dependence of W*(rc) is essentially due to that of Whs * (rc), the weak pressure dependence of W*(rc) results * (rc) and uattr(rc), from the opposite pressure dependences of Whs and the NaCl concentration dependence of W*(rc) comes from both of Whs * (rc) and uattr(rc). The osmotic second virial coefficients B were obtained as functions of the pressure and NaCl concentration. The results are consistent with those of W*(rc). We note that, as the mean-field approximation for liquids is not a good approximation near the critical point, the



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

K. Koga: 0000-0002-1153-5831 Present Address §

N.Y.: Faculty of Medicine, Medical School, Okayama University, Okayama 700-8558, Japan.

Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS One of us (K.K.) thanks Prof. Ben Widom for his insightful advice, fruitful collaboration, and warm hospitality at Cornell. K.K. also thanks Prof. Dor Ben-Amotz for helpful discussions of many aspects of the mean-field approximation. This work was supported by JSPS KAKENHI Grant No. JP26287099.



REFERENCES

(1) Battino, R. In IUPAC Solubility Data Series; Clever, H. L., Young, C. L., Eds.; Pergamon: Oxford, U.K., 1987; Vols. 27/28, pp 1−6. (2) Fernández Prini, R.; Crovetto, R. Evaluation of Data on Solubility of Simple Apolar Gases in Light and Heavy Water at High Temperature. J. Phys. Chem. Ref. Data 1989, 18, 1231−1243. (3) 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. (4) Kozak, J. J.; Knight, W. S.; Kauzmann, W. Solute-Solute Interactions in Aqueous Solutions. J. Chem. Phys. 1968, 48, 675−690. (5) Wood, R. H.; Thompson, P. T. Differences Between Pair and Bulk Hydrophobic Interactions. Proc. Natl. Acad. Sci. U. S. A. 1990, 87, 946−949. J

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B (6) Rankin, B. M.; Ben-Amotz, D.; van der Post, S. T.; Bakker, H. J. Contacts Between Alcohols in Water Are Random Rather than Hydrophobic. J. Phys. Chem. Lett. 2015, 6, 688−692. (7) Koga, K.; Bhimalapuram, P.; Widom, B. Correlation between Hydrophobic Attraction and the Free Energy of Hydrophobic Hydration. Mol. Phys. 2002, 100, 3795−3801. (8) Widom, B.; Bhimalapuram, P.; Koga, K. The Hydrophobic Effect. Phys. Chem. Chem. Phys. 2003, 5, 3085−3093. (9) Koga, K. Hydrophobic Effect in the Pressure-Temperature Plane. J. Chem. Phys. 2004, 121, 7304−7312. (10) Barbosa, M. A. A.; Widom, B. Molecular Correlations and Solvation in Simple Fluids. J. Chem. Phys. 2010, 132, 214506−1−8. (11) 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. (12) 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. (13) Koga, K. Osmotic Second Virial Coefficient of Methane in Water. J. Phys. Chem. B 2013, 117, 12619−12624. (14) Graziano, G. Hydrostatic Pressure Effect on Hydrophobic Hydration and Pairwise Hydrophobic Interaction of Methane. J. Chem. Phys. 2014, 140, 094503−1−8. (15) Graziano, G. Temperature Dependence of the Pairwise Association of Hard Spheres in Water. J. Phys. Soc. Jpn. 2016, 85, 024801−1−5. (16) Mochizuki, K.; Koga, K. Cononsolvency Behavior of Hydrophobes in Water + Methanol Mixtures. Phys. Chem. Chem. Phys. 2016, 18, 16188−16195. (17) Wilhelm, E.; Battino, R.; Wilcock, R. J. Low-Pressure Solubility of Gases in Liquid Water. Chem. Rev. 1977, 77, 219−262. (18) Widom, B. Some Topics in the Theory of Fluids. J. Chem. Phys. 1963, 39, 2808−2812. (19) Longuet-Higgins, H. C.; Widom, B. A Rigid Sphere Model for the Melting of Argon. Mol. Phys. 1964, 8, 549−556. (20) Widom, B. Intermolecular Forces and the Nature of the Liquid State: Liquids Reflect in Their Bulk Properties the Attractions and Repulsions of Their Constituent Molecules. Science 1967, 157, 375− 382. (21) 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. (22) Widom, B. Potential-Distribution Theory and the Statistical Mechanics of Fluids. J. Phys. Chem. 1982, 86, 869−872. (23) Chandler, D.; Weeks, J. D.; Andersen, H. C. Van der Waals Picture of Liquids, Solids, and Phase Transformations. Science 1983, 220, 787−794. (24) Weeks, J. D. Connecting Local Structure to Interface Formation: A Molecular Scale van der Waals Theory of Nonuniform Liquids. Annu. Rev. Phys. Chem. 2002, 53, 533−562. (25) Underwood, R.; Ben-Amotz, D. Length Scale Dependent OilWater Energy Fluctuations. J. Chem. Phys. 2011, 135, 201102−4. (26) Remsing, R. C.; Weeks, J. D. Dissecting Hydrophobic Hydration and Association. J. Phys. Chem. B 2013, 117, 15479−15491. (27) Rowlinson, J. S.; Widom, B. Molecular Theory of Capillarity; Clarendon Press: Oxford, U.K., 1982; p 130. (28) Guillot, B.; Guissani, Y. A Computer Simulation Study of the Temperature Dependence of the Hydrophobic Hydration. J. Chem. Phys. 1993, 99, 8075−8094. (29) Abascal, J. L. F.; Vega, C. A General Purpose Model for the Condensed Phases of Water: TIP4P/2005. J. Chem. Phys. 2005, 123, 234505−1−12. (30) Joung, I. S.; Cheatham, T. E. I. Determination of Alkali and Halide Monovalent Ion Parameters for Use in Explicitly Solvated Biomolecular Simulations. J. Phys. Chem. B 2008, 112, 9020−9041. (31) 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−1−9.

(32) Abe, K.; Sumi, T.; Koga, K. Mean-Field Approximation to the Hydrophobic Hydration in the Liquid−Vapor Interface of Water. J. Phys. Chem. B 2016, 120, 2012−2019. (33) Ashbaugh, H. S.; Weiss, K.; Williams, S. M.; Meng, B.; Surampudi, L. N. Temperature and Pressure Dependence of Methane Correlations and Osmotic Second Virial Coefficients in Water. J. Phys. Chem. B 2015, 119, 6280−6294. (34) Hummer, G.; Garde, S.; García, 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. (35) Ghosh, T.; García, A. E.; Garde, S. Molecular Dynamics Simulations of Pressure Effects on Hydrophobic Interactions. J. Am. Chem. Soc. 2001, 123, 10997−11003. (36) Ghosh, T.; Kalra, A.; Garde, S. On the Salt-Induced Stabilization of Pair and Many-body Hydrophobic Interactions. J. Phys. Chem. B 2005, 109, 642−651. (37) Widom, B.; Underwood, R. C. Second Osmotic Virial Coefficient from the 2-Component van der Waals Equation of State. J. Phys. Chem. B 2012, 116, 9492−9499. (38) Koga, K.; Holten, V.; Widom, B. Deriving Second Osmotic Virial Coefficients from Equations of State and from Experiment. J. Phys. Chem. B 2015, 119, 13391−13397. (39) Cerdeiriñ a, C. A.; Widom, B. Osmotic Second Virial Coefficients of Aqueous Solutions from Two-Component Equations of State. J. Phys. Chem. B 2016, 120, 13144−13151. (40) Ben-Naim, A. Solvation Thermodynamics; Plenum: New York, 1987. (41) Hirschfelder, J. O.; Curtiss, C. F.; Bird, R. B. Molecular Theory of Gases and Liquids; Wiley: New York, 1954. (42) 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. (43) Lebowitz, J.; Percus, J. Long-Range Correlations in a Closed System with Applications to Nonuniform Fluids. Phys. Rev. 1961, 122, 1675−1691. (44) Evans, R.; Leote de Carvalho, R. J. F.; Henderson, J. R.; Hoyle, D. C. Asymptotic Decay of Correlations in Liquids and Their Mixtures. J. Chem. Phys. 1994, 100, 591−603. (45) Ben-Amotz, D. Water-Mediated Hydrophobic Interactions. Annu. Rev. Phys. Chem. 2016, 67, 617−638.

K

DOI: 10.1021/acs.jpcb.7b12193 J. Phys. Chem. B XXXX, XXX, XXX−XXX